text:
# Discovery of Pulsed X-rays from the SMC Transient RX J0052.1-7319

## 1. Introduction

RX J0052.1-7319 is an X-ray source located in the the Small Magellanic Cloud, first detected with Einstein (1E 0050.3-7335, Wang & Wu  (1992)). It was classified by Kahabka and Pietsch (1996) as a transient X-ray binary candidate based on the ROSAT PSPC data from October 1991 and April 1992.

As part of a systematic search of the ROSAT data for pulsed sources, we have discovered pulsations from RX J0052.1-7319 at a frequency of 65.4 mHz (period = 15.3 s) using ROSAT HRI observations from 1996 November-December. Our search had previously detected another SMC pulsar, J0117.6-7330 (Macomb et al.  (1999)).

After a preliminary report of our discovery of the pulsar nature of RX J0052.1-7319 (Lamb et al.  (1999)), Kahabka estimated a new source position based on HRI observations of 1995 May (Kahabka 1999a ), and confirmed the detection of pulsations using an HRI observations from October 1996 (Kahabka & Pietsch  (1996)). Israel et al. (1999) identified a Be star as the likely optical counterpart to the source. Udalski (1999) then reported on the long term optical variability of this candidate, based on Optical Gravtational Lensing Experiment (OGLE) monitoring.

Here we report on the pulsations we detected in the ROSAT HRI data, detection of these pulses in hard X-rays using simultaneous BATSE data, the history of the outburst visible in the BATSE data, and optical observations of the candidate companion and nearby stars.

## 2. Observations

We initially discovered pulsations from RX J0052.1-7319 in observations made with the High Resolution Imager (HRI) detector of the ROSAT (Trümper  (1983)) in 1996 November and December. The HRI, which consisted of two cascaded microchannel plates (MCPs) with a crossed grid position readout system, was sensitive to X-rays in the $``$ 0.2 – 2 keV range. Thereafter we detected the pulsations in observations from the Burst and Transient Source Experiment (BATSE) (Fishman  (1989)) on the CGRO using data from the Large Area Detectors (LADs) which are NaI scintillation detectors sensitive to hard X-rays/soft gamma-rays in the 20 keV to 2 MeV range.

### 2.1. ROSAT Observations

Our discovery was made in the ROSAT HRI observations of 1996 November 10.71 – December 9.16 as part of a systematic search of the ROSAT data for previously undetected pulsars. These data were among 1365 data sets we selected from the catalog of 59911 ROSAT/HRI point source observations based on the potential for detection of significant pulsations. Each data set was processed using standard Ftools, with barycentered arrival times for the first 200 ks of each observation binned with 5 ms resolution, and Fourier transformed. The power spectra were searched for significant pulsed signals unassociated with the $``$5760 s spacecraft orbit or the 402 s period spacecraft orientation wobble. The RX J0052.1-7319 HRI observations were selected for further analysis because of a strong signal near 65 mHz.

These HRI observations, rh600811n00 and rh600812a01, are of two SMC fields which contain RX J0052.1-7319. The observations for these two fields are partially contemporaneous. Figure 1 shows the count rates in 1000 s bins for events within $`60\mathrm{}`$ of the source in these observations. The roughly 10-20% difference in counting rates between contemporaneous portions of the two observations may be due principally to differences in source vignetting. For observation rh600812a01, the source is $`17\mathrm{}`$ from the center of the HRI field and vignetting corrections are 10-20%; for observation rh600811n00 the source is $`8\stackrel{}{\mathrm{.}}4`$ from the HRI center and vignetting is substantially less. From observation rh600811n00 it is clear that the source flux is declining.

As can be seen in Figure 1, the initial 200 ks (2.3 days) interval of the observations, which was used for the initial search fft, contained only a portion of either of the two observations. We chose for further analysis the data of Nov 11.21-16.83 (all of rh600812a01 and the first portion of rh600811n00).

From initial analyses of this dataset we knew pulse frequency was changing significantly, and that there was substantial power at higher harmonics of the pulse frequency. To accurately estimate the pulse frequency and frequency derivitive we maximized the $`Z_3^2`$ statistic (Buccheri et al.  (1983)), where $`Z_3^2(f)=p(f)+p(2f)+p(3f)`$, and the Rayleigh statistic $`p(f)`$ is calculated as:

$$p(f)=\frac{2}{N}\left|\underset{k=1}{\overset{N}{}}exp(i2\pi f[t_k+\frac{1}{2}\alpha (t_k\tau )^2])\right|^2.$$
(1)
Here $`N`$ is the number of photons detected in the interval, $`f`$ is the analysis frequency, $`\alpha =\dot{\nu }/\nu `$, with $`\nu `$ the pulse frequency, $`t_k`$ is a barycentric arrival time of photon $`k`$, and $`\tau `$ is an epoch within the data interval.

For epoch MJD 50401.0 (1996 Nov. 14.0) we obtained $`\nu `$ = 0.06545850(7) Hz and $`\dot{\nu }=5.48(7)\times 10^{11}\mathrm{Hz}\mathrm{s}^1`$. The maximum value of $`Z_3^2`$ was 524, which is highly improbable by chance.111We estimate that the probability of obtaining a $`Z_3^2`$ of 524 or more due to Poisson noise is less than $`10^{98}`$, including the number of trials introduced by searching in frequency and frequency derivative. This calculation however neglects systematic signatures in the data, which would likely play a dominate role any false detection with a $`Z_3^2`$ this large. The Rayleigh statistic for this ratio of $`\dot{\nu }/\nu `$ is shown in figure 2.1. In addition to the pulse fundamental, four pulse harmonic overtones are clearly detected.

The HRI data from the 1996 November 11.21-16.83 interval epoch-folded with this ephemeris is shown in figure 2.1. The pulse profile has a single asymmetric peak, with a narrow valley at minimum. The pulsed fraction (\[mean-minimum\]/mean) is $`35\pm 3`$%.

For the HRI data from the interval 1996 December 5-9 we obtain a frequency $`\nu =0.0655330(3)`$ Hz at epoch MJD 50524.0, and $`\dot{\nu }=2.3(2)\times 10^{11}\mathrm{Hz}\mathrm{s}^1`$. Thus the decrease in counting rate noted above was accompanied by a decrease in spin-up rate, as would be expected from the intrinic correlation between mass accretion rate and angular momentum accretion rate in a disk fed accreting pulsar.

For the 1996 November 11-16 and December 5-9 intervals we find source count rates of 0.92 and 0.46 counts s-1 respectively. These rates are corrected for vignetting (David et al.  (1999)) and the $`15\mathrm{}`$ radius circle used to select the data, using the HRI off-axis encircled energy profile (Boese  (2000)).

To estimate fluxes from these count rates we have used the spectrum estimated for the source by Kahabka & Prietsch (1996) based on a fit of ROSAT PSPC data. This consists of thermal bremsstrahlung model with $`kT=17`$ keV, a galactic absorption column of $`3\times 10^{20}\mathrm{cm}^2`$, and an absorption column within the SMC with density of $`N_H=7.0\pm 4.0\times 10^{21}\mathrm{cm}^2`$ (95% confidence errors), with metallicities reduced by a factor of seven from solar system abundances. This model results in 0.1-2.0 keV unabsorbed flux estimates of $`5.9\pm 1.0\times 10^{11}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ and of $`3.0\pm 0.5\times 10^{11}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ for the 1996 November 11-16 and December 5-9 intervals respectively.

Using the 1996 November data of observation rh600811n00, we have determined the position of the source to be R. A. = $`0^\mathrm{h}52^\mathrm{m}13\stackrel{\mathrm{s}}{\mathrm{.}}65`$, Dec. = $`73\mathrm{°}19\mathrm{}19\stackrel{}{\mathrm{.}}5`$, with a statistical uncertainty ($`1\mathrm{}`$) which is negligible in comparison to systematic position errors, which result in an 68% confidence error circle of radius 7$`\mathrm{}`$ (Kürster  (1993), ROSAT Users Handbook222 http://heasearc.gsfc.nasa.gov/docs/rosat/ruh/handbook/node34.html). These systematic errors are due to uncertainty in the aspect of ROSAT. We note that this position differs by nearly 10” from the position quoted by Kahabka (1999b) and Kahabka (2000). However, since these latter determinations used either very low rate data ($`<`$0.01 c/s) or data in which the source was at the edge of the HRI field of view, we believe that the position given above may be the more accurate.

### 2.2. BATSE Observations

The Burst and Transient Source Experiment (BATSE) was an all sky monitor which flew onboard the Compton Gamma-Ray Observatory (CGRO). The data used in the analysis presented here were from the Large Area Detectors (LADs) which were NaI(Tl) scintillation counters 1.27 cm thick, with 2025 cm2 area, that were located at the eight corner of the CGRO spacecraft. Count rate data from the LADs was continuously available from 1991 April to 2000 May, except during SAA passages and occasional telemetry outages.

The 20-50 keV channel of the BATSE LAD discriminator rates (DISCLA channel 1) were analysed for pulsations from RX J0052.1-7319 using techniques discussed in Finger et al. (1999). Rates were combined from different detectors with coefficients optimal for a source with a spectrum of the form $`dN/dE=Aexp(E/kT)/E`$ with $`kT=20`$ keV. A large number of pulse profiles are obtained by fitting short segments of these combined rates with a model consisting of a quadratic spline background, plus a low order Fourier expansion pulse profile model. These profiles are then combined over multiple day intervals using trial frequencies and frequency rates. The resulting combined profiles are evaluated with the $`Y_n`$ statistic (Finger et al.  (1999)), which is simular to the $`Z_n^2`$ statistic, but accounts for possible non-Poisson noise.

Using data from the interval 1996 Nov 11.0–17.0, $`Y_3`$ was searched over a frequency range of width $`10^4`$ Hz centered on the ROSAT measurement, and a frequency rate range of $`10^{11}\mathrm{Hz}\mathrm{s}^1`$ to $`10^{10}\mathrm{Hz}\mathrm{s}^1`$ resulting in maximum value of $`Y_3=57.0`$ which we estimate333This estimate approximates the distribution of $`Y_3`$ with the chi-squared distribution with 6 degrees of freedom, which is valid in this case because of the large number of profiles being combined, and accounts for 1100 independent trials. has a probability of being exceeded by chance of less than $`2\times 10^7`$. The resulting frequency and frequency rate estimates were $`\nu =0.06545830(10)`$ Hz at epoch MJD 50401.0, and $`\dot{\nu }=5.62(12)\times 10^{11}`$ Hz s-1, in good agreement with the ROSAT results.

The BATSE pulse profile for this interval using the ROSAT HRI ephemeris is shown in figure 2.2. The profile is shown relative to the mean flux level, which cannot be determined from the data due to the weakness of the source and the high background level. The solid curve is the profile corresponding to six Fourier coefficients which have been estimated with the same fitting technique used with the frequency search. Error bars for the value of the curve at 13 approximately independent points are shown. Six harmonics were chosen because this approximately matches the 1.024 s resolution of the data.

Pulsations were then detected in nine additions six-day intervals, extending the total period the source was detected with BATSE to 60 days between 1996 September 18 and November 17. Beginning with the Nov 11-17 interval, the total period of detection was progressively extended outward by frequency and frequency rate searches in the adjacent six day intervals, using $`5\times 10^5Hz`$ frequency ranges centered on the extrapolation from the neighboring frequency and frequency rate estimate, and the frequency rate range given above. Searching in this manner at the boundaries of the known frequency history reduces the required frequency search range, providing the best sensitivity. Figure 2.3 shows the resulting measurements of frequency, frequency rate, and r. m. s. pulsed flux in the 20-50 keV energy band. The figure also shows for comparison the frequency and frequency rate measurements from the ROSAT HRI data.

Efforts were made to detect the source from 1996 August 1 to September 18, and from 1996 November 17 to 1997 January 8, but these failed. These attempts used searches with a frequency range of $`3\times 10^4`$Hz centered on the nearest detection, and a frequency rate range of $`10^{11}\mathrm{Hz}\mathrm{s}^1`$ to $`10^{10}\mathrm{Hz}\mathrm{s}^1`$. Prior to 1996 September 18, pulsations from the 130.4 mHz pulsar 4U 1626-67 interfered with first harmonic (i.e. $`2\times f`$) contribution to $`Y_3`$, in a narrow frequency band near 65.2 mHz. For this narrow band, this harmonic was left out of the search. We conclude that in these intervals prior to and following our 1996 September 18 to November 17 detections that the source had pulsed flux below our sensitivity level. We estimate an upper limit of 8 mCrab for the 20-50 keV r. m. s. pulsed flux for these intervals.

### 2.3. Optical Observations

In figure 2.3 we show a V band image taken from South African Astronomical Observatory (SAAO) 1.0m telescope on 20 Jan 1999. Marked on the figure is the X-ray error circle of Kahabka (2000) which has a 6” radius, and the error circle from this work, with a 7” radius. The object marked A is the source identified as the counterpart by Israel et al. (1999) based upon a red spectrum. The objects labeled A, B, C, and D are also identified in Figure 3.

If the optical counterpart of RX J0052.1-7319 is a Be star, then it should show excess $`H_\alpha `$ emission. In figure 3 we show $`H_\alpha `$ versus R band counts for 25 stars in the vicinity of Kahabka’s X-ray error circle. The line is a linear best fit through all the data points shown. Stars that lie above the line are ones that show an $`H_\alpha `$ excess – i.e. objects A, B and C as identified on figure 2.3. Object D, the one inside the Kahabka (2000) error circle, shows no excess. Objects B and C, though exhibiting an $`H_\alpha `$ excess, are too far away from the X-ray position to be the counterpart. If, as we might expect, the optical companion is a Be star, the only reasonable counterpart is object A, the star identified by Isreal et al. (1999). Object A is well contained with the error circle for RX J0052.1-7319 presented here, and thus the identification of the X-ray source as a member of a Be system seems reasonable.

The absolute V magnitude of object A may be determined from the value quoted by Udalski (1999) of $`m_v`$=14.67. The distance modulus to the SMC determined by Westerlund (1997) is $`(mM)_o`$ = 18.9. In addition, the average extinction is E(B-V) = 0.07-0.09 (Schwering & Israel  (1991)), though there are regions in the SMC where it can rise as high as 0.25. Combining these parameters leads to an estimate of $`M_v`$= -4.47$`\pm `$0.02. This magnitude corresponds to a star in the range B1III - B0V, very similar to the value obtained for Be/X-ray binary counterparts (e.g. that of RX J0117.6-7330 quoted in Coe et al. (1998)).

## 3. Discussion

The ROSAT HRI observations we have presented show that RX J0052.1-7319 is a pulsar with a 65 mHz rotation frequency. Variations in the flux by more than two orders of magnitude between ROSAT observations at different epochs demonstrates that the source is a transient (Kahabka & Pietsch  (1996),Lamb et al.  (1999)). The BATSE observations we have presented show that the November – December 1996 ROSAT HRI observations occurred on the tail of a outburst that began at least two months earlier, lasting more than 80 days.

The companions of transient pulsars with known spectral type are generally Be (or Oe) stars. There are now known a few transient pulsars in low-mass X-ray binaries, but these all have high spin frequencies ($`\nu >1`$ Hz). Be stars show both emission in the Balmer lines, and an excess in the IR, due to a disk of material shed from the equator of the rapidly rotating star. The accretion of this material, which is thought to form a quasi-keplerian disk, fuels the X-ray outbursts.

Two types of outbursting behavior are observed in Be/pulsar binaries: type I (“normal” outbursts), with a series of lower luminosity ($`L<10^{37}\mathrm{erg}\mathrm{s}^1`$ outbursts occurring once per orbit, generally near periastron; and type II (“giant” outbursts), single high luminosity ($`L10^{38}\mathrm{erg}\mathrm{s}^1`$ sometimes lasting several orbits (Stella, White & Rosner  (1986)). These outbursts are generally accompanied by rapid spin-up of the neutron star, indicating disk accretion (Finger, Wilson & Harmon  (1996); Bildsten et al.  (1997)).

The BATSE observations presented here are consistent with a giant outburst of a Be/pulsar transient in the SMC. From the observed spin-up we can show that a disk is present, and the the source is outside our galaxy. In disk accretion, the angular momentum accreted per unit accreted mass is $`l=(GMr_m)^{{\scriptscriptstyle \frac{1}{2}}}`$, were $`M`$ is the neutron star mass, and $`r_m`$ the radius of the magnetosphere.

The specific angular momentum for wind accretion is generally much less than this, since if it reaches this value, a disk forms. The magnetospheric radius is limited by the corotation radius $`r_c=(GM)^{{\scriptscriptstyle \frac{1}{3}}}(2\pi \nu )^{{\scriptscriptstyle \frac{2}{3}}}=1.0\times 10^9\mathrm{cm}`$, where we have used $`M=1.4M_{}`$. The specific angular momentum $`l`$ is proportional to the ratio $`\dot{\nu }/L`$ where $`L`$ is the pulsar’s luminosity. Since the rate at which angular momentum is accreted is $`l\dot{m}=2\pi I\dot{\nu }`$ where I is neutron star moment of inertia, and $`L=GM\dot{m}/R`$, where $`R`$ is the neutron star radius, we have the upper limit

$$\frac{\dot{\nu }}{L}<(2\pi )^{{\scriptscriptstyle \frac{4}{3}}}(GM)^{{\scriptscriptstyle \frac{1}{3}}}I^1R\nu ^{{\scriptscriptstyle \frac{1}{3}}}=3.8\times 10^{49}\mathrm{Hz}\mathrm{erg}^1$$
(2)
where we have used a moment of inertia $`I=10^{45}\mathrm{gm}\mathrm{cm}^2`$, and a radius $`R=10^6\mathrm{cm}`$. For the Nov 11-17 observation, for which $`\dot{\nu }=5.6\times 10^{11}\mathrm{Hz}\mathrm{s}^1`$, this implies $`L>1.5\times 10^{38}\mathrm{erg}\mathrm{s}^1`$, and for peak frequency rate observed we have $`L>2.0\times 10^{38}\mathrm{erg}\mathrm{s}^1`$. Using the 0.1-2.0 keV flux $`F=5.9\times 10^{11}\mathrm{erg}\mathrm{cm}^2\mathrm{s}^1`$ observed in the 1996 Nov 11-17 HRI observation, we find the lower limit on the source distance of

$$d>46\left(\frac{ϵ}{0.1}\right)^{{\scriptscriptstyle \frac{1}{2}}}\mathrm{kpc}$$
(3)
where $`ϵ`$ is the fraction of the luminosity in the 0.1-2.0 keV energy band. This eliminates the possibility that the RX J0052.1-7319 is in the foreground of the SMC.

The limit in equation 2 is reached only with disk accretion. Using the measured frequency rate, fluxes, and the distance of 60 kpc to the SMC we find a ratio of $`\dot{\nu }/L`$ comparable to this limit, suggesting that an accretion disk is present during the outburst.

The pulse periods of Be binary pulsars range from 69 ms to 1400 s, with determined orbital periods ranging from 17 to 250 days. Corbet (1986) showed that the orbital period is correlated with the pulse period. From the observed distribution of spin and orbital periods we would expect an orbital period in the range of 25 – 100 days. In giant outbursts of Be/X-ray pulsars we expect a strong correlation of flux and frequency rate (see e.g. Finger, Wilson & Harmon  (1996)). We note however that the history of the frequency rate in fig 2.3 is dissimilar in profile to that of the pulsed flux. This could be due to the doppler signature of a binary orbit. An alternate explanation would be changes in the spectra or pulse fraction. However, giant outbursts typically have simple flux and intrinsic spin-up rate profiles, with a steady rise to peak, and a somewhat slower fall (see e.g. Parmar et al.  (1989); Whitlock  (1989); Finger, Wilson & Harmon  (1996), but also, Negueruela et al.  (1997)).

The plausible identification of RX J0052.1-7310 as Be binary system in the SMC accentuates further the rather dramatic difference between the SMC and our Galaxy with regard to the population of high mass X-ray binaries. This fact has already been noted by several authors (Schmidtke et ak,  (1999); Jokogawa et al.  (2000)). A recent compilation of the known X-ray pulsars444http://gammaray.msfc.nasa.gov/batse/pulsar/asm\_pulsars.html gives within the SMC one supergiant system, three known Be systems, and 11 transients with uncertain companion class (likely to be Be systems), making 15 high mass X-ray pulsar binaries. For the Galaxy the corresponding number is 40. Therefore, using a mass ratio of the SMC to the Galaxy of 1/100, this suggests that high mass X-ray pulsar systems in the SMC are overabundant by roughly a factor of 30 relative to the Galaxy. This analysis ignores the important effects of obscuration within the Galaxy, the relative frequency of observations, and the low luminosity sensitivites obtained for the Magellantic clouds; nevertheless the apparent disparity is remarkable.

Since high mass X-ray binaries have lifetimes which are a very small fraction ($`10^3`$) of the age of the Galaxy, the dramatic difference between the SMC and the Galaxy points to a rather recent outburst of star-formation in the SMC within the last $`10^7`$ years. Further support of such an epoch of star formation comes from the radio observations of H1 by Stavely-Smith et al. (1997) and Putman et al. (1998), which show a strong bridge of material between the Magellanic Clouds and between them and our own galaxy. Furthermore, Stavely-Smith et al. have demonstrated the existence of a large number of supershells (created by multiple supernovae) of a similar age ($``$5 Myr), strongly suggesting enhanced starbirth has taken place as a result of tidal interactions between these component systems. Consequently it seems very likely that the previous closest approach of the SMC to the LMC $`10^8`$ years ago may have triggered the birth of many new massive stars which have given rise to the current population of HMXBs. In fact, other authors (e.g. Popov et al.  (1998)) claim that the presence of large numbers of HMXBs may be the best indication of starburst activity in a system.

M. H. F. acknowledges support through NASA grant NAG5-4238.


format: text

text:
# Scarring and the statistics of tunnelling

## I Introduction

Tunnelling rates from systems with complex internal behaviour are naturally described statistically. Such analyses have been used in particular to understand nuclear resonances , chemical reaction rates  and quantum dots . While the manner in which internal states couple to the continuum varies in these analyses, the internal quantum mechanics is usually modelled using one of the standard ensembles of random matrix theory (RMT). We show here that in systems where the complexity of the internal dynamics derives from low-dimensional chaos, sufficiently strong deviation from the standard RMT statistics is possible that it dominates the statistics of tunnelling. These deviations were pointed out in  and are explained in detail here using the theory of scarring developed by Heller, Kaplan and coworkers .

We consider tunnelling to and from regions of phase space associated with chaotic behaviour. In particular, our theory works for the calculation of level splittings in double well potentials and of resonance widths of quasibound states in metastable wells. With the development of experimental techniques capable of measuring state-specific reaction rates, such tunnelling from highly-excited chaotic states has become directly relevant to the analysis of chemical reactions, for example. Specific examples of reactions which have been treated using RMT-based statistical methods include the dissociation rates of $`\mathrm{D}_2\mathrm{CO}`$  and $`\mathrm{N}_2\mathrm{O}`$ . In , the distributions arising from such an analysis were shown to be determined in a simple way by the stability and action of a complex tunnelling orbit which crosses the potential barrier with minimum imaginary action. The resulting statistical distributions for the tunnelling rate agree well with numerically computed ensembles except when the real extension of the tunnelling orbit into the potential well is periodic; in that case, strong deviations from the RMT prediction are observed and it was proposed in  that these are due to the effect of scarring on wavefunction statistics as outlined in . Additional evidence in support of this has subsequently been provided in .

In this paper we derive a distribution describing scar-influenced tunnelling-rate statistics in two-dimensional potentials. It is determined completely by three dynamical parameters which are calculated from the monodromy matrices of the tunnelling orbit and of its real periodic extension. These three parameters can be understood as the stability of the scarred orbit, the stability of the complex tunnelling orbit and an angle relating stable and unstable manifolds. The root of this calculation is a conjecture governing the statistics of the components of chaotic wavefunctions in the eigenbasis of a tunnelling operator which was defined in . The low-lying states of this tunnelling operator are approximated by the eigenstates of a harmonic oscillator. Following the linear scarring theory of Kaplan  we show that the components of chaotic states in this basis have variances which deviate from those of RMT and which must furthermore be correlated. We conjecture that the eigenfunction distribution is a nonisotropic Gaussian. This is the simplest distribution consistent with the observed correlations and may also be obtained from maximum-entropy arguments. This conjecture is tested against quantum-map models and found to describe statistics of their eigenstates very accurately. The conjecture is then used to derive modified distributions for tunnelling rates in the presence of scarring and these are found to describe well the statistics of energy-level splittings in chaotic double well potentials.

We see the benefits of this work as being twofold. First, there is a solution to the primary problem of incorporating the effects of scarring into the statistics of tunnelling. This enables the detection and interpretation of system-specific dynamical detail in measurements of multidimensional tunnelling and furthermore provides a very direct manifestation of the phenomenon of scarring in quantities with real physical relevance. Second, we believe that from the point of view of scarring alone, even without reference to tunnelling, the correlated joint probability distribution proposed here for the statistics of eigenstates holds considerable promise as a theoretical tool. In particular, it suggests methods of quantifying the effects of scarring in a way that does not depend on a choice of test states. This second aspect is not fully exploited in the present paper, but we believe that the essential elements needed for a full development are put in place. We remark that the effect of correlations on the statistics of tunnelling rates was also considered in a somewhat different context in .

An outline of the paper is as follows. We begin in section II with a review of the existing theory of the statistics of tunnelling rates. This includes a definition of the tunnelling operator introduced in  and a summary of the derivation from it of nonscarred tunnelling-rate distributions as outlined previously in . In section III we use the linear theory of scarring to deduce that in the presence of scarring the components of chaotic eigenstates in an eigenbasis of the tunnelling operator must be nonRMT and conjecture a modified Gaussian distribution for them. This conjecture is tested on quantisations of perturbed cat maps and found to work well. In section IV we use the conjecture to deduce modified tunnelling-rate distributions in the presence of scarring and these are compared successfully with distributions of splitting in chaotic double-well potentials.

## II Tunnelling statistics without scarring

We begin with a brief description of existing theory of tunnelling-rate statistics in the absence of scarring, on which our present work is based. We first describe a tunnelling operator defined in  which relates tunnelling rates to the properties of wavefunctions in the allowed region. We then outline how this tunnelling operator is used to derive distributions for tunnelling rates in the absence of scarring.

### A The tunnelling operator

A full description of the semiclassical matrix element and complete definitions of the various objects needed to calculate it are described in . For the purposes of statistical analysis, an idealised model formulated in terms of quantum maps suffices, which we summarise here.
Our model starts with a finite-dimensional Hilbert space $``$. We take this space to be the quantum analog of a Poincaré section $`\mathrm{\Sigma }`$ in a potential well , with the dimension $`N`$ of $``$ being approximately proportional to the area of $`\mathrm{\Sigma }`$. A quantum mapping acts on the space $``$ as a unitary operator $`\widehat{U}`$, whose eigensolutions we denote by

$`\widehat{U}|n=e^{i\theta _n}|n`$

and whose classical analog is a symplectic map $`F`$ corresponding to chaotic motion within the well. In , the eigenstates $`|n`$ are represented as certain cross-sections of the eigenfunctions of the usual Hamiltonian in the full Hilbert space.

Tunnelling rates are calculated using a tunnelling operator $`\widehat{𝒯}`$ which also acts on $``$ and from which the scaled tunnelling rate associated with a chaotic eigenstate $`|n`$ is calculated from

$$y_n=\frac{n|\widehat{𝒯}|n}{\mathrm{Tr}\widehat{𝒯}}.$$
(1)
In this formula, we normalise the chaotic eigenstates according to $`n|n=N`$ so that, in particular, the overlaps with a second basis $`|\stackrel{~}{k}`$ of $``$ satisfy $`|\stackrel{~}{k}|n|^2=1`$. As a result, $`y_n=1`$. In resonance problems, this scaled tunnelling rate is $`y_n=\mathrm{\Gamma }_n/\mathrm{\Gamma }_n`$ where $`\mathrm{\Gamma }_n`$ is the resonance width associated with a particular metastable state labeled by $`n`$ and in double-well problems $`y_n=\mathrm{\Delta }E_n/\mathrm{\Delta }E_n`$ where $`\mathrm{\Delta }E_n`$ is the splitting of a doublet labeled by $`n`$.

The tunnelling operator $`\widehat{𝒯}`$ is interpreted as an evolution operator whose classical correspondent is a symplectic map $``$ which is complex. This map $``$ has a real fixed point $`\zeta _0`$ which corresponds to the most probable tunnelling route. That is, using $`\zeta _0`$ as an initial condition, a complex trajectory of imaginary period can be traced out which crosses a potential barrier with minimum imaginary action. Mappings of points in a neighbourhood of $`\zeta _0`$ in the Poincaré section $`\mathrm{\Sigma }`$ correspond to complex trajectories in full phase space which cross the potential barrier near this central tunnelling route. The matrix element in (1) then samples the state $`|n`$ in a small region of $`\mathrm{\Sigma }`$ surrounding $`\zeta _0`$ and with area of $`O(\mathrm{})`$ (as determined by say, a Wigner function on $`\mathrm{\Sigma }`$ ). This region is determined by the dynamics of $``$ near $`\zeta _0`$ and outside of it the symbol of $`\widehat{𝒯}`$ decays exponentially. For this reason it is consistent within semiclassical approximation to linearise dynamics around $`\zeta _0`$ and let $`\widehat{𝒯}`$ be the quantisation of the resulting complex symplectic matrix $`W`$. Note that in practical terms $`W`$ is simply the monodromy matrix of the tunnelling orbit.

We now describe the tunnelling operator for resonance problems of symmetric double wells where the symmetry is a reflection $`(x,y)(x,y)`$. If we choose a representation in which $`\zeta _0`$ is at the centre of coordinates $`\zeta =(q,p)`$ on $`\mathrm{\Sigma }`$, we may write

$`\widehat{𝒯}=e^{\alpha _0\widehat{h}/\mathrm{}},`$

for these cases, where $`\widehat{h}`$ is quadratic in the corresponding operators $`\widehat{\zeta }=(\widehat{q},\widehat{p})`$ and $`\alpha _0>0`$ . The corresponding symbol $`h(q,p)`$ is a quadratic form in $`(q,p)`$, expressed as

$`h(q,p)={\displaystyle \frac{1}{2}}\zeta ^TK\zeta ,`$

where $`K`$ is a positive-definite symmetric matrix. The matrix $`K`$ is determined from $`W`$ by writing

$`W=e^{i\alpha _0JK}.`$

where $`J`$ is the unit symplectic matrix. We normalise $`K`$ so that $`\mathrm{det}K=1`$ and then the parameter $`\alpha _0`$ is fixed by a calculation of $`W`$. That $`W`$ can be written in this way for real $`\alpha _0`$ and $`K`$ can be shown using the conjugate-time-reversal symmetry

$`W^{}=W^1,`$

which is a special case of a similar symmetry of $``$ . This symmetry is a reflection of the fact that the tunnelling orbit is periodic with an imaginary period and so complex conjugation is equivalent to time-reversal. Note that this is the case even if the problem does not have a proper time-reversal symmetry (such as in the presence of magnetic fields, for example).

We limit the detailed discussion in this paper to two-dimensional potentials, so that the Poincaré section $`\mathrm{\Sigma }`$ has one degree of freedom. In that case we may denote the eigensolutions of $`\widehat{h}`$ by

$`\widehat{h}|\stackrel{~}{k}=\left(k+{\displaystyle \frac{1}{2}}\right)\mathrm{}|\stackrel{~}{k},`$

using a tilde to distinguish them from the chaotic eigenstates $`|n`$. Note then that the eigensolutions of $`\widehat{𝒯}`$ are then

$`\widehat{𝒯}|\stackrel{~}{k}=e^{(k+1/2)\alpha _0}|\stackrel{~}{k}={\displaystyle \frac{1}{\mathrm{\Lambda }^{k+1/2}}}|\stackrel{~}{k},`$

where $`\mathrm{\Lambda }=e^{\alpha _0}`$. Notice that $`\mathrm{\Lambda }`$ is the larger of the two eigenvalues $`e^{\pm \alpha _0}`$ of $`W`$ (which are real even though $`W`$ is complex). Finally, we note that in double well problems for which the symmetry is an inversion $`(x,y)(x,y)`$, the eigenvalues of $`W`$ as defined in  are negative and the eigensolutions of the tunnelling operator are of a similar form except that the eigenvalues are $`\mathrm{\Lambda }^k\left|\mathrm{\Lambda }\right|^{1/2}`$.

### B Statistics using standard RMT

A standard statistical model of eigenstates of chaotic systems is that their components in a generic basis are Gaussian-distributed. This assumption has been combined with various models of coupling to the continuum  to produce statistical distributions for the scaled tunnelling rates in the absence of scarring. The discussion in  forms the basis for our treatment of scarred distributions and we therefore summarise it here. We assume for ease of presentation that the system is time-reversal symmetric so the assumptions of the GOE are adopted. Results for GUE systems are derived similarly and are summarised at the end.

We expand the eigenstate $`|n`$ in an eigenbasis of $`\widehat{𝒯}`$,

$`|n={\displaystyle \underset{k}{}}x_k|\stackrel{~}{k},`$

suppressing the dependence of the coefficients $`x_k`$ on the chaotic-state index $`n`$. The normalisation of $`|n`$ is such that

$`\left|x_k\right|^2=1`$

and in systems with time-reversal symmetry, the $`x_k`$’s are real. The standard RMT model is that for large $`N`$, the $`x_k`$’s are uncorrelated and distributed with a joint probability distribution

$$P(𝐱)=\underset{k=0}{\overset{N1}{}}\left[\frac{1}{\sqrt{2\pi }}e^{x_k^2/2}\right]=\frac{1}{(2\pi )^{N/2}}e^{𝐱^T𝐱/2},$$
(2)
where $`𝐱=(x_0,x_1,x_2,\mathrm{})`$. Given that (1) expresses each $`y`$ as a quadratic form

$`y=𝐱^TT𝐱`$

in $`𝐱`$, where $`T`$ is an $`N\times N`$ matrix representing $`\widehat{𝒯}/\mathrm{Tr}\widehat{𝒯}`$ in a $`|\stackrel{~}{k}`$-basis, we can write the distribution for $`y`$ in the form

$$p(y)=\delta \left(y𝐱^TT𝐱\right)P(𝐱)d𝐱.$$
(3)
Fourier-transforming with respect to $`y`$ gives the characteristic function,

$`\stackrel{~}{p}(q)`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}e^{iqy}p(y)dy`$ (4)
$`=`$ $`{\displaystyle \frac{1}{(2\pi )^{N/2}}}{\displaystyle \mathrm{exp}\left[\frac{1}{2}𝐱^T(I2iqT)𝐱\right]d𝐱}`$ (5)
$`=`$ $`{\displaystyle \frac{1}{\sqrt{\mathrm{det}\left(I2iqT\right)}}}.`$ (6)
One finds similarly that

$$\stackrel{~}{p}(q)=\frac{1}{\mathrm{det}\left(IiqT\right)}$$
(7)
in the absence of time-reversal symmetry.

These distributions have been shown in  to describe successfully the tunnelling rate statistics of chaotic double wells when the tunnelling route has a nonperiodic real extension. When the tunnelling route connects to a periodic orbit, however, strong deviations are found  which we explain in the next sections using the idea of scarring.

## III Wavefunction statistics and scarring

It has been pointed out by Kaplan and coworkers that the Gaussian distribution in (2) does not describe wavefunction statistics if a basis state $`|\stackrel{~}{k}`$ is localised in phase space near a periodic orbit of the real classical dynamics (corresponding to the Poincaré map $`F`$). In tunnelling-rate statistics, this can lead to strong deviations from the distribution described the previous section if $`\zeta _0`$, in addition to being a fixed point of the complex map $``$, is a fixed point or short periodic orbit of $`F`$. In this section we review elements of that argument and use them to calculate correlations among the components of $`𝐱`$. From this analysis there emerges a natural conjecture for generalising (2). The distributions we calculate on the basis of this conjecture are necessary to understand tunnelling but also, independently of the tunnelling problem, hold promise as a basis for understanding the statistics of scarring in a rather general way.

### A A Gaussian hypothesis

For a given ensemble of states with overlaps $`𝐱=(x_0,x_1,x_2,\mathrm{})`$, let us define a matrix of correlations $`C`$ whose elements are

$$C_{lk}=x_lx_k^{}.$$
(8)
The average is over the eigenstates $`|n`$ and will be made more precise below. We will in particular consider ensembles formed by the eigenstates $`|n`$ of $`\widehat{U}`$ whose eigenangles lie in a subset of the unit circle. By varying the quantum dimension $`N`$, for example, we can consider ensembles corresponding to fixed classical maps $`F`$ and $``$. If the joint probability distribution is of the form (2) or its GUE equivalent, we expect

$$C=I,$$
(9)
which says simply that $`|x_k|^2=1`$ and that $`x_k`$ and $`x_l`$ are uncorrelated unless $`k=l`$. In order to avoid confusion with a correlation function defined below, we will henceforth refer to $`C`$ as the covariance matrix. We now show that deviations of $`C`$ from the identity can be calculated on the basis of linearised dynamics around the point $`\zeta _0`$. The construction begins with a calculation of quantum recurrence for basis states localised near $`\zeta _0`$ and is finished by a Fourier transformation which relates these to the averages $`x_lx_k^{}`$.

Let $`\{|\stackrel{~}{k}|k=0,1,\mathrm{}\}`$ be an eigenbasis of $`\widehat{𝒯}`$ as defined in the previous section and suppose that $`\zeta _0`$ is a fixed point of the classical map $`F`$ corresponding to $`\widehat{U}`$. We now argue that for some ensembles (9) cannot hold. Following  we define the correlation function

$`A_{lk}(t)`$ $`=`$ $`\stackrel{~}{l}|\widehat{U}^t|\stackrel{~}{k}`$ (10)
$`=`$ $`{\displaystyle \frac{1}{N}}{\displaystyle \underset{n=1}{\overset{N}{}}}x_lx_k^{}e^{i\theta _nt}`$ (11)
and denote its Fourier transform by

$`S_{lk}(\theta )`$ $`=`$ $`{\displaystyle \underset{t=\mathrm{}}{\overset{\mathrm{}}{}}}e^{i\theta t}A_{lk}(t)`$ (12)
$`=`$ $`{\displaystyle \frac{2\pi }{N}}{\displaystyle \underset{n=1}{\overset{N}{}}}x_lx_k^{}\delta (\theta \theta _n).`$ (13)
The correlation function is approximated semiclassically by a sum of the form

$$A_{lk}(t)\underset{p}{}𝒜_pe^{2\pi NiS_p}$$
(14)
over classical orbits $`p`$, typically complex, which begin in and return to a small region around $`\zeta _0`$ after $`t`$ iterations of $`F`$, assuming $`k`$ and $`l`$ are $`O(1)`$. Each orbit contributes a term with an amplitude $`𝒜_p`$ which is independent of $`N`$ and an exponential $`\mathrm{exp}[2\pi NiS_p]`$ with a rapidly-varying complex phase which decays exponentially with the distance of the contributing orbit from $`\zeta _0`$. The primary such orbit corresponds to the fixed point $`\zeta _0`$ itself. We are free to choose the phase of $`\widehat{U}`$ so that the action of $`\zeta _0`$ vanishes (this amounts to a choice of the origin of the $`\theta _n`$ values). With that choice it is manifest that its contribution to $`A_{lk}(t)`$ is $`N`$-independent.

Consider averaging $`A_{lk}(t)`$ over a sequence of quantum systems with varying $`N`$ but identical classical limits. We find then that

$`A_{lk}(t){\displaystyle \underset{p}{}}𝒜_pe^{2\pi NiS_p}`$

and, with the exception of the fixed point $`\zeta _0`$ for which $`S_p=0`$, the contributions average to zero. The contribution from $`\zeta _0`$ can be obtained from linearised dynamics and we find that

$`A_{lk}(t)A_{lk}^{\mathrm{lin}}(t)\stackrel{~}{l}|\widehat{U}_{\mathrm{lin}}^t|\stackrel{~}{k}`$

where $`\widehat{U}_{\mathrm{lin}}`$ is the quantisation of a linearisation of the real map $`F`$ at $`\zeta _0`$. The linearised correlation function $`A_{lk}^{\mathrm{lin}}(t)`$ is calculated in the next section. It is $`N`$-independent and decays exponentially to zero with time when $`k`$ and $`l`$ are $`O(1)`$.

Performing the same averaging procedure on $`S_{lk}(\theta )`$ gives

$`S_{lk}(\theta )S_{lk}^{\mathrm{lin}}(\theta ){\displaystyle \underset{t=\mathrm{}}{\overset{\mathrm{}}{}}}e^{i\theta t}A_{lk}^{\mathrm{lin}}(t).`$

This is a real, $`2\pi `$-periodic function of $`\theta `$ which is also $`N`$-independent. Suppose we form an ensemble of chaotic states by varying $`N`$ and selecting states for which $`\theta _n`$ lies within a narrow window centred about $`\theta `$. The eigenfunction expansion in (12) indicates that

$$C_{lk}(\theta )=x_lx_k^{}=S_{lk}(\theta )S_{lk}^{\mathrm{lin}}(\theta ).$$
(15)
We find therefore that, if we produce an ensemble of overlaps from states within a narrow window of phase angle, the covariance matrix violates condition (9) and the joint probability distribution for $`𝐱`$ necessarily deviates from the standard RMT form.

The simplest distribution consistent with (15) is the Gaussian

$$P(𝐱;\theta )=\frac{1}{\left[(2\pi /\beta )^N\mathrm{det}C(\theta )\right]^{\beta /2}}e^{\beta 𝐱^TC^1(\theta )𝐱/2}$$
(16)
where $`\beta =1`$ in the GOE case and $`\beta =2`$ in the case of GUE. We conjecture that these distributions govern the wavefunction statistics in the presence of scarring when states are taken from a small window of eigenphases. A full justification of this would require a detailed analysis of long orbits contributing to (14), along the lines of the nonlinear theory of scarring in . We will not offer such an analysis here and assume (16) on the grounds of simplicity alone. This may be formalised by noting that (16) minimises the information content

$`\left[P(𝐱)\right]={\displaystyle P(𝐱)\mathrm{ln}P(𝐱)d𝐱},`$

subject to the constraint that the correlations are as written in (8) and that $`N`$ is large. Minimisation of information is not a proof, however, and ultimately we rely on detailed numerical testing to justify the choice of distribution. We now show in detail how $`C(\theta )`$ may be calculated and verify that the Gaussian hypothesis describes quite well the scarred wavefunction statistics of specific quantum-map models.

### B Calculation of the covariance matrix

Construction of the scarred distribution in (16) begins with the calculation of the linearised correlation function $`A_{lk}^{\mathrm{lin}}(t)`$, from which the covariance matrix $`C(\theta )`$ is obtained by Fourier transformation. We show in this subsection how $`A_{lk}^{\mathrm{lin}}(t)`$ may be calculated from the monodromy matrices linearising real and complex dynamics around the fixed point $`\zeta _0`$ (corresponding to the maps $`F`$ and $``$ respectively). In the main text we present enough detail that the method of calculation of $`A_{lk}^{\mathrm{lin}}(t)`$ should be clear and leave the derivation and other detailed discussion to appendices. Once $`A_{lk}^{\mathrm{lin}}(t)`$ has been constructed as outlined below, the covariance matrix is then easily calculated in practice using a FFT, for example.

It is shown in appendix A that $`A_{lk}^{\mathrm{lin}}(t)`$ vanishes if $`l`$ and $`k`$ are not both even or both odd and may otherwise be calculated from the polar form

$$A_{k+2n,k}^{\mathrm{lin}}(t)=G_{kn}(\psi (t))e^{i(k+1/2)(\varphi (t)\mu t\pi )+in\varphi (t)+in\pi /2},$$
(17)
where the Maslov index $`\mu `$ and the angles $`\psi (t)`$ and $`\varphi (t)`$ are defined by dynamics around $`\zeta _0`$ as described below. The amplitude is

$`G_{kn}(\psi )=\sqrt{{\displaystyle \frac{k!}{(k+2n)!}}}{\displaystyle \frac{(2n)!}{2^nn!}}\mathrm{sin}^n\psi \sqrt{\mathrm{cos}\psi }C_k^{n+1/2}(\mathrm{cos}\psi ),`$

where $`C_k^\alpha (x)`$ denotes a Gegenbauer polynomial. In the special case $`n=0`$ we get the autocorrelation function

$`A_{kk}^{\mathrm{lin}}(t)=\sqrt{\mathrm{cos}\psi (t)}P_k(\mathrm{cos}\psi (t))e^{i(k+1/2)(\varphi (t)\mu t\pi )},`$

where $`P_k(x)`$ is a Legendre polynomial. We assume $`n0`$ and $`t>0`$ in (17) and use $`A_{lk}^{\mathrm{lin}}{}_{}{}^{}(t)=A_{lk}^{\mathrm{lin}}(t)=A_{kl}^{\mathrm{lin}}(t)`$ to calculate the correlation function when $`t<0`$ or $`l<k`$. Note that this assumes a convention for the phases of the eigenstates $`|\stackrel{~}{k}`$ that is outlined at the end of appendix B. In particular, this convention leads to a covariance matrix $`C(\theta )`$ which is real-symmetric even in the case of GUE statistics; that is, while individual values of $`x_lx_k^{}`$ are complex, the phase convention is such that the averages are real.

To define the angles $`\varphi (t)`$ and $`\psi (t)`$, let $`M_0`$ be the real symplectic matrix linearising the unstable real dynamics around $`\zeta _0`$ and let

$$W(z)=e^{i(\mathrm{ln}z)JK}$$
(18)
be a generalisation of the complex symplectic matrix described in the the last section, reducing to that case when $`z=e^{\alpha _0}`$. Euler expansion of the exponential in (18) leads to the identity

$`\mathrm{Tr}W(z)M_0^t=(1)^{\mu t}\left(m(t)z+{\displaystyle \frac{m^{}(t)}{z}}\right),`$

where

$`m(t)=\mathrm{cosh}\rho t+iQ\mathrm{sinh}\rho t`$

and where $`(1)^\mu e^{\pm \rho }`$ are the eigenvalues of $`M_0`$. The real parameter $`Q`$ follows the notation of  and is given by $`Q=\mathrm{cot}\phi `$, where $`\phi `$ is the angle between the stable and unstable manifolds of $`M_0`$. Note that the parameter $`Q`$ depends implicitly on $`W`$ because the angle $`\phi `$ is calculated using the metric defined on phase space by $`K`$. Finally, we complete the explanation of (17) by defining the angles $`\psi (t)`$ and $`\varphi (t)`$ implicitly using the polar decomposition

$`m(t)=\mathrm{sec}\psi (t)e^{i\varphi (t)}`$

of $`m(t)`$. We may choose these angles to lie in the ranges $`0<\psi (t)<\pi /2`$ and $`\psi (t)<\varphi (t)<\psi (t)`$ respectively (note that $`\mathrm{sec}\psi \mathrm{cos}\varphi =\mathrm{cosh}\rho t>1`$) and, with this restriction on $`\varphi (t)`$, the half-angle in (17) is well-defined.

Note that, as $`t\mathrm{}`$, we find

$`\mathrm{cos}\psi (t){\displaystyle \frac{2}{\sqrt{1+Q^2}}}e^{\rho t}`$

and the correlation function in (17) therefore decays at the exponential rate $`e^{\rho t/2}`$ when $`k`$ is even and at the rate $`e^{3\rho t/2}`$ when $`k`$ is odd (since $`C_k^{n+1/2}`$ is then an odd polynomial). Less obvious from (17) is the fact that the correlation function also decays for fixed $`t`$ with increasing $`k`$ and $`l=k+2n`$. This is a reflection of the fact that more excited states of the oscillator $`\widehat{h}`$ are less localised at the fixed point $`\zeta _0`$ and therefore less affected by scarring.

When the linear correlation function is Fourier-transformed, we therefore get a covariance matrix $`C_{lk}(\theta )`$ which can be decomposed into two blocks corresponding to odd and even $`l`$ and $`k`$ and whose elements approach those of the identity matrix ($`C_{lk}(\theta )\delta _{lk}`$) as $`l`$ and $`k`$ increase. We expect scarring therefore to affect the statistics of $`x_k`$ significantly for relatively small values of $`k`$ and to have little effect when $`k`$ is large. This exponential convergence towards identity means that the matrix $`CI`$ is trace-class and the determinant of $`C`$ exists in the limit $`N\mathrm{}`$. In the next subsection we calculate the resulting covariance matrix for some model quantum maps and show that the Gaussian hypothesis provides a good model for the wavefunction statistics of scarred systems.

### C Testing the Gaussian hypothesis

We test the Gaussian hypothesis using a model system in which $`\widehat{U}`$ quantises the perturbed cat map,

$`F_ϵ=F_0M_ϵ,`$
where

$`F_0:\left(\begin{array}{c}q\\ p\end{array}\right)\left(\begin{array}{cc}1& 1\\ \\ 1& 2\end{array}\right)\left(\begin{array}{c}q\\ p\end{array}\right)\mathrm{mod}\mathrm{\hspace{0.33em}1}`$

is the standard cat map and

$`M_ϵ:\left(\begin{array}{c}q\\ p\end{array}\right)\left(\begin{array}{c}qϵ\mathrm{sin}2\pi p\\ p\end{array}\right)\mathrm{mod}\mathrm{\hspace{0.33em}1}`$

is a kicked Harper map. For sufficiently small values of the parameter $`ϵ`$ this map shares the hyperbolic structure of the unperturbed map but not its nongeneric degeneracies. Details of the quantisation of these maps for Hilbert spaces of arbitrary dimension $`N`$ may be found in .

The map $`F_ϵ`$ has a symmetry of inversion about the origin, which is always therefore a fixed point. It also has a less obvious time-reversal symmetry. One finds that

$`PF_ϵP=F_ϵ^1,`$

where

$`P:\left(\begin{array}{c}q\\ p\end{array}\right)\left(\begin{array}{cc}1& 0\\ \\ 1& 1\end{array}\right)\left(\begin{array}{c}q\\ p\end{array}\right)\mathrm{mod}\mathrm{\hspace{0.33em}1}.`$

Note that $`P`$ is antisymplectic and that $`P^2=I`$. The map therefore allows us to test statistics in the GOE scenario. (A systematic discussion of such time-reversing symmetries can be found in ).

To complete our model of a tunnelling system, we construct an analog of the tunnelling operator. Consider the Harper-like Hamiltonian

$`h(q,p)=2\mathrm{cos}\pi (2pq)\mathrm{cos}2\pi q.`$

This is invariant under the time-reversal symmetry $`P`$ and has a minimum at the origin in phase space. The function does not have the full periodicity of the torus in $`q`$ but, by placing the resulting discontinuity at $`q=1/2`$, may be quantised in such a way that the ground and first few excited states of its quantisation $`\widehat{h}`$ are localised near the origin. These are then approximated semiclassically by those of the harmonic oscillator $`(\widehat{p}\widehat{q}/2)^2+\widehat{q}^2`$. The statistics of the components $`x_k`$ of chaotic eigenstates of $`\widehat{U}`$ in an eigenbasis of $`\widehat{h}`$ therefore provide a test of the Gaussian hypothesis outlined in the previous subsection and since $`\widehat{U}`$ and $`\widehat{h}`$ share a time-reversal symmetry, the GOE hypothesis is appropriate.

We first consider the prediction for the covariance matrix $`C(\theta )`$. The element $`C_{00}(\theta )`$, which corresponds to the overlap statistics of a Gaussian wavepacket (defined as the ground state of $`\widehat{h}`$), has been considered extensively by Kaplan and coworkers. For comparison with later calculation it is plotted as a function of $`\theta `$ in Fig. 1$`(a)`$ for the map $`F_ϵ`$ with $`ϵ=0.1`$. Following , we refer to this curve as the spectral envelope. We expect larger than average overlaps where $`C_{00}(\theta )>1`$ and we will refer to this as the scarred region. The rest of the interval is the antiscarred region, where we expect smaller than average overlaps. The points in the figure give a representative sample of individual overlaps. There is also a curve representing a binned average of these overlaps which is obtained by accumulating quantum data for even values of $`N`$ between $`100`$ and $`700`$. It is also important for our analysis to understand the statistics of $`x_lx_k`$ when $`l`$ and $`k`$ lable excited states of $`\widehat{h}`$. We therefore make an analogous comparison with quantum data for $`C_{02}(\theta )`$ and $`C_{22}(\theta )`$ in Figs. 1$`(b)`$ and 1$`(c)`$ respectively. Note the different vertical scales in these cases. In all cases the average obtained by binning quantum data agrees well with the prediction $`C_{lk}(\theta )`$ obtained by Fourier-transforming the correlation function. This is to be expected since we have a relatively direct argument showing that this should be so. Less obvious is that the distribution of individual overlaps about these averages should be governed by the Gaussian hypothesis, since we do not have a proof for that. That this assumption can be verified numerically (below) is therefore nontrivial.

To compare the distributions about the mean $`C_{lk}(\theta )`$, lets us denote

$`\xi _{lk}=\mathrm{Re}x_lx_k^{}=𝐱^{}B_{lk}𝐱,`$

where $`B_{lk}`$ is a symmetric matrix whose $`ij^{\mathrm{th}}`$ entry is $`(\delta _{il}\delta _{jk}+\delta _{ik}\delta _{jl})/2`$. (For the moment we allow either GUE or GOE statistics even though all the numerical calculations we present are for the GOE case). Similarly to the calculation of $`p(y)`$ in section II, we find that the local distribution $`p(\xi _{lk};\theta )`$ is predicted by the Gaussian hypothesis to have the Fourier transform

$`\stackrel{~}{p}_{lk}(q;\theta )`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}e^{iq\xi _{lk}}p(\xi _{lk};\theta )d\xi _{lk}`$ (19)
$`=`$ $`{\displaystyle \frac{1}{\left[(2\pi /\beta )^N\mathrm{det}C(\theta )\right]^{\beta /2}}}{\displaystyle \mathrm{exp}\left[\frac{1}{2}𝐱^{}(C^1(\theta )\frac{2iq}{\beta }B_{lk})𝐱\right]d𝐱}`$ (20)
$`=`$ $`\left[\mathrm{det}\left(I{\displaystyle \frac{2iq}{\beta }}B_{lk}C(\theta )\right)\right]^{\beta /2}.`$ (21)
The determinant simplifies on expansion to

$`\mathrm{det}\left(I{\displaystyle \frac{2iq}{\beta }}B_{lk}C(\theta )\right)=1{\displaystyle \frac{2iq}{\beta }}C_{lk}+{\displaystyle \frac{q^2}{\beta ^2}}\left(C_{ll}C_{kk}C_{lk}^2\right)`$

and Fourier transformation of (19) then gives

$`p(\xi _{lk};\theta )=e^{\beta \xi _{lk}C_{lk}/D_{lk}}f(\xi _{lk};D_{lk},B_{lk}).`$

Here, $`D_{lk}=C_{ll}C_{kk}C_{lk}^2`$, $`B_{lk}=\sqrt{C_{ll}C_{kk}}`$ and $`f(\xi ;D,B)`$ is an even function of $`\xi `$, defined by

$`f(\xi ;D,B)=\{\begin{array}{cc}{\displaystyle \frac{1}{\pi D^{1/2}}}K_0\left(B\xi /D\right),\hfill & \text{for GOE and}\hfill \\ \multicolumn{2}{c}{}\\ {\displaystyle \frac{1}{B}}\mathrm{exp}\left[2B|\xi |/D\right],\hfill & \text{for GUE,}\hfill \end{array}`$

(in which $`K_0(z)`$ is a modified Bessel function).

When $`l=k`$ we find in the GOE case that

$$p(\xi _{kk};\theta )=\frac{\xi _{kk}^{1/2}}{\sqrt{2\pi C_{kk}(\theta )}}e^{\xi _{kk}/2C_{kk}(\theta )},$$
(22)
which is just the Porter-Thomas distribution scaled so that $`\xi _{kk}=C_{kk}(\theta )`$, consistent with .

In order to test these predictions with reasonable statistics, we create ensembles by selecting states with phase angles restricted to subsets $``$ of the unit circle and compare their statistics with the accumulated distribution

$$p_{}(\xi _{lk})=\frac{1}{\mu ()}_{}p(\xi _{lk};\theta )d\theta ,$$
(23)
where $`\mu ()=_{}d\theta `$. These ensembles provide averages that are local in $`\theta `$ while allowing sufficient data to make meaningful numerical tests. The resulting statistical distributions are compared in Fig. 2 with the statistics of $`\xi _{02}=x_0x_2`$ and $`\xi _{22}=x_2^2`$. In each case the statistics in two subintervals are shown: one covering a region where there is scarring (larger than average overlaps) and one in a region where there is antiscarring (smaller than average overlaps). In all cases the Gaussian hypothesis works well in describing the whole distribution. Note that similar agreement has been found for models with GUE statistics (which are not shown here in the interests of brevity).

We also show, in Fig. 3, distributions for the individual overlaps $`x_k`$. The local distributions $`p(x_k;\theta )`$ are Gaussians scaled to that the variances are $`C_{kk}(\theta )`$, analogously to (3). We then define accumulated distributions by integrating over an interval $``$ with respect to $`\theta `$, as in (23). In practice, the most natural ensemble is formed by averaging over the whole unit circle, corresponding to taking $`=[0,2\pi )`$. A sequence of such distributions is shown in Fig. 3 corresponding to $`k=0,2,4,6`$. One sees that the distribution approaches the normal GOE prediction as $`k`$ increases but that for smaller values of $`k`$ the deviation from GOE is quite strong. Note that the Gaussian hypothesis describes quite accurately the statistics for moderate values of $`x_k`$ as well as in the tails which correspond to the strong scarring limit. There is some underestimation in the antiscarring limit where $`x_k`$ is small, but the overall agreement is good over the whole range of values of $`x_k`$. Similar or better agreement has been found in all the models we have examined.

All the numerical evidence we have collected indicates strongly that the Gaussian hypothesis accurately describes the complete statistics of chaotic eigenstates in a neighbourhood of a periodic orbit. While motivation for extending the linear theory of scarring to a complete basis has initially come from a desire to treat tunnelling, we note that the results in this section are interesting from the point of view of scarring alone. The covariance matrix $`C(\theta )`$ provides a means of characterising the eigenstate statistics in the neighbourhood of a periodic orbit in a way that is independent of the choice of test state. We can accomodate the characterisation of statistics using a basis other than $`\{|\stackrel{~}{k}\}`$ simply by altering the covariance matrix using a similarity transformation. Ideally, any basis we use should in the interests of compactness cover the immediate neighbourhood of the periodic orbit (as measured by Wigner functions, for example) with as few members as possible, but we are otherwise then free to manipulate those basis elements at will. In addition to being essentially basis-independent in this way, such a characterisation is complete, in the sense if we choose a subset of the basis which covers a semiclassical neighbourhood of the periodic orbit, the statistics of any other measure of the eigenstate near the periodic orbit could in principle then be determined from the corresponding sub-block of the covariance matrix.

Note finally that the deviations from RMT of the statistics of low-lying components of the eigenstate persist, and in fact remain fixed, in the semiclassical limit. This can be understood intuitively to derive from the fact that as $`\mathrm{}0`$ (or $`N\mathrm{}`$) the support of a harmonic oscillator basis state shrinks at the same rate as the semiclassical footprint of the periodic orbit. In $`\mathrm{\Sigma }`$, the area occupied by each decays as $`O(\mathrm{})`$ in the semiclassical limit. This semiclassical shrinking of the test states we use allows the scarring effect to remain fixed while other measures of scarring (such as the value of a wavefunction at a fixed point in space) lead to a vanishing effect in the semiclassical limit.

## IV Tunnelling statistics with scarring

The scarred wavefunction statistics described in the previous section are now used to derive distributions for the anomalous tunnelling rate statistics of . The first step is to outline how the discussion for maps in the preceding section may be adapted to the statistics of states in potential-well problems and then we can describe how these ideas are used to construct distributions for tunnelling rates.

While the formalism of tunnelling is presented in  in terms of quantum maps and their eigenstates, quantisation of a potential-well problem leads more naturally to wavefunctions in a full Hilbert space and at first sight it is not obvious how to relate the two descriptions. Faced with a numerically-produced wavefunction, for example, there is no exact, purely quantum-mechanical way of defining a corresponding state in a finite-dimensional Hilbert space; all existing descriptions are semiclassical in nature. Our task is made easier by the fact that tunnelling depends only on properties of the state in a small region of phase space (or of a Poincaré section in the reduced description) surrounding a tunnelling orbit. Within that region, the dynamical properties of the tunnelling orbit have been used in  to give an explicit prescription for the definition of a Wigner function on a Poincaré section. We will not repeat the details of that prescription here but state simply that a recipe exists for turning an eigenfunction $`\psi _n(x,y)`$ of a Hamiltonian with potential $`V(x,y)`$ into a Wigner function $`𝒲_n(y,p_y)`$ on a Poincaré section $`\mathrm{\Sigma }`$ defined by $`x=x_0`$ and on which $`(y,p_y)`$ are canonical coordinates. This Wigner function is invariant near the tunnelling orbit under iteration of a semiclassical quantisation of the Poincaré map . With appropriate normalisation this construction allows us to compute overlaps, using

$$\left|n|\stackrel{~}{k}\right|^2=_\mathrm{\Sigma }dydp_y𝒲_n(y,p_y)\stackrel{~}{𝒲}_k(y,p_y),$$
(24)
between a reduced chaotic state $`|n`$ and an eigenstate $`|\stackrel{~}{k}`$ of a harmonic oscillator Hamiltonian $`h(y,p_y)`$ defined on the Poincaré section $`\mathrm{\Sigma }`$.

We now check that these overlaps are described statistically by the distributions derived in the preceding section. Until now the discussion has been given entirely in terms of simple quantum maps for which the chaotic spectrum is described in terms of eigenangles on the unit circle. In physical problems the spectrum is more naturally described in terms of energy levels on the real line and we must connect the two pictures before proceeding. A detailed connection can be made using the transfer-operator approach of , but for present purposes it is sufficient to identify the part of the spectrum corresponding to a narrow segment of the spectral envelope (so that we may partition states into ensembles with simple Gaussian overlap statistics). A continuous-time analog of the phase angle $`\theta `$ is provided by the action $`S(E)`$ of the periodic orbit responsible for scarring . We therefore partition the quantum states into ensembles for which $`\theta =S(E)/\mathrm{}`$ lies in narrow intervals, mod $`2\pi `$. We assert that the statistics of the overlaps in (24) for each such ensemble should be a nonisotropic Gaussian such as described in (16).

We test this assertion using the two-dimensional double-well potential

$`V(x,y)=(x^21)^4+x^2y^2+\mu y^2+\nu y+\sigma x^2y`$

whose dynamics are predominantly chaotic near the saddle energy for appropriate parameter values and which we will later use to test tunnelling rate statistics. In practice we form statistical ensembles by fixing the energy and quantising $`\mathrm{}`$ and in this way classical-dynamical parameters are held fixed throughout the spectrum. Our system does not allow the accumulation of sufficient data to test the distributions meaningfully for a fixed value of $`\theta =S(E)/\mathrm{}`$ and we must instead pass straight to a distribution in which all states are combined, of the kind illustrated in Fig. 3. We show in Fig. 4$`(a)`$ the accumulated distribution $`p(x_0)`$ of overlaps with a harmonic oscillator ground state (with $`k=0`$) for the parameter values $`(\mu ,\nu ,\sigma )=(0.1,0,0)`$. The potential in this case has a symmetry of reflection in $`y`$ and therefore has a real periodic orbit connecting to a tunnelling orbit in the $`x`$-axis. The harmonic oscillator used is the one which generates the monodromy matrix $`W`$ of the tunnelling orbit and from which the tunnelling operator is constructed. The distribution derived from the Gaussian hypothesis describes the quantum data well and indicates that it is a reasonable basis from which to tackle tunnelling-rate statistics.

To describe tunnelling-rate statistics, it is convenient to begin once again with an ensemble taken from a fixed part of the spectral envelope (with $`\theta =S(E)/\mathrm{}`$ fixed mod $`2\pi `$). Repeating the derivation of (19), with $`B_{lk}`$ relaced by the matrix $`T`$ representing the scaled tunnelling operator, we obtain

$`p(y;\theta )={\displaystyle \frac{1}{2\pi }}{\displaystyle _{\mathrm{}}^{\mathrm{}}}e^{iqy}\stackrel{~}{p}(q,\theta )dq,`$

where

$$\stackrel{~}{p}(q,\theta )=\left[\mathrm{det}\left(I\frac{2iq}{\beta }TC(\theta )\right)\right]^{\beta /2}.$$
(25)
Our numerical calculations are made for GOE systems, in which case $`\beta =1`$. (Note that while the mechanics of the calculation are similar to those leading to (19), the quantities under investigation are very different in each case.) This distribution describes an ensemble of tunnelling rates taken with a fixed value of $`\theta `$. As was the case for overlap statistics, it is not possible with our present model to produce enough data to test the distribution for narrow ranges of $`\theta `$ and instead we make a comparison in which states from all parts of the spectral envelope are combined. Theoretically, we produce a combined distribution describing such collective statistics by averaging $`p(y,\theta )`$ over $`0\theta <2\pi `$.

Comparison of the resulting distribution with quantum data for the potential $`V(x,y)`$, with the same parameter values as in Fig. 4$`(a)`$, is shown in Fig. 4$`(b)`$ as a solid curve. Agreement with the scarred distribution is good and there is marked deviation from the prediction of standard RMT (short-dashed curve). Also shown for comparison is the Porter-Thomas distribution  (long-dashed curve).

Systems with symmetry such as that used for Figs. 4$`(a)`$ and 4$`(b)`$ are expected generically to show scarred statistics of the kind described in this paper. Nonsymmetric systems may also have scarred statistics, however, if parameters are tuned so that the real extension of the tunnelling orbit is periodic. We show tunnelling rate statistics in Figs. 4$`(c)`$ and 4$`(d)`$ for one such potential, corresponding to the parameter values $`(\mu ,\nu ,\sigma )=(0.25,0.40,0.254)`$. Quantum data for this potential were used in  to illustrate that dynamics could affect tunnelling statistics and the present theory gives a quantitative explanation of the qualitative observations made there. In Fig. 4$`(c)`$ we compare the scarred $`p(y)`$ (solid curve) with the standard RMT prediction (short-dashed curve) and the Porter-Thomas distribution (long-dashed curve).

Note that while the scarred distribution appears not to differ strongly from the Porter-Thomas distribution in Figs. 4$`(b)`$ and 4$`(c)`$ , the difference is significant when plotted on a logarithmic scale. This is done for the nonsymmetric potential in Fig. 4$`(d)`$ (and omitted for the symmetric potential in the interests of space). We also note that in the case of the symmetric potential the first hundred or so states (omitted from Figs. 4$`(a)`$ and 4$`(b)`$) give a distribution that seems better described by Porter-Thomas, though the scarred distribution unambiguously describes better the more highly excited states corresponding. This semiclassical limit seems to be approached more rapidly for the nonsymetric potential.

We note that while the detailed discussion has been presented here for two-dimensional systems, the general theory will also work in higher dimensions. In particular, for example, (25) holds except that the matrices $`T`$ and $`C(\theta )`$ have additional structure. Generalisation of $`T`$ to higher dimensions is straighforward (see ) and, for $`C(\theta )`$, the generation of the autocorrelation function as described in the appendices may also be extended, though there are now more than the three classical parameters $`(\mathrm{Tr}W,\mathrm{Tr}M_0,Q)`$, controlling quantum recurrence.

The success of scarring theory in describing the tunnelling-rate distribution therefore confirms the hypothesis of  that dynamical details of such low-dimensional systems can dominate the statistics of tunnelling and be responsible for strong deviation from RMT. The deviations we have seen are expected to be common in systems with discrete symmetry, such as the double well potential we have treated in Figs. 4$`(a)`$ and 4$`(b)`$. Deviation will also be seen in more generic systems, however, as parameters are tuned to make the real extension of the tunnelling orbit periodic . Such sensitive dependence of experimentally-accessible data on system parameters and dynamics might be useful, for example, as a means of probing the internal dynamics of relatively complex systems.

## V Conclusion

We have shown that wavefunction statistics in the neighbourhood of a fixed point or short periodic orbit are well described by a normal distribution which allows for correlations between wavefunction components along different basis elements. We have also shown that this wavefunction distribution can be used to describe scarring anomalies previously detected in tunnelling-rate statistics.

Scarring effects lead to strong deviations from the standard distributions of RMT in the neighbourhood of short periodic orbits. To characterise these deviations fully, we examined the components of the eigenstates of a chaotic map in the eigenbasis of a harmonic oscillator centred on a periodic orbit. The low-lying states of such an eigenbasis are localised around the corresponding periodic orbit and overlaps with a chaotic wavefunction are strongly affected by scarring. In appropriate subintervals of the chaotic spectrum, we have proposed that the joint probability distribution describing the corresponding components of chaotic eigenstates remain normal as in standard RMT but has a nontrivial covariance matrix (whose elements are completely determined by linearised dynamics around the periodic orbit). This hypothesis suffices to describe accurately all the measures of deviation from RMT that we have examined. The collective statistics of chaotic wavefunctions taken from a complete chaotic spectrum are described by superimposing these distributions. Deviations from RMT remain strong (and are effectively reproduced by the theory) in these combined ensembles, which are more natural to work with from a practical point of view.

We believe that examining the wavefunction statistics in a complete basis such as this provides a promising technique for understanding scarring in general terms. The description is essentially basis-independent and also completely determines the nature of wavefunction statistics near a periodic orbit.

The assumption of normal deviation from RMT also allows the effect of scarring to be incorporated very naturally into the calculation of tunnelling rate statistics. The theory was tested on a two-dimensional double-well potential with chaotic dynamics and found to describe the statistics of its eigenfcuntions and energy-level splittings well. Scarring can account for quite strong deviation from standard RMT in this system. In this way we might use tunnelling and the statistics of tunnelling to probe sensitively the internal dynamics of unstable systems.

Acknowledgements
SC Creagh and SY Lee are supported by the EPSRC under the Fast Stream scheme. SY Lee also acknowledges support from KOSEF.

## A generating the correlation function

We calculate the correlation function $`A_{lk}(t)`$ defined in section III using a generating function

$`G(w,z,t)=\mathrm{Tr}e^{w\widehat{a}}\widehat{𝒯}(z)\widehat{U}_{\mathrm{lin}}^t,`$ (A1)
where $`\widehat{a}`$ is an annihilation operator for the harmonic oscillator $`\widehat{h}`$ and

$`\widehat{𝒯}(z)=e^{(\mathrm{ln}z)\widehat{h}/\mathrm{}},`$

reduces to the tunnelling operator when $`z=e^{\alpha _0}`$. Each of the three operators in this sequence is the quantisation of either a linear symplectic transformation or a phase-space translation, possibly complex. The unitary operator $`\widehat{U}_{\mathrm{lin}}^t`$ quantises a real symplectic matrix $`M_0^t`$ which linearises the real unstable motion near $`\zeta _0`$. The operator $`\widehat{𝒯}(z)`$, which is nonunitary, quantises the complex symplectic matrix $`W(z)`$ defined in (18). Finally, the operator $`e^{w\widehat{a}}`$ can be interpreted as a translation along the complex phase-space displacement

$$\delta \zeta =w\sqrt{\frac{\mathrm{}}{2}}\left(\begin{array}{c}1\\ i\end{array}\right).$$
(A2)
The combined operator therefore quantises the affine transformation

$`\zeta W(z)M_0^t\zeta +\delta \zeta .`$

Because the transformation is affine, the trace formula may be used to evaluate the right hand side (A1) exactly. It is shown in appendix B that this gives

$$G(w,z,t)=\frac{1}{\sqrt{\mathrm{\Delta }(z,t)}}\mathrm{exp}\left[(1)^{\mu t}\frac{iw^2}{2}\frac{z}{\mathrm{\Delta }(z,t)}\mathrm{tan}\psi (t)\right]$$
(A3)
where

$`\mathrm{\Delta }(z,t)=\mathrm{Tr}W(z)M_0^t2=(1)^{\mu t}\left(m(t)z+{\displaystyle \frac{m^{}(t)}{z}}\right)2`$

and $`m(t)=e^{i\varphi (t)}\mathrm{sec}\psi (t)`$ as described in Sec. III. Note that $`G(w,z,t)`$ is independent of $`\mathrm{}`$. Note also that since $`\mathrm{\Delta }(z,t)`$ is complex, care must be taken in deciding the branch of the square root in (A3) — the branch is chosen using the Maslov index (computed from real dynamics), so that $`\mathrm{\Delta }(z,t)`$ reduces to the standard formula for real maps when $`z1`$.

Expanding the right hand side of (A3) gives

$`G(w,z,t)={\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}(1)^{n\mu t}{\displaystyle \frac{i^nw^{2n}}{2^nn!}}{\displaystyle \underset{k=0}{\overset{\mathrm{}}{}}}z^{k+2n+1/2}\mathrm{sin}^n\psi \sqrt{\mathrm{cos}\psi }C_k^{n+1/2}(\mathrm{cos}\psi )e^{i(n+k+1/2)(\varphi \mu t\pi )},`$

where $`C_k^{n+1/2}(\mathrm{cos}\psi )`$ are Gegenbauer polynomials defined by the generating function

$`{\displaystyle \frac{1}{(12xz+z^2)^{n+1/2}}}={\displaystyle \underset{k=0}{\overset{\mathrm{}}{}}}C_k^{n+1/2}(x)z^k.`$

Comparing this with the expansion

$$G(w,z,t)=\underset{k=0}{\overset{\mathrm{}}{}}\underset{m=0}{\overset{\mathrm{}}{}}\frac{w^m}{m!}z^{k+m+1/2}\stackrel{~}{k}|\widehat{a}^m\widehat{U}_{\mathrm{lin}}^t|\stackrel{~}{k}$$
(A4)
obtained by evaluating the trace (A1) in an eigenbasis of $`\widehat{h}`$ leads us to the result

$`\stackrel{~}{k}|\widehat{a}^{2n}\widehat{U}_{\mathrm{lin}}^t|\stackrel{~}{k}={\displaystyle \frac{(2n)!}{2^nn!}}\mathrm{sin}^n\psi \sqrt{\mathrm{cos}\psi }C_k^{n+1/2}(\mathrm{cos}\psi )e^{i(k+1/2)(\varphi \mu t\pi )+in\varphi +in\pi /2}`$

from which we determine $`A_{k+2n,k}^{\mathrm{lin}}(t)`$ as given in (17). Note that we find $`A_{lk}^{\mathrm{lin}}(t)=0`$ unless $`l`$ and $`k`$ are both even or both odd. This is because the linearised system enjoys an inversion symmetry even when the original system does not.

## B tracing the affine transformation

In this appendix we derive (A3) using the trace formula . We consider first the general affine transformation

$$\zeta M\zeta +Z$$
(B1)
where $`M`$ and $`Z=(a,b)`$ are an arbitrary symplectic matrix and phase-space displacement respectively. For convenience we may initially suppose that $`M`$ and $`Z`$ are real and generalise to complex transformations later using analytic continuation.

The transformation in (B1) has a single fixed point

$`\zeta =(MI)^1Z.`$

The trace formula then gives

$`\mathrm{Tr}\widehat{T}(Z)\widehat{U}(M)={\displaystyle \frac{e^{iS_0/\mathrm{}}}{\sqrt{\mathrm{Tr}M2}}},`$

where $`\widehat{T}(Z)`$ and $`\widehat{U}(M)`$ are the operators which quantise $`\zeta \zeta +Z`$ and $`\zeta M\zeta `$ respectively and $`S_0`$ is the action of the fixed point. By evaluating the generating function for $`M`$ at the fixed point one can show that

$`S_0`$ $`=`$ $`{\displaystyle \frac{1}{2}}(bqap)`$
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{\Omega }(\zeta ,Z)`$
$`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{\Omega }(Z,(MI)^1Z),`$
where $`\mathrm{\Omega }(u,v)=u^TJv`$ is the symplectic form and $`\zeta =(q,p)`$.

We specialise to the trace in (A1) by substituting the $`\delta \zeta `$ defined in (A2) for $`Z`$ and $`W(z)M_0^t`$ for $`M`$, giving

$`G(w,z,t)={\displaystyle \frac{e^{iw^2\mathrm{\Xi }(z,t)/4}}{\sqrt{\mathrm{Tr}W(z)M_0^t2}}},`$

where

$$\mathrm{\Xi }(z,t)=\mathrm{\Omega }(\zeta _c,\frac{1}{W(z)M_0^tI}\zeta _c)\text{and}\zeta _c=\left(\begin{array}{c}1\\ i\end{array}\right).$$
(B2)
Notice that, because $`\delta \zeta `$ scales as $`\mathrm{}^{1/2}`$, there is no $`\mathrm{}`$-dependence in the final result.

We are free to work in a coordinate system in which $`K=I`$ and then

$`W(z)=e^{i\mathrm{ln}zJ}=z(I+iJ)+{\displaystyle \frac{1}{z}}(IiJ).`$

Consider first the case

$$M_0^t=(1)^{\mu t}\left(\begin{array}{cc}\mathrm{cosh}\rho t& e^\sigma \mathrm{sinh}\rho t\\ \\ e^\sigma \mathrm{sinh}\rho t& \mathrm{cosh}\rho t\end{array}\right).$$
(B3)
where $`\sigma `$ and $`\rho `$ are constants. Substitution in (B2) gives then

$`\mathrm{\Xi }(z,t)=2(1)^{\mu t}\mathrm{cosh}\sigma \mathrm{sinh}\rho t{\displaystyle \frac{z}{\mathrm{\Delta }(z,t)}},`$

where $`\mathrm{\Delta }(z,t)`$ was defined in appendix A. Explicit calculation of $`\mathrm{\Delta }(z,t)`$ for this $`M_0`$ gives

$`Q=\mathrm{sinh}\sigma .`$

Using the identity

$`\mathrm{sec}^2\psi (t)=\mathrm{cosh}^2\rho t+\mathrm{sinh}^2\sigma \mathrm{sinh}^2\rho t=1+\mathrm{cosh}^2\sigma \mathrm{sinh}^2\rho t`$

we then get $`\mathrm{tan}\psi (t)=\mathrm{cosh}\sigma \mathrm{sinh}\rho t`$ and hence (A3).

Any hyperbolic matrix can be put in the form (B3) following an orthogonal transformation which rotates coordinates by an angle $`\theta _0`$. The matrix $`W(z)`$ is unchanged by such a transformation. However, the annihilation operator is multiplied by a phase

$`\widehat{a}e^{i\theta _0}\widehat{a}`$

and the exponent of the generating function by

$`\mathrm{\Xi }(z,t)e^{2i\theta _0}\mathrm{\Xi }(z,t).`$

In systems with time-reversal symmetry we can demand that the eigenvectors $`|\stackrel{~}{k}`$ be real and this allows us to fix the possible rotation angles $`\theta _0`$ to multiples of $`\pi /2`$. Among these, only $`\theta _0=0`$ and $`\theta _0=\pi /2`$ give distinct results for $`\mathrm{\Xi }(z,t)`$. The formula we have given applies to the case $`\theta _0=0`$. In the case $`\theta _0=\pi /2`$, the formula holds as given if we rephase the eigenvectors according to $`|\stackrel{~}{k}(1)^k|\stackrel{~}{k}`$.

We can similarly use a convention for $`|\stackrel{~}{k}`$ such that (17) holds in the GUE case. This amounts in particular to adopting a convention in which the covariance matrix $`C_{lk}(E)=x_lx_k^{}`$ is real even though individual values of $`x_lx_k^{}`$ are complex. This reflects the fact that even in systems without global time-reversal symmetry, the linearised dynamics from which $`A_{lk}^{\mathrm{lin}}(t)`$ is calculated is time-reversal invariant and we are free to choose the eigenstates $`|\stackrel{~}{k}`$ in a way which reflects this. This is in direct analogy to the spatial symmetry mentioned in the previous appendix which the linearised dynamics also possesses even in the absence of an equivalent symmetry for the global dynamics.


text:
# Can spectral and local information decide separability?

## I Introduction

Entanglement has always been a key issue in the ongoing debate about the foundations and interpretation of quantum mechanics since Einstein  and Schrödinger  expressed their deep dissatisfaction about this astonishing part of quantum theory. Whereas for the long period from 1935 to 1964, until Bell  published his famous work, discussions about entanglement were purely meta-theoretical, nowadays quantum information theory has established entanglement as a physical resource and key ingredient for quantum computation and quantum information processing. This new point of view led to a dramatic increase of general structural knowledge about entanglement in the last few years. However, there are still many open problems and deciding them often fails at what we try to profit from in quantum computation: the rapid increase of dimensions. One of these unsolved problems is in fact the question whether a given quantum state is entangled or merely classically correlated, i.e. separable.

The present paper is devoted to the question to what extent the combination of spectral and local information can decide separability. We will show that the discriminating power of any separability criterion based on this kind of information is rather limited or even null for some otherwise simple classes of quantum states.

The paper is organized as follows. In Sec. II we will discuss entropic inequalities as common examples for the kind of criteria we want to address. We give a small counterexample for the conjectured monotonicity of the conditional Tsallis entropy and prove that the criteria based on conditional Tsallis and Rényi entropies are implied by the reduction criterion for any positive value of the entropic parameters. In Sec. III, IV we will construct pairs of isospectral states having the same spectra and the same reductions such that one state is entangled and the other is separable. Sec. III establishes such an isospectral and separable counterpart for any entangled Werner state in odd dimensions. Sec. IV then provides a similar example even for the case of two qubits by utilizing a phase gate operation showing clearly how weak the discriminating power of the considered class of separability criteria is.

To fix ideas we will start by recalling some of the basic notions and previous results.

A bipartite quantum state described by its density matrix $`\rho `$ acting on a Hilbert space $`=^{(1)}^{(2)}`$ is said to be separable, unentangled or classically correlated if it can be written as a convex combination of tensor product states

$$\rho =\underset{j}{}p_j\rho _j^{(1)}\rho _j^{(2)},$$
(1)
where the positive weights $`p_j`$ sum up to one and $`\rho ^{(i)}`$ describes a state on $`^{(i)}`$. This means in particular that pure states are separable if and only if they are product states. Moreover, all entanglement properties of pure states, which can always be written in their Schmidt form (cf. ) as $`|\mathrm{\Psi }=_i\sqrt{\lambda _i}|i|i,`$ are completely determined by the eigenvalues $`\{\lambda _i\}`$ of the reduced state $`\rho _A=\mathrm{tr}_B|\mathrm{\Psi }\mathrm{\Psi }|`$.

For mixed quantum states however, the situation is much more difficult and deciding whether a state is entangled or separable is not yet feasible in general. Currently, the most efficient necessary criterion for separability is the positivity of the partial transpose (PPT), i.e., the condition that $`\rho ^{T_1}`$ has to be a positive semi-definite operator . The partial transpose of the state is thereby defined in terms of its matrix elements with respect to some basis by $`kl|\rho ^{T_1}|mn=ml|\rho |kn`$. For the smallest non-trivial systems with $`2\times 2`$ resp. $`2\times 3`$ dimensional Hilbert spaces and a few other special cases the PPT-criterion also turned out to be sufficient . In higher dimensional systems however, so called bound entangled states exist, which satisfy the PPT-condition without being separable .

Another well known condition is given by the reduction criterion

$$\rho _A\mathrm{𝟏}\rho 0,$$
(2)
which is implied by the PPT-criterion but nevertheless an important condition since its violation implies the possibility of recovering entanglement by distillation (which is yet unclear for PPT violating states). For the case of two qubits the reduction criterion is also known to be sufficient . The general line of implication is

$$\rho \text{ sep.}\rho ^{T_1}0\rho \text{ undistillable}\rho _A\mathrm{𝟏}\rho 0.$$

## II Entropic inequalities

The idea to use entropic inequalities as separability criteria for mixed states goes back to the mid nineties when Cerf and Adami  and the Horodecki family  recognized that certain conditional Rényi entropies are non-negative for separable states, and it was recently resurrected by several groups  in the form of conditional Tsallis entropies. Characteristic for these criteria is, that they just use the spectra of the state and its reductions.

The quantum Rényi entropy depending on the parameter $`\alpha \text{ }\text{R}`$ is given by

$$S_\alpha (\rho )=\frac{\mathrm{log}\mathrm{tr}(\rho ^\alpha )}{1\alpha }.$$
(3)
For the case of separable states it was shown in :

$$S_\alpha (\rho )S_\alpha (\rho _A)0$$
(4)
for $`\alpha =0,\mathrm{}`$ and $`\alpha [1,2]`$, where $`S_0,S_1,S_{\mathrm{}}`$ reduces to the logarithm of the rank, the von Neumann entropy and the negative logarithm of the operator norm respectively.

In Ref.  essentially the same criterion was expressed in terms of the conditional Tsallis entropy

$$T_\alpha (\rho )=\frac{\mathrm{tr}(\rho _A^\alpha )\mathrm{tr}(\rho ^\alpha )}{(\alpha 1)\mathrm{tr}(\rho _A^\alpha )},$$
(5)
which was conjectured  to be monotone in $`\alpha `$ and non-negative for all positive values of $`\alpha `$ if and only if $`\rho `$ is separable.

However, it is definitively neither monotone nor sufficient for separability. Monotonicity can most easily be ruled out by low rank examples like

$$\rho =\frac{1}{2}\left(|\mathrm{\Phi }_+\mathrm{\Phi }_+|+|0101|\right),|\mathrm{\Phi }_+=\frac{1}{\sqrt{2}}\left(|00+|11\right),$$

for which the reduced state has eigenvalues $`\frac{1}{4},\frac{3}{4}`$ and therefore $`T_0=T_{\mathrm{}}=0T_2=\frac{1}{5}`$. Fortunately however, monotonicity is not necessary for proving the validity of the Tsallis/Rényi separability criteria for other values than $`\alpha =0,\mathrm{}`$, $`\alpha [1,2]`$ :

Proposition 1 Any state for which the reduction criterion (2) holds, in particular any separable state, satisfies

$$\mathrm{tr}(\rho _A^\alpha )\mathrm{tr}(\rho ^\alpha )\{\begin{array}{cc}0,& \alpha >1\\ 0,& \mathrm{\hspace{0.33em}0}\alpha <1\end{array}.$$
(6)
This implies the validity of the Tsallis/Rényi separability criteria for any positive value of the entropic parameter.

Proof :  For $`\alpha >1`$ the proof is essentially based on the Golden-Thompson inequality (cf.) stating that $`\mathrm{tr}\left(e^Ae^B\right)\mathrm{tr}\left(e^{A+B}\right)`$ for hermitian matrices $`A,B`$. Utilizing the definition of the reduced state, i.e., $`P0:\mathrm{tr}\left(\rho (P\mathrm{𝟏})\right)\mathrm{tr}\left(\rho _AP\right)`$ this leads to:

$`\mathrm{tr}\left(\rho _A^\alpha \right)`$ $`=`$ $`\mathrm{tr}\left[\rho (\rho _A^{\alpha 1}\mathrm{𝟏})\right]`$ (7)
$`=`$ $`\mathrm{tr}\left[\mathrm{exp}\left(\mathrm{ln}\rho \right)\mathrm{exp}\left((\alpha 1)\mathrm{ln}(\rho _A\mathrm{𝟏})\right)\right]`$ (8)
$``$ $`\mathrm{tr}\left[\mathrm{exp}\left(\mathrm{ln}(\rho )+(\alpha 1)\mathrm{ln}(\rho _a\mathrm{𝟏})\right)\right]`$ (9)
$``$ $`\mathrm{tr}\left[\mathrm{exp}\left(\alpha \mathrm{ln}\rho \right)\right]=\mathrm{tr}\left(\rho ^\alpha \right),`$ (10)
where Eq.(9)-(10) is implied by the reduction criterion (2) since the logarithm is operator monotone  and the exponential function is monotone under the trace. The latter can be seen by noting that for any $`A`$ hermitian, $`P0`$ and $`B=(A+ϵP)`$ with $`ϵ0`$:

$$\frac{}{ϵ}\mathrm{tr}\left(e^B\right)=\mathrm{tr}\left(e^BP\right)0.$$
(11)
Hence $`\mathrm{tr}\left(e^B\right)\mathrm{tr}\left(e^A\right)`$ is implied by $`BA`$.

For $`0\alpha <1`$ the reduction criterion can immediately be applied since $`f(A)=A^r`$ is an operator decreasing function for $`1r0,A0`$ (cf.) and thus

$$\mathrm{tr}\left(\rho _A^\alpha \right)=\mathrm{tr}\left[\rho (\rho _A^{\alpha 1}\mathrm{𝟏})\right]\mathrm{tr}\left(\rho ^\alpha \right),$$
(12)
which completes the proof.

For negative values of $`\alpha `$, $`T_\alpha (\rho )0`$ holds for any state no matter if it is entangled or not. This can be shown by noting that for $`\{|a\}`$ being an eigenbasis of $`\rho _A`$

$`\mathrm{tr}\left(\rho _A^\alpha \right)`$ $`=`$ $`{\displaystyle \underset{a}{}}a|\rho _A|a^\alpha ={\displaystyle \underset{a}{}}\left[{\displaystyle \underset{i}{}}ai|\rho |ai\right]^\alpha `$ (13)
$``$ $`{\displaystyle \underset{a,i}{}}ai|\rho |ai^\alpha \mathrm{tr}\left(\rho ^\alpha \right),`$ (14)
where Eq.(13-14) uses that $`\left(_ib_i\right)^\alpha _ib_i^\alpha `$ holds for $`b_i0,\alpha 0`$, and the last inequality is implied by the convexity of negative powers.

We note that the von Neumann case $`S_1(\rho _A)S_1(\rho )`$ resp. $`T_1(\rho )0`$ also follows from the reduction criterion due to the operator monotonicity of the logarithm. Hence, separability criteria based on conditional Tsallis-Rényi entropies for all positive values of the entropic parameter are implied by the reduction criterion and therefore certainly weaker than the latter. In particular they cannot be sufficient.

This raises the question about the resolution that separability criteria, which are based on the spectrum of a state and its reductions, can have in general. The next two sections are devoted to the somewhat disillusioning answer.

## III Separable counterparts for Werner states

Werner states  have always played an important and paradigmatic role in quantum information theory. Their characteristic property is that they commute with all unitaries of the form $`UU`$ and they can be expressed as

$$\rho (p)=(1p)\frac{P_+}{r_+}+p\frac{P_{}}{r_{}},0p1,$$
(15)
where $`P_+`$ ($`P_{}`$) is the projector onto the symmetric (antisymmetric) subspace of $`\text{ }\text{C}^d\text{ }\text{C}^d`$ and $`r_\pm =\mathrm{tr}[P_\pm ]=\frac{d^2\pm d}{2}`$ are the respective dimensions. Werner showed that these states are entangled iff $`p>\frac{1}{2}`$ independent of the dimension $`d`$. The following proposition shows however, that none of these entangled states for odd dimension $`d`$ can be detected by any separability criterion, which is based on the spectrum of the state and its reductions.

Proposition 2 Any entangled state in $`\text{ }\text{C}^d\text{ }\text{C}^d`$ with maximal chaotic reductions and eigenvalues having multiplicities which are multiples of $`d`$, has a separable isospectral counterpart, which is locally undistinguishable as it has the same reductions.

Proof:  Let us consider a special basis of maximally entangled states in $`\text{ }\text{C}^d\text{ }\text{C}^d`$ :

$$|\mathrm{\Psi }_{jk}=\frac{1}{\sqrt{d}}\underset{n=1}{\overset{d}{}}\mathrm{exp}\left(\frac{2\pi i}{d}jn\right)|n,nk,$$
(16)
where $`j,k=1,\mathrm{},d`$ and $``$ means addition modulo $`d`$. Any equal weight combination of all states of the form (16), which belong to the same value of $`k`$, is then a projector onto a separable state since

$`P_k`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{d}{}}}|\mathrm{\Psi }_{jk}\mathrm{\Psi }_{jk}|`$ (17)
$`=`$ $`{\displaystyle \frac{1}{d}}{\displaystyle \underset{j,n,m=1}{\overset{d}{}}}\mathrm{exp}\left[{\displaystyle \frac{2\pi i}{d}}j(nm)\right]|n,nkm,mk|`$ (18)
$`=`$ $`{\displaystyle \underset{n=1}{\overset{d}{}}}|nn||nknk|`$ (19)
is an equal weight combination of product states. Here we have used that $`\frac{1}{d}_{j=1}^d\mathrm{exp}\left[\frac{2\pi i}{d}j(nm)\right]=\delta _{n,m}`$. Moreover, the reductions of the respective states $`P_k/d`$ are maximally chaotic, i.e. $`\rho _A=\mathrm{𝟏}/d`$, just as the reductions of any maximally entangled state.

If we now have a state with multiplicities being multiples of $`d`$ we can replace the projectors onto its eigenspaces with sufficiently many projectors of the form $`P_k`$. The resulting state will then be again a convex combination of product states, i.e., separable, having the same spectrum and maximal chaotic reductions.

For the case of Werner states we note that the unitary invariance of the state $`\rho (p)`$ in Eq. (15) implies that its reductions are $`\rho _A=\mathrm{𝟏}/d`$. Moreover $`\rho (p)`$ has two eigenvalues $`(1p)/r_+`$ and $`p/r_{}`$ with multiplicities $`r_+`$, $`r_{}`$ which are indeed multiples of $`d`$ in odd dimensions.

Following Proposition 2 we can now construct a state

$$\rho ^{}(p)=\frac{(1p)}{r_+}\underset{k=1}{\overset{r_+/d}{}}P_k+\frac{p}{r_{}}\underset{l=1}{\overset{r_{}/d}{}}P_{l+r_+/d},$$
(20)
which has then both, the same spectrum and the same reductions as $`\rho (p)`$. However, as convex combination of separable states it is itself separable for any $`0p1`$.

## IV The two qubit case

For the case of two qubits the antisymmetric subspace is one dimensional and the Werner state $`\rho (p)`$ is therefore entangled iff it has an eigenvalue larger than one half. It was essentially shown by Bennett et al.  that this is indeed sufficient to decide separability for two qubit states with maximal chaotic reductions, such that these states cannot have separable counterparts.

However, the following will show that separability is in general not encoded in the spectrum of a two qubit state and its reductions.

It is often advantageous to express two qubit states in terms of the matrix

$$R_{ij}=\mathrm{tr}\left(\rho \sigma _i\sigma _j\right),i,j=0,\mathrm{},3$$
(21)
with $`\sigma _0=\mathrm{𝟏}`$ and $`\sigma _{1,2,3}`$ being the Pauli matrices. Normalization then requires $`R_{0,0}=1`$ and the reductions are determined by the vectors $`R_{0,i}`$ and $`R_{i,0}`$. However, positivity of the respective state has to be verified separately.

Let us now consider a controlled phase gate operation described by the unitary

$`U`$ $`=`$ $`|00|\mathrm{𝟏}+|11|\sigma _3`$ (22)
$`=`$ $`{\displaystyle \frac{1}{2}}(\mathrm{𝟏}+\sigma _3)\mathrm{𝟏}+{\displaystyle \frac{1}{2}}(\mathrm{𝟏}\sigma _3)\sigma _3`$ (23)
which applies a phase operation on the target bit if the source bit is in the state $`1`$ and leaves it unchanged if the source bit is in the state $`0`$. Controlled phase gates as well as controlled NOT gates can be considered as universal basic building blocks for future quantum computers since together with single qubit gates they can be used to implement an arbitrary quantum circuit (cf. ).

The operation $`U`$ interchanges the matrix elements $`R_{1,0}R_{1,3}`$, $`R_{0,1}R_{3,1}`$ and $`R_{1,1}R_{2,2}`$ since

$`U(\sigma _1\mathrm{𝟏})U^{}`$ $`=`$ $`(\sigma _1\sigma _3),`$ (24)
$`U(\mathrm{𝟏}\sigma _1)U^{}`$ $`=`$ $`(\sigma _3\sigma _1),`$ (25)
$`U(\sigma _1\sigma _1)U^{}`$ $`=`$ $`(\sigma _2\sigma _2).`$ (26)
If we now choose the state such that $`R_{0,1}=R_{1,0}=R_{3,1}=R_{1,3}=r`$ take on the same values, $`R_{1,1}=\frac{1}{2}`$ and the others vanish, then the controlled phase gate will leave the reductions unchanged and simply interchange $`R_{1,1}R_{2,2}`$:

$$R=\left(\begin{array}{c}\text{1}\text{r}\text{0}\text{0}\\ \text{r}& \frac{1}{2}& \text{0}& \text{r}\\ \text{0}& \text{0}& \text{0}& \text{0}\\ \text{0}& \text{r}& \text{0}& \text{0}\end{array}\right)\stackrel{U}{}R^{}=\left(\begin{array}{c}\text{1}\text{r}\text{0}\text{0}\\ \text{r}& \text{0}& \text{0}& \text{r}\\ \text{0}& \text{0}& \frac{1}{2}& \text{0}\\ \text{0}& \text{r}& \text{0}& \text{0}\end{array}\right)$$
(27)
These matrices correspond to positive states for $`|r|\frac{3}{8}`$. Moreover, since the partial transposition only changes the signs of the $`\sigma _2`$ components, the initial state is equal to its partial transpose and therefore separable on the whole interval. This argumentation fails however for the state $`\rho ^{}=U\rho U^{}`$, which corresponds to $`R^{}`$, and in fact the determinant of its partial transpose is negative and the state is thus entangled for $`\frac{\sqrt{3}}{8}<|r|\frac{3}{8}`$.

Hence, the controlled phase gate operation can produce an entangled state out of a separable one, which has the same reductions and due to the unitarity of the operation the same spectrum.

It is worth mentioning that although $`\rho ^{}`$ is entangled for $`\frac{\sqrt{3}}{8}<|r|\frac{3}{8}`$ it does not violate the CHSH Bell inequality . It can be shown that $`T_20`$ is indeed stronger than the Bell inequality for two qubit systems . Hence, a similar example with one state being separable and the other violating the CHSH inequality cannot exist. This is however no longer true for higher dimensional systems.

## V Conclusion

Deciding whether a given quantum state is entangled or not is still one of the big open problems in the theory of quantum information. We discussed to what extent criteria based on the spectrum of a state and its reductions can decide separability. We showed that criteria, which are based on entropic inequalities using conditional Tsallis and Rényi entropies, are definitely weaker than the reduction criterion but nevertheless valid for any value of the entropic parameter. Moreover, the fact that there exists a separable counterpart for any entangled Werner state in odd dimensions shows that there is only a rather loose relation between the spectra of a state and its reductions on the one hand and its entanglement properties on the other.

Since the latter must not depend on the choice of local bases, it is sufficient to look at local unitary invariants of a state in order to decide its separability. Surprisingly, the most efficient separability criterion is still based on just one rather simple invariant: the smallest eigenvalue of the partial transpose. One remaining question in this context would therefore be: how can other (easy calculable) invariants provide information about the separability of a state, which is not yet encoded in the smallest eigenvalue of its partial transpose?

## Acknowledgement

The authors would like to thank R.F. Werner for useful discussions. Funding by the European Union project EQUIP (contract IST-1999-11053) and financial support from the DFG (Bonn) is gratefully acknowledged.