text:
# Probing the statistics of the temperature-density relation of the IGM

## 1 Introduction

In the last few years our understanding of small-scale structure at high redshift has increased enormously through the study of the Ly-$`\alpha `$ forest, the series of absorption features in the spectra of distant QSOs. Based in large part upon insights gained from numerical simulations (see, e.g., Cen, Miralda-Escude, Ostriker, & Rauch, 1994; Zhang, Anninos, & Norman, 1995; Hernquist, Katz, Weinberg, & Jordi, 1996) we now know that the physics governing the Ly-$`\alpha `$ forest is relatively simple. At the relevant redshifts the gas making up the intergalactic medium (IGM) is in photoionization equilibrium, which results in a tight density–temperature relation for the absorbing material with the neutral hydrogen density proportional to a power of the baryon density. Since pressure forces are sub-dominant, the neutral hydrogen density closely traces the total matter density on the scales relevant to the forest ($`0.110h^1`$Mpc). The structure in QSO absorption thus traces, in a calculable way, slight fluctuations in the matter density of the universe back along the line of sight to the QSO, with most of the Ly-$`\alpha `$ forest arising from over-densities of a few times the mean density (see, e.g., Croft, Weinberg, Katz, & Hernquist, 1998; Croft et al., 2002).

While most attention has focused on the low-order statistics of the Ly-$`\alpha `$ forest, the statistical properties of the flux are significantly non-Gaussian. Viel et al. (2004) studied the bispectrum of a large sample of QSO Ly-$`\alpha `$ flux spectrum. While their results from observations and simulations and theoretical predictions reasonably agree with each other at large scales, they found the errors on the bispectrum are too large to discriminate between models with different Ly-$`\alpha `$ forest distribution based on the bispectrum method.

However, while in general probing high-$`z`$ physics with the bispectrum is complicated, in some special cases patterns of higher order correlations can be visible in the moments of the flux. One particular pattern is that due to the action of gravity in hierarchical structure formation. The presence or absence of this pattern provides us with a powerful probe of ‘extra’ physical effects, beyond gravity, in the formation of the Ly-$`\alpha `$ forest, as first emphasized by Zaldarriaga, Seljak & Hui (2001; hereafter ZSH01). These authors noted that it was possible to constrain fluctuations in the continuum by appealing to higher order moments of the flux. In this paper we wish to address another physical effect: scatter in the density-temperature relation induced for example by different thermodynamical histories for the gas, fluctuations in the UV background radiation field or radiative transfer effects. We shall see that it is possible to put tight constraints on the scatter in the $`\rho T`$ relation by considering the higher order moments of the spectra.

The outline of the paper is as follows. In §2 we briefly discuss the simulations involved in this work and how we extract “mock” Ly-$`\alpha `$ spectra. The statistical method introduced by ZSH01 is reviewed in §3. Our main results are presented in §4, and we discuss the application to high resolution Keck spectra in section §5. In §6 we discuss several related issues.

## 2 Simulation

Since we are interested in the higher order correlations in the Ly-$`\alpha `$ flux which are induced by gravitational processes, we first need a model of gravitational clustering. For this purpose we use a particle mesh code to follow the evolution of dark matter clustering, including only gravity. Details of the PM simulation can be found in Meiksin, White, & Peacock (1999), White (1999) and Meiksin & White (2003). Given the particles’ positions and velocities, the Ly-$`\alpha `$ forest spectrum can be constructed based on the so-called “fluctuating Gunn-Peterson approximation”, in which Ly-$`\alpha `$ absorbing gas is assumed to trace the dark matter distribution, and be in a state of photoionization equilibrium with the UV background radiation (see, e.g., Cen, Miralda-Escude, Ostriker, & Rauch, 1994; Zhang, Anninos, & Norman, 1995; Hernquist, Katz, Weinberg, & Jordi, 1996; Hui & Gnedin, 1997). Under these assumptions the gas temperature can be obtained through

$$T=T_0\left(\frac{\rho }{\overline{\rho }}\right)^{\gamma 1}$$
(1)
where $`\rho `$ is the matter density and $`\overline{\rho }`$ the cosmic mean. The index $`\gamma `$ is expected to vary from near unity right after reionization, to $`1.5`$ at low-$`z`$ (Hui & Gnedin, 1997). We will set $`T_0=2\times 10^4`$K throughout, but in the next section we will introduce scatter around the values given by Eq. 1.

The neutral hydrogen number density, based on our ionization equilibrium assumption, is $`n_{\mathrm{HI}}\left(\rho /\overline{\rho }\right)^2T^{0.7}`$, with the constant of proportionality dependent on the baryon density, redshift, and hydrogen photoionization rate. We shall fix it using the observed value of the mean flux (see below). With the gas density, temperature and velocity data, the optical depth $`\tau `$, at a given velocity $`u`$ is

$`\tau `$ $`=`$ $`{\displaystyle n_{\mathrm{HI}}(x)\sigma (x)𝑑x}`$ (2)
$``$ $`{\displaystyle 𝑑x\left(\frac{\rho }{\overline{\rho }}\right)^2T^{0.7}b^1\mathrm{exp}\left[\frac{(uu_0)^2}{b^2}\right]}`$
where $`u_0`$ is the sum of Hubble expansion and peculiar velocity at position $`x`$, and $`b=\sqrt{2k_BT/m_p}`$ is the Doppler parameter. Here $`m_p`$ is the proton mass, and $`\sigma `$ is the photoionization cross section. The flux at velocity $`u`$ is then $`\mathrm{exp}(\tau )`$.

For definiteness we concentrate on the $`z=3`$ output of a $`512^3`$ particle and $`1024^3`$ mesh PM simulation in a box of side $`60h^1`$Mpc, which gives a simulation resolution of $`6\mathrm{km}\mathrm{s}^1`$ at $`z3`$. The cosmology simulated was of the $`\mathrm{\Lambda }`$CDM family, with $`\mathrm{\Omega }_{\mathrm{mat}}=0.3`$, $`h=0.67`$ and $`\sigma _8=0.9`$. The density and velocity fields were Gaussian smoothed by $`100h^1`$kpc to simulate the effects of gas pressure and $`32^2`$ spectra, each of $`1024`$ pixels, were computed on a regular grid parallel to the box axes. The normalization constant was fixed to reproduced the mean flux, $`\overline{F}=0.684`$, observed by McDonald et al. (2000) at $`z=3`$ (see also Appendix B of Meiksin & White, 2004).

## 3 Statistics

To constrain scatter in the temperature density relationship on small scales we adopt the high-order statistical method proposed by ZSH01. We briefly describe this technique here and refer the reader to their paper and Mandelbaum, McDonald, Seljak, & Cen (2003) for details.

We define the relative fluctuation of the flux, $`\delta _FF/\overline{F}1`$, and its Fourier transform

$$F(k)=𝑑xe^{ikx}\delta _F(x).$$
(3)
We bandpass filter $`\delta _F`$ with a window function $`W(k)`$: $`F_f(k)=F(k)W(k)`$, taking the window function $`W(k)`$ to be either a square window,

$$W_{k_1,k_2}(k)=\{\begin{array}{cc}1\hfill & \mathrm{if}k_1kk_2;\hfill \\ 0\hfill & \mathrm{otherwise}\hfill \end{array}$$
(4)
or a Gaussian window function. Our results are not sensitive to this choice, so in what follows we shall concentrate on the square window.

The filtered field $`F_f(k)`$ in Fourier space is then inverse transformed back into real space to define

$$\delta _H(x)=\frac{dk}{2\pi }e^{ikx}F_f(k).$$
(5)
and this is squared to produce $`h(x)\delta _H^2(x)`$. Following the conventions in Mandelbaum, McDonald, Seljak, & Cen (2003), we define $`P`$ as the power spectrum of $`\delta _F`$, $`T`$ as the power spectrum of $`h`$, and $`B`$ as the cross-spectrum between $`\delta _F`$ and $`h`$. Since $`h`$ depends quadratically on $`\delta _F`$, $`B`$ and $`T`$ can be related to the flux bispectrum and trispectrum, respectively. We also introduce the same cross-correlation coefficient $`C(k)`$ to study the effect of large scale density fluctuation on small scale power, defining $`CB/\sqrt{PT}.`$

## 4 Scatter in $`\rho T`$

We model the scatter in the temperature-density relation as log-normal in $`T`$ at fixed $`\rho `$, with mean value determined from Eq. 1. We choose the width, $`\sigma `$, to mimic the results seen in hydrodynamic simulations. Specifically we choose $`\sigma \delta \mathrm{ln}T=0.1`$ and $`\sigma =0.2`$ to study the effect. For comparison, by checking $`\rho T`$ scatter plot, we find a scatter of $`\sigma =0.1`$(0.2) is consistent with results from low (high) resolution hydrodynamic simulations in Bryan, Machacek, Anninos, & Norman (1999), respectively (see their Fig. 3). Figure 1a shows the power spectrum $`P`$ of our simulation under the assumptions of no scatter, $`\sigma =0.1`$ or $`\sigma =0.2`$. At $`k0.1\mathrm{s}\mathrm{km}^1`$ (scales $`10h^1`$Mpc), scatter in the temperature-density relation adds power.

We show in Fig. 1 the power spectrum, bispectrum, trispectrum and cross-correlation coefficient from the simulations with a top-hat filter with modes between $`0.2k0.3\mathrm{s}\mathrm{km}^1`$. We chose that range to maximize the effect of $`\rho T`$ scatter on $`C(k)`$. As pointed out by ZSH01, a positive correlation between the large scale density fluctuations and the power on small scales is induced by gravitational instability. Since, within the “fluctuating Gunn-Peterson approximation”, the fluctuations in the dark matter distribution are directly traced by fluctuations in the Ly-$`\alpha `$ forest, the positive correlation induced by gravitational instability is manifest in the Ly-$`\alpha `$ flux. In the Ly-$`\alpha `$ forest the correlation is negative since the flux decreases with increasing density.

We see that our model for scatter in the $`\rho T`$ relation does not change the large scale fluctuations of the Ly-$`\alpha `$ flux spectrum; however, by locally varying the optical depth $`\tau `$ through Eq. 2, it effectively adds noise at small scales. Figure 1d shows $`C(k)`$ at low-$`k`$ becomes steadily less negative as the scatter is increased. Note also that the bispectrum, tri-spectrum and cross-correlation coefficient go quickly to zero for $`k0.1\mathrm{s}\mathrm{km}^1`$. This is because after filtering the highest mode in the squared field, $`h(k)`$, is $`k_2k_1`$, or $`0.1\mathrm{s}\mathrm{km}^1`$ in this case.

Finally, we have investigated the requirements on signal-to-noise and resolution of spectra which would be suitable for measuring $`C(k)`$ on the appropriate scales. By adding Gaussian random pixel noise to our mock spectra we found that a $`S/N`$ approaching 50 – 100 is required in order not to wash out the small scale signal. We also notice that the level of scatter is closely related to the simulation resolution. The flux spectrum is obtained by averaging out fluctuations in simulations cells over the thermal broadening width. The smaller the cells are, the more cells will be averaging out, and the smaller the fluctuation will be. In our simulation, the typical temperature of the IGM is $`2\times 10^4`$ K, which corresponds to a size of $`18\mathrm{km}\mathrm{s}^1`$. So the flux spectrum in each pixel is roughly averaged over three neighboring cells. To test this effect, we check simulations with two other resolutions that are close to the ones in hydrodynamic simulations by Bryan, Machacek, Anninos, & Norman (1999), which are two and four times higher than the one used in the paper. We found that there are no significant differences among these cases.

## 5 An Application to Quasar Data

In this section, we study the cross-correlation coefficient $`C(k)`$ with two, high resolution, high $`S/N`$ quasar spectra, Q 1422+231 at $`z=3.62`$ (Rauch, Sargent, Barlow, & Carswell, 2001) and HS 1946+7658 at $`z=3.05`$ (Kirkman & Tytler, 1997). These two quasars were observed with Keck High Resolution Echelle Spectrometer (HIRES), with a signal-to-noise ratio ($`S/N`$) of $`140`$ and $`50100`$, respectively. We split the two spectra into a total of 13 separate pieces between Ly-$`\alpha `$ and Ly-$`\beta `$ emission to avoid contamination from Ly-$`\beta `$ absorption. Each piece has a comoving length of 60$`h^1`$ Mpc to match our PM simulation. We calculate $`C(k)`$ among these 13 samples following the method described in section §3, and compare these results against our PM simulation. We made no attempt to correct for metal lines. In order to account for evolution in the mean flux, we adjusted each spectral segment to $`F=0.684`$ by scaling the implied optical depths as in the simulations.

In Figure 2 we show $`C(k)`$ from both Keck data and our PM simulation. To study the filter effect, we select three filters: (a): $`0.02k0.06\mathrm{s}\mathrm{km}^1`$; (b): $`0.1k0.2\mathrm{s}\mathrm{km}^1`$; and (c): $`0.2k0.3\mathrm{s}\mathrm{km}^1`$. With each filter we plot the no scatter case (solid lines), scattering with $`\sigma =0.1`$ (dotted lines) and scattering with $`\sigma =0.2`$ (dashed lines). Keck data are labeled as solid circles, with error bars computed from the simulations (the sample size in the data is too small to give constraints on variances).

Clearly, with large fluctuations, negative signals present in cases with filter (a) & (b), which agrees with simulation quite well. This is consistent with what ZSH01 found in Q 1422+231 data, and confirms the correlation between the large and small scale fluctuations of the Ly-$`\alpha `$ flux spectrum induced by gravitational instability. Also in case (a) & (b) , there is almost no difference between no-scatter and scatter because their power spectra (see Figure 1a) show no difference between $`0.01`$$`0.2\mathrm{s}\mathrm{km}^1`$. At smaller scales (larger $`k`$ modes) the scattering increases the power spectrum and starts to differentiate between the scatter and no-scatter cases. We also find, with filter (c), the $`C(k)`$ from Keck data is essentially zero, consistent with $`20\%`$ scatter in $`\rho T`$ relation case. Though other effects (such as random noise) may also drive $`C(k)`$ to zero, increasing the $`\rho T`$ scatter above 20% drives $`C(k)`$ positive, so we can conservatively set an upper limit of $`20\%`$ on the the scatter in the $`\rho T`$ relation.

In the above we adopted an independent random Gaussian-distributed $`T`$ for each pixel. In reality, however, pixels tend to correlate with neighbor pixels: pixels around a high $`T`$ pixel appear to also have a high temperature and vice versa. Such a correlation among pixels tends to produce a “correlated” noise from $`\rho T`$ scatter and enhance the signals in cross-correlation coefficient $`C(k)`$: the larger the correlated pixel length scale, the stronger the $`C(k)`$. So with filter (c), we can ask: what is the largest length scale that pixels can be correlated with each other, without showing a strong signal in $`C(k)`$ that conflicts with the Keck data? We can parameterize this by adopting a correlation parameter $`r`$: $`T_i=rT_{i1}+(1r)T_{\mathrm{new}}`$, while $`T_i`$ is the temperature of pixel $`i`$, $`T_{\mathrm{new}}`$ is the temperature obtained from an independent random draw, and $`r`$ is a number between 0 and 1: $`r=0`$ means pixels are independent of each other, while $`r=1`$ mean the maximum correlation. We found, for the 20% scatter case, $`r`$ must be smaller than $`0.3`$ – 0.5 not to be ruled out by the Keck data. This roughly corresponds to a distance of $`5`$ pixels if we define two pixels correlated at the 1% level as “correlated”. This means, with a 20% scatter, the correlated scale length must be smaller than $`0.3h^1`$ Mpc, or $`30\mathrm{km}\mathrm{s}^1`$.

The observational evidence for the gravitationally induced correlation, and for small-scale decorrelation, currently rests on two high resolution spectra. It would obviously be of some interest to improve the analysis presented here and to increase the sample of QSO spectra for which the analysis has been performed.

## 6 Discussion

ZSH01 showed that gravitational instability, and the assumption of the “fluctuating Gunn-Peterson approximation”, predict a specific correlation between large- and small-scale power in the Ly-$`\alpha `$ forest. They also showed, based on high resolution Keck data, that observations were consistent with the correlation predicted by gravity. Here we point out that this observation allows one to limit the deviations from the “fluctuating Gunn-Peterson approximation” and in particular to constrain the fluctuations in the density-temperature relation of the IGM. Such fluctuations are expected to arise from varying thermodynamical histories for the gas and from fluctuations in the ultra-violet background radiation field or radiative transfer effects. We find that these fluctuations must be $`20\%`$ in $`T`$ at fixed $`\rho `$ in order to be consistent with high resolution Keck data.

Our method provides a sensitive way to probe for scatter in the temperature density relation that is expected from the complex physics of the IGM and seen in cosmological hydrodynamic simulations. Persistence of gravitationally induced higher order correlations in the spectra requires that the scatter be no more than twice as large as already seen in hydrodynamic simulations of structure formation. This limits the role of “additional” sources of fluctuations such as an inhomogeneous UV background field or radiative transfer effects. To further strengthen the case requires this correlation to be computed from more high $`S/N`$, high resolution quasar spectra such as could be obtained with Keck, Magellan, or similar telescopes. We hope that when more high quality spectra along different lines of sight become available, this method will be a powerful tool to study any “extra” physical process, beyond gravity, in the IGM.

Acknowledgments: T. Fang thanks Rupert Croft, Patrick Mcdonald, Ravi Sheth and Matias Zaldarriaga for helpful conversations. We thank M. Rauch and W. Salgent for providing Keck data of QSO 1422+231, and A. Meiksin for Keck data of HS 1946+7658. We also thank the referee for comments on our manuscript. T. Fang was supported by the NASA through Chandra Postdoctoral Fellowship Award Number PF3-40030 issued by the Chandra X-ray Observatory Center, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of the NASA under contract NAS8-39073. M. White was supported by the NSF and NASA. Parts of this work were done on the IBM-SP at the National Energy Research Scientific Computing Center.


format: text

text:
# Superstatistics, escort distributions, and applications

## 1 Introduction

In this paper we will briefly review the so-called superstatistics concept. It yields a simple and plausible argument why nonextensive methods, as developed by Tsallis and many others , are relevant for nonequilibrium systems with large-scale spatio-temporal fluctuations of an intensive parameter.

Superstatistics is a superposition of two (or even more) different statistics: One given by ordinary Boltzmann factors, and another one given by large-scale fluctuations of one (or several) intensive parameters (e.g. the inverse temperature). The corresponding stationary probability distributions arise as a convolution of the various statistics. Tsallis statistics is a special case that arises naturally in this approach for a $`\chi ^2`$-distributed parameter, if the marginal distributions are formed. But other generalized statistics are possible as well. Recently it has become clear that the concept is not only a theoretical construct but of practical physical relevance. For example, it is possible to develop superstatistical models that well describe the statistical properties of measured accelerations in Lagrangian turbulence experiments , or the statistics of defect velocities in pattern forming systems , or the measured energy distribution of cosmic rays .

## 2 Basic idea of superstatistics

Consider a driven nonequilibrium systems with spatio-temporal fluctuations of an intensive parameter $`\beta `$. These fluctuations are externally produced, by constantly putting energy into the system which is dissipated. The intensive parameter may be the inverse temperature, or an effective chemical potential, or a function of the fluctuating energy dissipation in the flow (for the turbulence application). Locally, i.e. in spatial regions (cells) where $`\beta `$ is approximately constant, the system may be described by ordinary statistical mechanics, i.e. ordinary Boltzmann factors $`e^{\beta E}`$, where $`E`$ is an effective energy in each cell. In the long-term run, the system is described by a spatio-temporal average over the fluctuating $`\beta `$.

Suppose the probability density of $`\beta `$ is $`f(\beta )`$. By integrating over $`\beta `$ one obtains marginal distributions of the form

$$p(E)=_0^{\mathrm{}}f(\beta )\frac{1}{Z(\beta )}\rho (E)e^{\beta E}𝑑\beta ,$$
(1)
where $`\rho (E)`$ is the density of states and $`Z(\beta )`$ is the normalization constant of $`\rho (E)e^{\beta E}`$ for a given $`\beta `$. Eq. (1) describes the long-term behaviour of the system (we have chosen type-B superstatistics ).

Let us consider three important examples (many more are possible).

1. $`\beta `$ is distributed according to the $`\chi ^2`$-distribution

$$f(\beta )=\frac{1}{\mathrm{\Gamma }\left(\frac{n}{2}\right)}\left\{\frac{n}{2\beta _0}\right\}^{\frac{n}{2}}\beta ^{\frac{n}{2}1}\mathrm{exp}\left\{\frac{n\beta }{2\beta _0}\right\}.$$
(2)
$`n`$ is a parameter. The average of the fluctuating $`\beta `$ is given by

$$\beta =_0^{\mathrm{}}\beta f(\beta )𝑑\beta =\beta _0$$
(3)
and the second moment of $`\beta `$ is given by

$$\frac{\beta ^2}{\beta _0^2}=1+\frac{2}{n}.$$
(4)
Let $`E=\frac{1}{2}u^2`$ be a kinetic energy. The marginal distribution (1), obtained by integration over all $`\beta `$, yields the generalized canonical distributions of nonextensive statistical mechanics

$$p(u)=\frac{1}{Z_q}\frac{1}{\left(1+\frac{1}{2}\stackrel{~}{\beta }(q1)u^2\right)^{\frac{1}{q1}}}$$
(5)
where $`Z_q`$ is a normalization constant and

$`q`$ $`=`$ $`1+{\displaystyle \frac{2}{n+1}}`$ (6)
$`\stackrel{~}{\beta }`$ $`=`$ $`{\displaystyle \frac{2}{3q}}\beta _0.`$ (7)
The distributions (5) maximize the Tsallis entropies $`S_q`$ subject to suitable constraints .

2. $`\beta `$ is log-normally distributed:

$$f(\beta )=\frac{1}{\beta s\sqrt{2\pi }}\mathrm{exp}\left\{\frac{(\mathrm{log}\frac{\beta }{m})^2}{2s^2}\right\}$$
(8)
The average $`\beta _0`$ of the above log-normal distribution is given by $`\beta _0=m\sqrt{w}`$ and the variance by $`\sigma ^2=m^2w(w1)`$, where $`w:=e^{s^2}`$. One obtains the superstatistics distribution

$$p(u)=\frac{1}{2\pi s}_0^{\mathrm{}}𝑑\beta \beta ^{1/2}\mathrm{exp}\left\{\frac{(\mathrm{log}\frac{\beta }{m})^2}{2s^2}\right\}e^{\frac{1}{2}\beta u^2}.$$
(9)
The integral cannot be evaluated in closed form, but it can be easily numerically evaluated and compared with experimental data .

3. Yet another interesting example of a superstatistics is obtained if we assume that the temperature $`T=\beta ^1`$ (rather than $`\beta `$ itself) is $`\chi ^2`$-distributed. In this case the corresponding marginal distribution $`p(u)`$ has exponential tails .

## 3 Escort distributions in superstatistical systems

For a given probability distribution $`\{p_i\}`$ of microstates $`i`$ and a given parameter $`q`$ the escort distribution $`\{P_i\}`$ is defined as

$$P_i=\frac{p_i^q}{_ip_i^q}.$$
(10)
Both types of distributions, $`\{p_i\}`$ and $`\{P_i\}`$, notoriously occur in the formulation of nonextensive statistical mechanics , like two spin-degrees of freedom which are not there in classical statistical mechanics ($`q=1`$). We will now show that within the superstatistical approach these two degrees of freedom can be easily understood.

Generally, it makes sense to $`\beta `$-average typical thermodynamic relations if the local cells are sufficiently large and if they are in local equilibrium. In this case one has for each local cell

$$S=\beta (UF),$$
(11)
where $`S`$ is the ordinary Boltzmann-Gibbs-Shannon entropy, $`U`$ is the internal energy and $`F`$ the free energy associated with each cell. The averaged version reads

$$_0^{\mathrm{}}f(\beta )S𝑑\beta =_0^{\mathrm{}}\beta f(\beta )U𝑑\beta _0^{\mathrm{}}\beta f(\beta )F𝑑\beta .$$
(12)
Note that on the right-hand side we obtain averages formed with $`\beta f(\beta )`$. For $`\chi ^2`$-distributed $`\beta `$ this is like transforming $`nn+2`$ in the distribution (2). Since

$$\frac{1}{q1}=\frac{n+1}{2}$$
(13)
and

$$\frac{q}{q1}=\frac{n+3}{2}$$
(14)
this means one is forming averages with respect to the escort distributions on the right-hand side of eq. (12), as required in the modern versions of nonextensive statistical mechanics .

Our argument is valid for arbitrary energies $`E`$. We see that within the superstatistical $`\beta `$-averaging approach entropic quantities most naturally correspond to the original distribution $`p_i`$, energetic quantities such as $`U`$ or $`F`$ most naturally to the escort distribution $`P_i`$, since there is an additional factor $`\beta `$ on the right-hand side of eq. (12).

In general, the averaged Shannon entropy on the left-hand side does not coincide with the Tsallis entropy $`S_q`$, which satisfies

$$S_q=\stackrel{~}{\beta }(U_qF_q).$$
(15)
However, it might be interesting to construct special situations where eq. (12) and (15) are equivalent.

## 4 Applications

Langevin equations where the parameters fluctuate on a large spatio-temporal scale yield a concrete dynamical realization of superstatistics . For example, to model the motion of a single test particle in a turbulent flow one considers a superstatistical extension of the Sawford model . The Sawford model describes the joint stochastic process $`(a(t),u(t),x(t))`$ of acceleration, velocity and position of a Lagrangian test particle in a turbulent flow by the stochastic differential equation

$`\dot{a}`$ $`=`$ $`(T_L^1+t_\eta ^1)aT_L^1t_\eta ^1u`$ (16)
$`+\sqrt{2\sigma _u^2(T_L^1+t_\eta ^1)T_L^1t_\eta ^1}L(t)`$
$`\dot{u}`$ $`=`$ $`a`$ (17)
$`\dot{x}`$ $`=`$ $`u,`$ (18)
where $`L(t)`$ is Gaussian white noise. Note that in this model the acceleration $`a`$ and not the velocity $`u`$ is driven by white noise, so the meaning of the variables is different as compared to an ordinary Brownian particle. $`T_L`$ and $`t_\eta `$ are two time scales, with $`T_L>>t_\eta `$ and

$`T_L`$ $`=`$ $`{\displaystyle \frac{2\sigma _u^2}{C_0ϵ}}`$ (19)
$`t_\eta `$ $`=`$ $`{\displaystyle \frac{2a_0\nu ^{1/2}}{C_0ϵ^{1/2}}}.`$ (20)
$`ϵ`$ is the energy dissipation, $`\nu `$ the kinematic viscosity, $`C_0,a_0`$ are Lagrangian structure function constants, and $`\sigma _u^2`$ is the variance of the velocity distribution. The Sawford model with constant coefficients predicts Gaussian stationary distributions for $`a`$ and $`u`$, and is thus at variance with recent measurements . However, a superstatistical generalization of the Sawford model  well fits the data.

For $`T_L\mathrm{}`$ one can derive  that the fluctuating parameter $`\beta `$ is given by

$$\beta =\frac{2a_0}{C_0^2}\nu ^{1/2}ϵ^{3/2},$$
(22)
i.e. the fluctuating energy dissipation $`ϵ`$ is proportional to $`\beta ^{2/3}`$. Fig. 1 shows the measured probability density of the acceleration of a Lagrangian test particle in a turbulent flow as obtained in the experiment of Bodenschatz et al. . Log-normal superstatistics with $`s^2=3.0`$ yields an excellent fit .

Similar superstatistical models can be formulated for pattern forming systems, just that the meaning of the variables in the Langevin equation is different. The probability density of defect velocities in inclined layer convection experiments can be quite precisely measured. As shown in Fig. 2, it quite precisely coincides with a Tsallis distribution with $`q1.46`$ . A superstatistical model well describes various stochastic properties of these defects.

Superstatistical models also have applications in high energy physics, for example as simple models to explain the measured statistics of cosmic rays .


text:
# EXPLORE/OC: A Search for Planetary Transits in the Field of NGC 2660

## 1 Introduction

As part of the EXPLORE111http://www.ciw.edu/seager/EXPLORE/explore.htm Project Mallén-Ornelas et al. (2003), we have recently begun a survey of southern open clusters (OCs) with the aim of detecting planetary transits around cluster member stars (EXPLORE/OC222http://www.ciw.edu/seager/EXPLORE/open\_clusters\_survey.htm). Probing cluster populations provides a complement to our ongoing deep monitoring studies of rich Galactic fields Mallén-Ornelas et al. (2003).

Open cluster monitoring provides the following advantages and incentives:

* In general, metallicity, age, distance, and foreground reddening are either known or may be determined for cluster members (more easily than for random field stars). Thus, planets detected around cluster stars will readily represent data points for any statistic correlating planet frequency with age or metallicity of the parent star.
* The planet-formation processes, and hence planet frequencies, may differ between the open cluster and Galactic field populations. This study allows the EXPLORE Project to compare these two different environments.
* Specific masses and radii for cluster stars may be targeted in the search by the choice of cluster and by adjusting exposure times for the target. In general, smaller stars offer better chances to detect the low-amplitude transit signal.

The difficulties and challenges involving open cluster surveys are:

* The number of monitored stars is typically lower than in rich Galactic fields, reducing the statistical chance of detecting planets.
* Determining cluster membership of stars in the open cluster fields without spectroscopic data is difficult due to the contamination by Galactic field stars. Since the clusters are typically concentrated toward the Galactic disk, this contamination may be significant.
* Significant differential reddening across the cluster field and along the line of sight can make isochrone fitting (and subsequent determination of physical parameters such as age, distance, and metallicity) difficult.

Note that both of the latter two difficulties may at least in part be circumvented by obtaining spectra. In this work, we illustrate some of the points mentioned above, describe our observing and data-reduction strategies, and show some of our preliminary results of the southern open cluster NGC 2660. Preliminary results of the study of our second target, NGC 6208, are presented in a companion paper Lee et al. (2004) in this volume.

## 2 Data

Our observations of NGC 2660 ($`l=265.9^{};b=3.01^{}`$) were obtained with the Carnegie Institution’s Swope 1m Telescope at the Las Campanas Observatory in Chile during the nights of 2003 Feb. 10–28 which were partially hampered by poor weather conditions. Our field of view was approximately 23’ $`\times `$ 15’ in size. Only $`I`$-band data with a cadence of around 7 minutes were obtained for monitoring. In that filter, a transit will be most easily distinguishable from contaminants such as grazing binaries because color-dependent limb-darkening results in steeper ingress and egress as well as a flatter eclipse bottom in $`I`$ than in other bands Seager & Mallén-Ornelas (2003). We show in Fig. 1 our observing window function which indicates the likelihood of our detecting existing transits and measuring their periods to within an aliasing factor (we require at least two transits). A necessary condition for detecting existing transits is, of course, high-precision photometry Mallén-Ornelas et al. (2003). In Fig. 2, we show our photometric precision as a function of magnitude.

## 3 Results

Around 7000 unphased photometric light curves of the stars with the best photometric precision were visually inspected to detect transit-like signals. We show in Fig. 3 three examples of light curves with the low-amplitude signature typically expected for transiting hot Jupiters. Although the eclipses visible in the light curves are most likely not due to planets, we are currently analyzing spectral data on them to determine spectral types of the parent stars as well as their variations in radial velocities. Regardless of whether or not we are able to find a planetary transit in NGC 2660, we have illustrated the potential of our methods to detect (unphased) amplitude variations of less than 1%, the most fundamental requirement to detect planetary transits.

We would like to express our gratitude to the Las Campanas Staff for their unparalleled helpfulness and dedication to optimizing every little aspect of our observing runs.


text:
# Research Problem: The Completion Number of a Graph

## 1. Introduction

A partial matrix is a matrix in which some entries are specified and the others are left as free complex variables. A completion of a partial matrix is a conventional matrix obtained by the specification of all free variables. We call a matrix $`A=\{a_{ij}\}_{i,j=1}^n`$ partial Hermitian if $`a_{ii}`$ is a specified real number for every $`i=1,2,\mathrm{},n`$, $`a_{ij}`$ is specified if and only if $`a_{ji}`$ is specified, and $`\overline{a}_{ji}=a_{ij}`$. In this paper we are only interested in Hermitian completions of Hermitian partial matrices. A partial Hermitian matrix $`A`$ is called partial positive (semidefinite) if all specified principal submatrices of $`A`$ are positive.

Let $`G=(V,E)`$ be an undirected graph in which $`V`$ is the finite set of vertices and $`E\{(x,y)|x,yV,xy\}`$ is the set of edges. We shall assume throughout that $`V=\{1,2,\mathrm{},n\}`$ for a positive integer $`n`$. A subset $`K\{1,2,\mathrm{},n\}`$ is called a clique if $`(x,y)E`$ for every $`x,yK`$, $`xy`$. The graph $`G=(\{1,2,\mathrm{},n\},E)`$ is said to be the graph of the Hermitian partial matrix $`A=\{a_{ij}\}_{i,j=1}^n`$ if $`a_{ij}`$, $`ij`$ is specified if and only if $`(i,j)E`$. If $`K`$ is a clique of $`G`$, we denote by $`A(K)`$ the (specified) principal submatrix $`\{a_{ij}\}_{i,jK}`$.

The inertia of an $`n\times n`$ Hermitian matrix $`B`$ is the triple $`(i_+(B),i_{}(B),i_0(B))`$, in which $`i_+(B)`$ (resp. $`i_{}(B)`$) is the number of positive (resp. negative) eigenvalues of $`B`$ counting multiplicities, and $`i_0(B)=ni_+(B)i_{}(B)`$ is the number of zero eigenvalues of $`B`$. The following is an immediate consequence of the so-called interlacing inequalities.

###### Lemma 1.1.

Let $`B`$ be an $`(n1)\times (n1)`$ principal submatrix of an $`n\times n`$ Hermitian matrix $`M`$. Then

$$i_\pm (B)i_\pm (M)i_\pm (B)+1.$$

Let $`A`$ be a partial Hermitian matrix with graph $`G`$ and let $`K`$ be a clique of $`G`$. By Lemma 1.1, for every Hermitian completion $`M`$ of $`A`$ we have that $`i_\pm (M)i_\pm (A(K))`$. We define

$$i(A,)=\mathrm{max}i_{}(A(K))$$

and

$$i(A;0,)=\mathrm{max}[i_0(A(K))+i_{}(A(K))],$$

in which the maximum is taken over all maximal cliques $`K`$ of $`G`$.

An undirected graph $`G`$ is called chordal if it has no minimal simple circuits with $`4`$ or more edges. The following is the main result of .

###### Theorem 1.2.

Let $`G=(\{1,2,\mathrm{},n\},E)`$ be a chordal graph, and let $`A`$ be a partial Hermitian matrix with graph $`G`$. Then:

(a) There exists a Hermitian completion $`M`$ of $`A`$ such that

$$i_{}(M)+i_0(M)=i(A;0,).$$

(b) If, in addition, $`A(K)`$ is nonsingular for every maximal clique $`K`$ of $`A`$, then there exists a nonsingular Hermitian completion $`M`$ for which $`i_{}(M)=i(A,)`$. If $`A`$ is real, then in either statement (a) or (b) M can be taken real as well.

The following well-known result in  is a corollary of Theorem 1.2 for the case when $`i(A,)=0`$

###### Theorem 1.3.

Every partial positive matrix with a chordal graph admits a positive completion.

In  it was also shown that given any nonchordal graph $`G`$ there exists a partial positive matrix $`A`$ with graph $`G`$ such that $`A`$ does not admit a positive completion.

It appears as a potentially interesting question to ask what is the minimum number of negative eigenvalues one can always obtain in case of a given nonchordal graph.

Definition. Let $`G`$ be an undirected graph. The completion number of $`G`$ is the smallest number $`k`$ such that every partial positive matrix $`A`$ with graph $`G`$ has a completion $`M`$ with $`i_{}(M)k`$.

The following is our proposed Research Problem.

Question: Given a graph $`G`$, how can one find its completion number?

More specifically, it would be of interest to find graph theoretic characterizations of those graphs with given completion number. By the results in , the completion number of a graph $`G`$ is $`0`$ if and only if $`G`$ is chordal. The completion number does not distinguish the lenght of a chordless cycle with $`n4`$ vertices. More preciselly, if $`G`$ is the chordless $`4`$ cycle with vertices $`\{1,2,3,4\}`$, by adding the edge $`13`$ to $`G`$ we obtain the graph $`G^{}`$ which has the maximal cliques $`\{1,2,3\}`$ and $`\{1,3,4\}`$ that are not cliques in $`G`$. If $`A`$ is a partial positive matrix with graph $`G`$ then we can specify the $`13`$ and $`31`$ entries to obtain a partial Hermitian matrix $`A^{}`$ such that $`A^{}(\{1,2,3\})`$ is positive while by Lemma 1.1, $`i_{}(A^{}(\{1,3,4\}))1`$. Since $`G^{}`$ is chordal, by Theorem 1.2 $`A^{}`$ has a Hermitian completion $`M`$ with $`i_{}(M)1`$. This shows that the completion number of $`G`$ is $`1`$. Similarly, the completion number of any chordless cycle with $`n4`$ vertices is $`1`$.

In fact, we can show that any positive integer is the completion number of a certain graph. First, we show this for $`n=2`$.

###### Proposition 1.4.

Let $`G_2=(V_2,E_2)`$ be the graph with $`V_2=\{1,2,\mathrm{},8\}`$ and $`E_2`$ defined by its complement as $`E_2^c=\{13,24,57,68\}`$. Then the completion number of $`G_2`$ is $`2`$.

###### Proof.

Let

$$A=\left(\begin{array}{cccccccc}1& 1& \mathrm{?}& 1& 0& 0& 0& 0\\ 1& 1& 1& \mathrm{?}& 0& 0& 0& 0\\ \mathrm{?}& 1& 1& 1& 0& 0& 0& 0\\ 1& \mathrm{?}& 1& 1& 0& 0& 0& 0\\ 0& 0& 0& 0& 1& 1& \mathrm{?}& 1\\ 0& 0& 0& 0& 1& 1& 1& \mathrm{?}\\ 0& 0& 0& 0& \mathrm{?}& 1& 1& 1\\ 0& 0& 0& 0& 1& \mathrm{?}& 1& 1\end{array}\right),$$

where ? denotes an unspecified entry. The graph of $`A`$ is $`G_2`$ and $`A`$ is partial positive. Since every Hermitian completion of both the upper-left and lower-right $`4\times 4`$ principal submatrices of $`A`$ (these matrices were first used as counterexamples in ) have at least one negative eigenvalue, it follows that the completion number of $`G_2`$ is at least $`2`$.

To prove that the completion number is at most $`2`$, let $`A`$ be a partial positive matrix with graph $`G_2`$. Consider $`G_2^{(1)}`$ be the graph obtained by adding the edge $`13`$ to $`G_2`$. The following is the list of the maximal cliques of $`G_2^{(1)}`$ that are not cliques in $`G_2`$: $`\{1,2,3,5,6\}`$, $`\{1,2,3,6,7\}`$, $`\{1,2,3,7,8\}`$, $`\{1,2,3,5,8\}`$, $`\{1,3,4,5,6\}`$, $`\{1,3,4,6,7\}`$, $`\{1,3,4,7,8\}`$, and $`\{1,3,4,5,8\}`$. By Theorem 1.3, we can specify the $`13`$ and the $`31`$ entires of $`A`$ to obtain a partial Hermitian matrix $`A_1`$ with graph $`G_2^{(1)}`$ such that $`A_1(\{1,2,3,5,6\})`$ is positive definite. By Lemma 1.1 the principal submatrices of $`A_1`$ corresponding to the next $`7`$ maximal cliques have at most one negative eigenvalue.

Let $`G_2^{(2)}`$ be the graph obtained by adding the edge $`24`$ to $`G_2^{(1)}`$. The following is the list of the maximal cliques of $`G_2^{(2)}`$ that are not cliques in $`G_2^{(1)}`$: $`\{1,2,3,4,5,6\}`$, $`\{1,2,3,4,6,7\}`$, $`\{1,2,3,4,7,8\}`$, and $`\{1,2,3,4,5,8\}`$. By Theorem 1.2 we can specify the $`24`$ and $`42`$ entries of $`A_1`$ to obtain a partial Hermitian matrix $`A_2`$ with graph $`G_2^{(2)}`$ such that $`i_{}(A(\{1,2,3,4,6,7\}))1`$. By Lemma 1.1, $`i_{}(A_2(\{1,2,3,4,5,6\}))1`$, while the principal submatrices of $`A_2`$ corresponding to the cliques $`\{1,2,3,4,7,8\}`$ and $`\{1,2,3,4,5,8\}`$ have at most two negative eigenvalues.

Let $`G_2^{(3)}`$ be the graph obtained by adding the edge $`57`$ to $`G_2^{(2)}`$. Then $`\{1,2,3,4,5,6,7\}`$ and $`\{1,2,3,4,5,7,8\}`$ the maximal cliques of $`G_2^{(3)}`$ that are not cliques in $`G_2^{(2)}`$. By Theorem 1.2, we can specify the $`57`$ and $`75`$ entries of $`A_2`$ to obtain a partial Hermitian matrix $`A_3`$ with graph $`G_2^{(3)}`$ such that $`i_{}(A_3(\{1,2,3,4,5,7,8\}))2`$. By Lemma 1.1, $`i_{}(A_3(\{1,2,3,4,5,6,7\}))2`$. Since $`G_2^{(3)}`$ is chordal, we can finally specify the $`68`$ and $`86`$ entries of $`A_3`$ to obtain a Hermitian matrix with at most two negative eigenvalues. This completes the proof of the proposition. ∎

The general case is settled by a similar construction.

###### Proposition 1.5.

Let $`G_n=(V_n,E_n)`$ be the graph with $`V_n=\{1,2,\mathrm{},4n\}`$ and $`E_n`$ defined by its complement as $`E_n=\{(4k+1,4k+3),(4k+2,4k+4):k=0,1,\mathrm{},n1\}`$. Then the completion number of $`G_n`$ is $`n`$.

###### Proof.

A straightforward generalization of the proof of Proposition 1.4. ∎

There are many other questions that can be asked in connection with the completion number. For instance, is there any connection between the completion number of a graph and its order as defined in ? Any connection to ?


text:
# Diffractive Vector Meson Production with a Large Momentum Transfer

## Acknowledgements

We thank the organizers of the Aachen EPS conference for providing the opportunity to present this work. This research was funded in part by the UK Particle Physics and Astronomy Research Council (PPARC), by the Swedish Research Council, and by the Polish Committee for Scientific Research (KBN) grant no. 5P03B 14420.


text:
# Quantum Process Tomography on vibrational states of atoms in an Optical Lattice

## Abstract

Quantum process tomography is used to fully characterize the evolution of the quantum vibrational state of atoms. Rubidium atoms are trapped in a shallow optical lattice supporting only two vibrational states, which we charcterize by reconstructing the 2x2 density matrix. Repeating this process for a complete set of inputs allows us to completely characterize both the system’s intrinsic decoherence and resonant coupling.

In recent years, there have been remarkable advances in directly controlling and observing the dynamics of individual quantum systems in a variety of domains. This degree of control of microscopic systems is one of the technological advances underlying the myriad proposals for realistic quantum-information processing systems. For instance, a number of quantum computing proposals rely on atoms trapped in an optical lattice lattice\_comp99Brennen ; Jaksch99-QC , a system in which a great deal of work has investigated coherent centre-of-mass motionWavepacket-Ertmer00 ; tunnel\_resonance\_Raithel02 ; Raithel02feedback , full characterisation of spin statesKlose-recon\_dm01 , loading of individual atoms into lattice sites patternedlattice ; Bloch-collapse02 ; Weiss03-initializelattice , and coherent interactions between atomslattcomp\_Zoller00 ; Cphase\_Deutsch02 ; Collisions03 . Here we demonstrate a technique for completely reconstructing the quantum state of motion of an atom trapped in a lattice well. By performing this density-matrix reconstruction for a complete set of input states, we are able to completely characterize the quantum evolution of the system (the “superoperator”), including decoherence.

The development of a quantum computer will rely on the success of error-correction to reduce errors to an acceptable levelQCThreshold1 ; QCThreshold2 ; QCThreshold3 ; QCThreshold4 . A quantum computer is extremely vulnerable to errors and decoherence. Experimental characterization of both the decoherence and the operations will be required in order to implement quantum error correction Shor-QEC95 ; Steane96-QEC ; Zanardi97 ; LidarDFS98 . An arbitrary operation may be characterized using quantum process tomographyChuang-Nielson01 ; Poyatos\_quantum-process97 , or QPT. The result of QPT is the superoperator, a positive, linear map from density matrices to density matrices, which governs the evolution of the density matrix for the operation. Unlike a propagator, the superoperator allows for non-unitary evolution of the system, thoroughly characterizing decoherence, relaxation and loss in a system. From the superoperator one can determine the types of errors which occur and develop procedures to reduce or eliminate errors, without requiring a priori assumptions about the underlying physical mechanism causing the errorsKnill97 .

QPT has recently been demonstrated using spins in a NMR systemChilds-Leung , the polarization of single photonsancilla-assisted\_QPT  and a singlet state filter for photon pairsMitchell03 . QPT is performed by preparing a complete set of input density matrices, subjecting each to the operation being tested, and measuring the resultant output density matrices. Due to the linearity of quantum mechanics, QPT of a process on an N dimensional system can be achieved by sending in $`N^2`$ linearly independent density matrices (alternatively, it has recently been shown that one can use a single state in a larger Hilbert space as the inputancilla-assisted\_QPT ; Tomo-operations01 .)

In this experiment we perform quantum process tomography using the motional states of atoms trapped in the potential wells of a 1-D optical lattice. We examine processes which are independent of the electronic state of the atom and are only dependent upon the motional states of the atom. The measurements are insensitive to the long-range degrees of freedom and effectively trace over the quasimomentum in the Bloch state picture (or equivalently the well index in the Wannier state picture). We use a shallow 1-D lattice which only supports 2 bound bands, which we label as ground ($`|0`$) and excited ($`|1`$). The lattice is vertically oriented, causing all atoms in higher energy, classically unbound states to quickly fall out of the lattice and become spatially separated from atoms which remain in bound states of the lattice (the Landau-Zener tunneling rates from the 3 lowest energy bands are 3$`10^7`$, 14.5 and 1150 per second in increasing order). A typical sequence, from state preparation in the lattice to measurement, lasts 20 ms.

We begin by cooling and trapping Rubidium-85 atoms in a standard vapour cell MOT to a temperature of 7 $`\mu K`$ with an rms radius of approximately 1 mm. During the optical molasses stage we turn on an optical lattice in a vertical orientation. The optical lattice is created by interfering two laser beams which are detuned 30 GHz below the Rb D2 resonance at 780.03 nm. Each beam travels through an acousto-optic modulator, each of which is driven by a function generator, providing control of the relative phase between the lattice beams. By modulating the phase of one of the lasers the lattice may be displaced by up to a lattice spacing within one microsecond. The beams are superposed after the acousto-optic modulators, and co-propagate to the vacuum chamber with orthogonal polarizations, to reduce phase fluctuations between the beams which would impress noise on the lattice. The beams are separated on a polarizing beam splitter near the MOT. The polarization of one beam is then rotated such that the beams have parallel polarizations in the MOT. The beams have an angle of 50 degrees between them, creating an optical lattice with a lattice constant of $`L=0.93`$ microns. The depth of the lattice is controlled by the intensity and detuning of the beams, and is chosen to be $`18`$ $`E_r`$ (where $`E_r=h^2/8L^2m=h690`$ Hz is the effective recoil energy of the lattice) at which depth it contains two bound states. The energy separation between the states is $`\mathrm{}2\pi 5.0`$ kHz. The scattering rate from the lattice beams is on the order of 4 Hz, which is insignificant on the timescale of the experiment.

The population in each band of the lattice is determined by adiabatically decreasing the lattice potentialJessen-adiab95 ; Mlynek97 . In order to satisfy the adiabatic criterion we must decrease the potential slower than $`h`$ $`10^8/s^2`$, whereas the fastest decrease we use is $`h4.110^6/s^2`$ with non-adiabatic effects appearing if the turn-off is faster than $`h1.410^7/s^2`$. As the depth of the wells decreases the energy bands gradually get closer to the top of the potential. Once an energy band becomes classically unbound, then all the atoms in that band accelerate downwards due to gravity. Since each state becomes unbound at a different time, each band becomes mapped into a different location in space. The spatial distribution can be recorded on a CCD camera by fluorescence imaging. Figure 1a shows a sample spatial distribution after ramping down the potential over a time of 45 milliseconds. Alternatively, the lattice depth may be quickly, but still adiabatically, lowered to a depth of 9 $`E_r`$. At this depth the ground state has a Landau-Zener-limited lifetime of 250 milliseconds while the excited state has a lifetime of 0.5 milliseconds. The excited state atoms quickly escape the lattice while ground state atoms remain trapped. After holding the depth constant for some time (typically about 20 ms), the relative populations can be determined by fluorescence imaging, as shown in Figure 1b. The beams creating the lattice have a Gaussian shape, with rms width of 3 mm, causing the lattice to be shallower as we move farther away from the center. The spatial variation of the lattice depth causes atoms at the edge of the lattice to tunnel out sooner than atoms near the center, resulting in the curved clouds seen in Figure 1a. To reduce broadening effects due to inhomogeneous well depths we integrate the flourescence signal only over the central 600 $`\mu `$m of the cloud where the potential varies slowly.

A sample of ground state atoms is prepared by filtering out the excited state atoms from the lattice. This is accomplished by reducing the well depth to a depth of approximately 9 recoil energies, at which point only the ground state is bound. The potential depth is held at this value for 3 ms, long enough for almost all the excited state atoms to escape. The well depth is then adiabatically increased back to the original depth, preparing a sample of atoms with up to 95% occupation of the ground state.

To prepare a variety of initial states, we make use of our ability to displace the lattice, and of the atoms’ free evolution. Displacement of the lattice is equivalent to a spatial translation of the atom cloud in the lattice’s reference frame, constituting a coherent coupling between the energy eigenstates. In addition to coupling the two bound states, this induces some transitions to unbound states, which are lost from the lattice. Spatial translations change the coefficients of the states as described by the following equations.

(1)
The coefficients $`c_{00},c_{10}`$ and $`c_{11}`$ are determined by displacing the lattice with different initial state populations. A typical displacement used during the measurement is 116 nm. For this displacement we measure $`c_{00}`$,$`c_{10}`$ and $`c_{11}`$ to be 0.86(2), 0.50(2) and 0.53(10) respectively, close to the theoretical values of 0.87,0.45 and 0.63. During a time period $`\mathrm{\Delta }t`$ of free evolution in the lattice, the ground and excited states acquire a relative phase shift of $`\omega \mathrm{\Delta }t`$, where $`\omega `$ is the oscillation frequency in the lattice. Using a combination of displacement and time delay we can prepare superposition states with arbitrary relative phase.

State tomographyBigelow-staterecon  is performed by projecting the unknown state onto a set of known states. We use a set of non-orthogonal states $`\{\mathrm{\Phi }_1,\mathrm{},\mathrm{\Phi }_4\}=`$ $`\{|0,|1,|\theta _x,|\theta _y\}`$ whereThew-Qudit02  $`|\theta _x=\mathrm{cos}\theta |0+\mathrm{sin}\theta |1`$ and $`|\theta _y=\mathrm{cos}\theta |0i\mathrm{sin}\theta |1.`$ We project onto states of the form $`\mathrm{cos}\theta |0+\mathrm{sin}\theta |1`$ by spatially displacing the trapping potential before separating the resulting energy eigenstates. We choose our displacements to be small in order to have negligible coupling to the higher energy, unbound states. We use a displacement of L/8=116 nm, generating states with an experimental value of $`\theta 0.5`$ radians. The state $`|\theta _x`$ can be changed into state $`|\theta _y`$ by a time delay of a quarter period after displacement. State tomography is performed by measuring the projections $`m_i=\mathrm{\Phi }_i\left|\rho \right|\mathrm{\Phi }_i`$ for all states $`|\mathrm{\Phi }_i`$. The resulting measurements $`\{m_1,\mathrm{},m_4\}`$ are used to reconstruct the density matrix.

$$\rho =\left(\begin{array}{cc}m_1& \frac{(m3+im4)m_1\mathrm{cos}^2\theta m_2\mathrm{sin}^2\theta }{2\mathrm{sin}\theta \mathrm{cos}\theta }\\ \frac{(m3im4)m_1\mathrm{cos}^2\theta m_2\mathrm{sin}^2\theta }{2\mathrm{sin}\theta \mathrm{cos}\theta }& m_2\end{array}\right)$$
(2)
Process tomography is first performed on the free evolution of a quantum state in the lattice for one period. This allows us to completely characterize decoherence intrinsic to the lattice. To perform process tomography we prepare a complete set of density matrices as input states. The four linearly independent density matrices we used correspond to a ground state, prepared by filtering out the excited state population; a coherent state with real coherence prepared by displacing a ground state; a coherent state with imaginary coherence prepared by adding a quarter-period time delay after displacement; and a mixed state. The mixed state can be prepared by either skipping the filtering step, or by preparing a coherent state and waiting 3 ms for it to decohere (see discussion below).

QPT proceeds as follows: one of the input density matrices is prepared, and characterized with state tomography as outlined above; the same state is again prepared and allowed to freely evolve in the lattice for 200$`\mu s`$ (one oscillation period), and then state tomography is performed on the resulting state. Figure 2 shows the projection of each input density matrix onto the projection states, before and after the ‘operation’.

The super-operator $``$ resulting from QPT can be expressed in a number of ways. One common form is the operator sum representation,

$$\rho _{out}=\left(\rho _{in}\right)=\underset{i}{}\widehat{A}_i\rho _{in}\widehat{A}_i^{}$$
(3)
where $`\widehat{A}_i`$ are operational elements, often called Kraus operatorsKraus\_Ops , subject to the constraint $`\underset{𝑖}{}\widehat{A}_i^{}\widehat{A}_i=I`$ . The Choi matrixChoi75 ; Choi\_Proof-Leung  provides a straightforward procedure to obtain experimental Kraus operators. The Choi matrix is defined as

$$\underset{¯}{C}=\underset{i,j}{}|ij|(|ij|)$$
(4)
where $`|ij|`$ is an outer product of basis states and $`(|ij|)`$ is the super-operator acting on the matrix given by the outer product $`|ij|`$. Then

$$\left(\rho \right)=\underset{i,j}{}C_{i,j}\rho _{i,j}$$

where $`C_{i,j}=(|ij|)`$ is the i,jth 2$`\times `$2 submatrix of $`C`$ and $`\rho _{i,j}`$ is the i,jth element of the density matrix. The eigenvalues and eigenvectors of $`C`$ can then be used to determine the canonical Kraus operators, given by $`\widehat{A}_i=\sqrt{\kappa _i}\widehat{k_i}`$ where $`\kappa _i`$ is the ith eigenvalue and $`\widehat{k_i}`$ is the corresponding eigenvector written in matrix form. The matrix $`|ij|`$ is not necessarily a density matrix, but can be written as a linear combination of the measured input density matrices. Using a maximum-likelihood technique, we find the Choi matrix which best predicts the measured output states given the measured input states. This search is limited to physical, i.e., completely positive, Choi matrices. The resulting Choi matrix is shown in Figure 3. From the Choi matrix we find that the populations are preserved to within experimental uncertainty while the coherences decay by 36 percent.

The same data may be visually displayed using a Bloch sphere representationQuantuComp\_Nielson ; ancilla-assisted\_QPT , which has the advantage of showing how any state on the surface of the Bloch sphere evolves into a new state. Figure 4a shows the initial, undisturbed Bloch sphere before evolution in the lattice, and Figure 4b shows the Bloch sphere after one oscillation. The sphere becomes prolate, by contracting toward the z-axis by 36%, as expected. The Kraus operators are determined from the Choi matrix. The most significant Kraus operators are found to be $`\widehat{A}_1=0.90\widehat{I}+\widehat{R}_1`$ and $`\widehat{A}_2=0.41\widehat{\sigma }_z+\widehat{R}_2`$ where $`\widehat{I}`$ is the identity matrix, $`\widehat{\sigma }_z`$ is the z Pauli matrix and $`\widehat{R}_1,\widehat{R}_2`$ are small remainders with magnitudes bound by $`Tr[\widehat{R}_i^{}\widehat{R}_i]`$ $`0.03`$. The other two Kraus operators are insignificant on the scale of our experimental resolution, also satisfying a similar bound. Kraus operators of the form $`\widehat{I}`$ and $`\widehat{\sigma }_z`$ are consistent with pure dephasing, as expected for from either inter-well tunneling or inhomogeneous broadening.

The coherence time is thus found to be 555 $`\mu s`$ (2.78 periods); this is shorter than the coherence time of 2 ms expected based on the width of the excited band (the inter-well tunneling rate) and on the variation in well depth across the finite Gaussian lattice beams. The number of oscillations is however consistent with oscillations observed in other worktunnel\_resonance\_Raithel02 ; Haycock00 ; Raithel98 . It is believed that the observed decoherence is caused by small-scale inhomogeneity, such as fringes, in the lattice beams. The lattice beams have been spatially filtered before travelling to the MOT. Unlike earlier experiments, we use a two-level system, and anharmonicity is therefore not a factor in the observed decoherence.

An example of an operation necessary for quantum information processing is a single-qubit rotation. To demonstrate the applicability of process tomography to characterizing such operations we attempt to perform a rotation of the Bloch sphere by a method equivalent to a Rabi oscillation. We sinusoidally drive the displacement of the lattice at the trap frequency, thereby coupling neighbouring states coherently. Process tomography is performed after a single period of this drive. The lattice displacement is driven by applying a sinusoidal phase modulation with an amplitude of $`\pi /9`$ radians to one of the lattice beams. The displacement is kept small to ensure that coupling is predominantly between neighbouring states. We test both a sine drive, $`\mathrm{\Delta }x(t)=x_m\mathrm{sin}(\omega _0t)`$, and a cosine drive, $`\mathrm{\Delta }x(t)=x_m(\mathrm{cos}(\omega _0t)1)`$, where $`x_m=`$ $`26`$ $`nm`$ and $`\omega _0`$ is the oscillation frequency in the lattice and the pulse lasts from $`t=0`$ to $`t=2\pi /\omega _0`$. We again find the Choi matrix from a maximum-likelihood model, but find the Bloch sphere representation to be the most intuitive. Figures 4c and 4d show the resulting Bloch spheres, which have rotations of 35.5 degrees and 36.4 degrees about the y-axis and x-axis respectively. The Bloch spheres are rotated 90 degrees out of phase with one another as expected for driving fields 90 degrees out of phase. The resulting shape is close to a sphere, but the radius has decreased in all dimensions. In particular, the length of the semi-minor axis’ for the sinusoidal drive is 0.69 (while it should be noted that in the absence of the coherent drive, it decayed to a value of 0.64). A simulation using a truncated harmonic-oscillator model predicts a rotation of 35.0 degrees about the y-axis.

We have introduced a new technique for determining the motional quantum state of atoms in an optical lattice, and applied it to a demonstration of quantum process tomography. We have extracted the “superoperators” fully characterizing the action of several different operations on an arbitrary state of atoms in the lattice, specifically, free evolution for one period and two different resonant-coupling protocols. In this way, we have characterized the intrinsic dephasing of the system over time, and the effectiveness of single-qubit rotations induced by resonant modulation of the lattice phase. We plan to extend these techniques to test the Markovian approximation; to characterize and optimize bang-bang methodsLloyd-bangbang98 ; LidarBB  for removing inhomogeneous-broadening effects; and to study the well-depth-dependence of the decoherence, investigating the role of inter-well tunneling, Wannier-Stark transitions, and Bloch oscillations. The procedure can be extended to higher-dimensional Hilbert spaces, although the number of measurements required grows exponentially with the dimensionality. Process tomography should prove essential for tailoring error-correction protocols to the observed behaviour of particular physical realisations of quantum-information systemsCory-incoher03 ; Cory-incoher03b ; it is reasonable to expect that such system-by-system tailoring may prove necessary if error thresholds on the order of $`10^4`$ or $`10^5`$ are ever to be reachedQCThreshold1 ; QCThreshold2 ; QCThreshold3 ; QCThreshold4 . More generally, it is the only method to permit complete characterisation of the evolution of open quantum systems, and should play a central role in the toolbox for control and study of individual quantum systems, whether for quantum-information processing or for other applications.

We would like to thank Daniel Lidar and Samansa Maneshi for helpful discussions, and we thank Sara Schneider for assistance with the maximum-likelihood methods. We would also like to thank the DARPA-QuIST program (managed by AFOSR under agreement No. F49620-01-1-0468), the National Science and Engineering Research Council of Canada, the Canadian Institute for Photonics Inovation, and eMPOWR for financial support of this project.

Figures:

Figure 1) a) A fluorescence image of the state populations in a lattice obtained by adiabatic decrease of the lattice potential is shown . A stepwise decrease of the potential leads to a clearer separation of the states as shown in b).

Figure 2) Matrix of measured projections. Input density matrices (reading left to right: ground state; mixed state; superposition with real coherence and superposition with imaginary coherence) are shown along the top while the post-selected states are listed on the side. The table on the left shows the projections for the input states while the table on the right shows the corresponding projections after one period. Note the decreased contrast in the $`\theta _x`$ and $`\theta _y`$ projections for the coherent states. All populations are unchanged to within experimental error.

Figure 3) The Choi matrix after one oscillation in the lattice characterizing decoherence. The left graph shows the real part of the Choi matrix and the right graph shows the imaginary. The matrix is dominated by real components at the corners. The diagonal corners represent the mapping of populations into populations. The off-diagonal corners, which map coherences into coherences, are significantly less than one, showing a loss of coherence. The dotted lines separates the $`2\times 2`$ submatrices of the Choi matrix.

Figure 4) Bloch sphere representation of process tomography. a) The initial Bloch sphere representing the space of pure input states. b) after 1 period of free evolution the sphere contracts horizontally due to a loss of coherence. c) B.S. after sine drive showing rotation about the y-axis. d) B.S. after cosine drive, showing rotation about the x-axis.