text:
# A two dimensional position sensitive gas chamber with scanned charge transfer readout

## 1 Introduction

To monitor high rate beam fluency gas ionization chambers based devices  can be used. In order to obtain a 2d profile of radiation dose or fluency it is necessary to have a highly segmented anode  . This represents a major problem when the size of the chamber is large compared with the anode segmentation pitch. In the present work we describe a new and simple solution developed to give real two dimensional readout on a gas ionization chamber. The method consists on scanning electronically the detector active area in the perpendicular direction to the readout strips orientation (strips instrumented with the measuring electro-meters). Gated charge transfer is a well known method which has been used in different radiation detectors like, for example, gated gas-filled Time Projection Chambers  in High Energy Physics.

The transfer of electrons drifting from the detector volume to the collector electrode strips can be blocked by means of a local electric field that dominates the drift field in the proximity of the gate electrode strips (see simulation in Figure 1). In this way we provide a single stage control of the charge transfer from the gas drift region to the collector electrodes. In this work we profit from the great development of Micro–Pattern Gas Detectors (MPGD)  in the last years, specially the advances on kapton etching techniques. In our case we integrated a gating grid using a two metal layer circuit based in a kapton foil glued on a FR4 epoxy substrate to provide mechanical rigidity. The use of a 25$`\mu `$m thick kapton foil to produce a two level readout circuit allows to use a small gating voltage.

## 2 Prototype description

To design and optimize the geometry of the gate and collector electrodes we simulated the electric field map of a transfer gap cell through MAXWELL3D111Ansoft Corporate, Pittsburgh USA and MATLAB222The MathWorks, Inc. 3 Apple Hill Drive, Natick MA(USA) 01760-2098 programs. For the final design it was chosen a 1:2:1 ratio corresponding to the transfer gap height:gate–spacing:gate–width. With this geometry we expected to have a low electron transfer transparency (well bellow 1%) using a relatively low positive (+5V) ”closed gate” voltage333In the present work collector electrodes are always grounded. and moderate drift fields ($``$ $`10^5`$ V/m). A wider gate strip spacing would increment the electron transparency for a given negative ”open gate” voltage, but also higher positive voltages would be needed to block efficiently the electron drift. On the other hand the local electric field present around the gate strip either in open or closed mode should not affect significantly the electric field in the drift region. This is another important reason to use a small transfer gap between the gate and collector strips.

The readout circuit layout of this prototype is described in Figure 2. We used a 25$`\mu `$m thick kapton foil metal coated on both sides with 5$`\mu `$m copper. After etching of the copper, we obtained 35$`\mu `$m wide and 1 mm wide strips on the gate and collector plane respectively. The foil was glued on a 1mm FR4–epoxy substrate and subsequently kapton in the gap between gate strips was removed obtaining the Figure 2 micro-pattern. In order to provide a simple connection procedure a 64 pin 1/10 inch pitch standard connector was included in the board edge design for both gate and collector strips. This choice implies a 1.27 mm pitch in the detector active area. The collector strips were thus made 1mm wide with 1.27 pitch while 17 gate strips were joined in a group in the detector edge to achieve the previous effective pitch. The drift gas gap was built using a G10 6mm height frame glued on the readout board, and this volume was closed with a drift electrode made of 200$`\mu `$m thick G10 copper clad on one side.

## 3 Charge transfer

The experimental transparency for electrons obtained as a function of the gate voltage is shown in Figure 4. This measurements were done using pure Ar at atmospheric pressure, and corresponds to the ratio of the ionization current under X–ray irradiation measured at the collector and at the drift electrode. The plot shows that very small transparency values ($`10^3`$) can be obtained for low positive voltages at ”closed mode”. The transparency increases when applying negative voltages at the gating electrodes, showing that there is a transition region centered on -15V, and can reach values very close to 1 (”open mode”) when the gate voltage is below -20V. If the drift field is increased, the gate voltage has also to be increased proportionally to maintain the same transparency value. The fraction of current produced by electrons arriving to the gate electrode is the complementary value of the transparency calculated for the collector. We have chosen -17V as open gate voltage as this value can be commutated using standard CMOS analogue switches. If we consider the width of the gate strips $`w`$ (35 $`\mu `$m in Figure 2) and the corresponding space between gate strips $`c`$ (40 $`\mu `$m in Figure 2), the expected electron transparency $`t`$ is equal to the optical limit defined by

$$t=\frac{c}{c+w}=0.53\mathrm{at}V_g=g\frac{V_d}{h}=0.5V$$

being $`V_g`$ and $`V_d`$ the gate and drift voltages respectively, $`g`$ the height of the gate strips over the collector electrodes and $`h`$ the gas gap below the drift electrode. The measured maximum transparency significantly deviates from this expectation due to the trapezoidal etching of the kapton layer that partially covers the collector electrodes in the space between gate strips.

## 4 Two dimensional readout

To demonstrate the working principle of our position sensitive charge transfer scan readout method, a small area of the prototype was instrumented with charge integration electronics. To drive the detector, two separate electronic boards were used, one with the commutation circuit to sequentially commutate from closed to open gate voltage on the gate strips and another board with the necessary electronic readout channels to integrate the charge collected by each individual collector strip (integrators board). To integrate the collector current we used a precision switched IVC102  integrating amplifier (with an internal feedback capacitor set to 10 pF).

To sequentially change the voltage value applied to each gate strip a shift register was used, built by means of daisy chained D–type flip–flops (HEF40174B) whose outputs were connected to double analogue switches (AD7592DI). These analogue switches commute from closed gate voltage value (that keeps the electron transparency at lowest values), to a negative open gate voltage value giving a high electron transparency for the selected gate strip (typically $`>`$70% in our setup). A scheme of the readout circuit is shown on Figure 5. For the synchronous control of the two electronic boards we used a personal computer PCI embedded card444National Instruments PCI-6023E DAQ card..

The first prototype was instrumented with a total number of 96 effective pixels (on an active area of 15.2 $`\times `$ 10.2 mm2), corresponding to 12 gate channels and 8 collector channels. The total capacitance at the detector between gate and collector plane is 7nF, meaning a 2pF capacity per pixel. This gate–collector capacity gives a 2$`\mu `$s transient time during commutation, considering the 3k$`\mathrm{\Omega }`$ impedance of the IVC102 (with switches S1+S2 closed). This dead time is small compared with the typical integration times used (of the order of milliseconds). The total time required to obtain an image is equal to the integration time needed to integrate the charge transfered by each gate strip times the number of gate strips.

Figure 6 and 7 show two X–ray images obtained with this prototype using a Chromium X–ray tube: a 1.5 mm slit between two 5mm thick aluminum plates and a 5mm screw nut were illuminated. The $`closed`$ and $`open`$ gate voltage values used were +5V and -17V respectively, and the drift field applied in the image was 1.7$`\times `$104 V/m in pure Ar. To obtain this images a scan readout cycle time of less than 10 seconds was required, with a collector current value of 7 pA in the pixels with maximum signal.

The charge transparency value when the gate electrode is set to the closed  state should be ideally zero. But the real value differs from zero by a small amount ($`i.e.`$ 0.1%). This causes a leakage current, of a value proportional to the detector area, that can seriously distort the image. Nevertheless, if the leakage current per effective pixel is small we can correct this effect by what we have called differential readout. The correction term is measured by integrating the detector current during the same time interval used for the standard readout with all the gate electrodes closed.

This value will then be subtracted as a pedestal from the values obtained with the gate electrodes on opened state. In this way the image will not be dramatically affected by the non zero transparency value at the closed  gate state.

## 5 Conclusions

We have proved the working principle of a simple and reliable readout method for 2d position sensitive gas ionization detectors. This two dimensional charge transfer readout solution allows to cover large detector areas giving a high number of effective pixels and minimizing the number of readout electronic channels. By using a kapton–insulated two layer readout circuit, charge transfer does not affect the drift field and the control of the charge transfer process can be done with very low voltage values. The detector can be used not only with gases but with other photo-conductive media like non–polar liquids ($`i.e.`$ isooctane or tetrametilsilane). For low beam intensity or fast readout cycle applications, an intermediate gas avalanche device, like a Gas Electron Multiplier , can be added to improve its sensitivity. Considering moderate gas gains around 100, an image could be readout with a cycle period 100 times faster than an equivalent ionization chamber, leading to a fast beam imaging device.

## 6 Acknowledgments

We are grateful to Manuel Sanchez from CERN EST/DEM group for his permanent technical support and collaboration.


format: text

text:
# Chapter 1 Energetic Particles in Clusters of Galaxies and Starburst Galaxies

## Chapter 1 Energetic Particles in Clusters of Galaxies and Starburst Galaxies

### Abstract

The nonthermal characteristics of starburst galaxies and galaxy clusters are reviewed and the perspectives for their observation in high energy gamma rays discussed. It is argued that starburst galaxies should have physically similar and even simpler nonthermal properties as our Galaxy, the particle transport being convective in their strong galactic winds. They should be measurable in TeV gamma rays by large detectors, and in fact the recent detection of NGC 253 is a first demonstration. Galaxy clusters are equally expected to produce a nonthermal energetic nuclear particle component that gives unambiguous information about the entropy production history in clusters. Both types of objects should be main sources of and at the same time also major visible representatives for the ubiquitous nonthermal component of matter in the Universe.

### 1. Introduction

In this review I will first discuss the production of energetic particles in starburst (SB) galaxies and galaxy clusters. The resulting particle density and the energy spectrum are not determined by the sources alone, but depend to the same extent on the transport and confinement properties in such objects; especially in clusters of galaxies the question of confinement plays a decisive role. Therefore it is necessary to address the combined problem involving both the sources and the propagation (e.g. Völk et al. 1996). I shall in fact argue that, at least for strong sources and extended objects, energetic particle transport corresponds to simple limits. Using multi-wavelength information from the radio to the X-ray region, as well as recent $`\gamma `$-ray observations and theory for Supernova Remnants (SNRs) in our Galaxy, the expected energetic particle densities are very high, with a Cosmic Ray energy density in SB galaxies that is about two orders of magnitude higher than in our Galaxy. In clusters, the average nonthermal energy density corresponds roughly to that in the Interstellar Medium of the Galaxy. If we consider the enormous gas masses and physical sizes of rich galaxy clusters, we expect a sizeable overall $`\gamma `$-ray emission in the sub-TeV and TeV range which may be observable with the new generation of instruments, both on the ground and in space, despite the large distances and angular sizes of even nearby clusters in excess of hundred Mpc.

With this overall picture in mind we can attempt to interpret the observed nonthermal radiation. For clusters the physical characteristics of the emission in the synchrotron and Inverse Compton channels depend strongly on the strength of the magnetic field. Therefore I shall consider separately the cases of ’low’ ($`10^1\mu `$G) and ’high’ ($`1\mu `$G) field strengths, respectively. Even though radio observations indicate high Intracluster (IC) fields, this may not be true everywhere in a nonuniform IC medium. Finally, we shall ask ourselves about the perspectives for high energy $`\gamma `$-ray astronomy. The result is that nearby SB galaxies should have been marginally observable with the previous instruments like Whipple, CANGAROO I, HEGRA, or CAT. However, the next generation of instruments should be able to identify both SB galaxies and clusters of galaxies as TeV sources, even though, in the case of clusters, not by a large margin. Although not yet confirmed by another instrument, it is therefore gratifying that the CANGAROO II telescope appears to have successfully detected the nearby Southern Hemisphere SB galaxy NGC 253 in the TeV range.

### 2. Particle production and propagation in starburst galaxies

External galaxies should also produce Cosmic Rays (CRs), from sources similar to those in our Galaxy – plus possibly through other mechanisms, like overall accretion or an AGN. SB galaxies are assumed to be analogous to the Milky Way, although with a significantly enhanced star formation rate (SFR). This should lead to a higher CR production rate as well. However, the form and the strength of the galactic mass loss is probably different. Whereas in our own and similar spiral galaxies, the CR pressure drives a moderate mass loss (Breitschwerdt et al. 1991; Zirakashvili et al. 1996) which ultimately, at large distances, implies convective removal of all energetic particles produced in the disk (Ptuskin et al. 1997), the wind from SB galaxies is driven by the pressure of an extremely hot ($`10^8`$ K) component of the interstellar gas (e.g. Schaaf et al. 1989). Like the nuclear CR component, it can basically only cool by adiabatic expansion, resulting in wind velocities of thousands of km/sec. Therefore one must assume that beyond the dense (molecular) gas disk, CR propagation is essentially convective at all particle energies of interest for $`\gamma `$-ray astronomy.

The comparison with the Galaxy is interesting and important for several reasons:

1. To ascertain the contribution of recent star formation to CR production, compared to other processes (e.g. the action of the central Black Hole, of SNRs from old stars, and of X-ray binaries from the old stellar population)
2. To study the CR propagation process, i.e. the relative role of convection and diffusion, as indicated above
3. To test the theory of an universal CR origin (Burbidge and Brecher 1971) which would imply a universal CR intensity.

For definiteness I shall discuss two concrete examples, M82 in the Northern Hemisphere, and NGC 253, its ’twin’ in the Southern Sky. They are the closest SB galaxies at distances of $`3.2`$ Mpc (M 82) and $`2.5`$ Mpc (NGC 253).

#### 2.1. M82

M82 is the small sister of the large spiral M81 in the Northern Hemisphere, at a declination of 70 degrees. The high radio synchrotron intensity in the central region (galactocentric radius $`300`$ pc; height above disk midplane $`100`$ pc) of M82 suggests very high CR intensities. From minimum energy arguments we have a total CR energy density $`u_{CR}80`$ eV/cm3 and a magnetic field strength $`B50\mu `$G, assuming an electron to proton ratio of $`10^2`$, similar to that in the Galaxy. Thus, $`u_{CR}`$ and the magnetic field energy density $`u_B`$ may exceed their values in the Galaxy ($`1`$ eV/cm3) by about two orders of magnitude.

Also the far infrared (FIR) emission is large with an energy density $`u_{rad}200\mathrm{e}\mathrm{V}/\mathrm{cm}^3u_{rad}u_B`$ in the disk, presumably due to the absorption of the UV radiation from the many young and massive stars by interstellar dust in the high mean (molecular) gas density $`140`$ H-atoms/cm3 of the SB nucleus. This is consistent with the fact that M82 also obeys the radio/FIR-correlation that characterizes ”normal” star forming galaxies, i.e. galaxies whose emission is dominated by star formation and not by the radiation from an AGN (e.g. Wunderlich et al. 1987; Völk 1989). Observations suggest that the average magnetic field strength increases about as fast as the radiation field in starbursts (Lisenfeld et al. 1996).

The extended soft X-ray emission as well as the spectral break in the radio synchrotron spectrum at a close distance from the CR source region imply a fast, thermally driven galactic wind with a velocity of several 1000 km/sec (Chevalier and Clegg 1985). As mentioned before, this strongly suggests that beyond the immediate neighborhood of the SB nucleus CR transport is dominated by convection in this wind rather than by diffusion, in contrast to the situation in the lower halo of our Galaxy (Ptuskin et al. 1997). With an observed supernova rate $`\nu _{\mathrm{M82}}(SN)0.3/\mathrm{yr}`$, which is $`10`$ times the Galactic SN rate $`\nu _{\mathrm{Gal}}(SN)`$, it follows that this convective escape is consistent with SNRs being the sources for the CRs, both for M82 and the Galaxy, i.e. $`dE_{CR}/dt=A\nu (SN)`$, where the proportionality constant $`A`$ is the same for both galaxies, despite their enormously different characteristics (Völk et al. 1989). The question is then, whether $`A`$ is a universal constant for ”normal” galaxies. An answer can only be obtained by an investigation of other objects.

##### The role of Gamma-ray astronomy

However plausible such results may be, they can only furnish a necessary condition for the proof of the supposition that the CR sources are principally the same in all galaxies. Since the arguments presented depend on the assumption that the energy density of nuclear CRs exceeds that of synchrotron radiating CR electrons by a factor of 100, it is necessary to directly measure the $`\pi ^0`$-decay hadronic $`\gamma `$-rays in those objects.

The closest to such a measurement is a detection in TeV $`\gamma `$-rays which shows at least the presence of particles in roughly the same energy range, even if a separation between hadronic and leptonic $`\gamma `$-ray emission is still required. In fact, the leptonic Inverse Compton $`\gamma `$-ray spectrum at very high energies should be softer than the hadronic spectrum. The two contributions can thus be distinguished by a broad-band measurement of the $`\gamma `$-ray energy spectrum. The reason behind this is as follows: Assuming the CRs to be the result of diffusive shock acceleration in SNRs, then the CR nuclei will have a spectral index $`\alpha _p^{SNR}2`$. The primary CR electrons of high energies, on the other hand, will be radiatively cooled inside the SNR with a corresponding index $`\alpha _e^{SNR}3`$. This leads to an almost identical photon spectral index $`1`$ for $`\pi ^0`$-decay and Inverse Compton $`\gamma `$-rays from the sources themselves (Berezhko et al. 2002; Völk et al. 2002). After release from the sources and upon transport out of the galaxy the electrons will continue to cool radiatively in the strong radiation field and the interstellar magnetic field, so that their spatially integrated energy spectrum will have an even softer index $`\alpha _e^{diff}4`$ leading to a spectral index of their diffuse Inverse Compton $`\gamma `$-ray emission $`\alpha _{IC}^{diff}1.5`$. However, the spectral index of the diffuse hadronic $`\gamma `$-ray emission will still be $`1`$, assuming strictly convective CR transport. A second distinguishing feature should be the lower cutoff of the Inverse Compton $`\gamma `$-ray spectrum, due to the radiation loss-dominated cutoff of the electron acceleration process in the sources.

I should add that this last consideration is modified when we go to sufficiently low electron energies, where radiative cooling inside the SNR sources does not occur during their life time. This should be the case at $`\gamma `$-ray energies of some GeV. Below those energies the diffuse Inverse Compton spectrum should flatten to an index $`\alpha _{IC}^{diff}1`$, or even somewhat below, making detection with space instruments at these energies more difficult than one might naively expect.

##### Gamma-ray detectability

What is then the $`\gamma `$-ray emission from M82 and what are the prospects for its detection? The SB nucleus has an angular size of about 40 arcsec. Such a region would clearly be a $`\gamma `$-ray point source. Even for the high molecular gas density Völk et al. (1996) have estimated a rather low hadronic $`\gamma `$-ray flux:

$`F_\gamma (E>100\mathrm{MeV})`$ $`=`$ $`1.3\times 10^9\mathrm{ph}\mathrm{cm}^2\mathrm{sec}^1`$
$`F_\gamma (>E)`$ $`=`$ $`2.6\times 10^{13}(E/1\mathrm{TeV})^{1.1}\mathrm{ph}\mathrm{cm}^2\mathrm{sec}^1,`$
where $`\alpha _p^{SNR}2.1`$ has been assumed. A slightly harder spectral index $`\alpha _p^{SNR}2.0`$ will not significantly change these numbers. Even though the gas density and the target photon field energy density in M82 are both about equally strongly increased relative to their Galactic values, the Inverse Compton flux might be somewhat higher, because its emission volume should be larger. However M82 has neither been detected by EGRET nor up to now in TeV $`\gamma `$-rays, in agreement with the above estimate. The HEGRA upper limit after 35 hours of observations with the stereoscopic system corresponds to 2 percent of the Crab flux. However, the future instruments in the Northern Hemisphere like VERITAS and MAGIC, as well as GLAST in the GeV range, should have a good chance to detect this galaxy.

#### 2.2. NGC 253

NGC 253 at a declination of -25 degrees has a projected optical size of $`27.5\times 6.8\mathrm{arcmin}^2`$. As such it would therefore expected to be an extended $`\gamma `$-ray source for the 0.1 degree angular resolution stereoscopic arrays CANGAROO III and H.E.S.S. in the Southern Hemisphere. The SB nucleus is very much smaller than this geometrical size of the stellar disk. The supernova rate is estimated to be between 0.1/yr and 0.3/yr (Bhattacharya et al. 1994, i.e. at least 3 to 10 times larger than for the Galaxy. The CR energy density $`u_{CR}`$ exceeds the Galactic value by a factor up to 100. A large radio synchrotron halo rises to a height of at least an optical disk radius (Carilli et al. 1992), consistent with a correlation of the soft X-ray with the H$`\alpha `$ emission over the same region above the disk (Strickland et al. 2002). The same phenomenon exists for M82 as well. Dynamically this again implies a thermally driven galactic wind, with an estimated mass loss rate of about $`24M_{}/yr`$ (Cappi et al. 1999). The far infrared (FIR) flux, the gas density in the SB nucleus, and the magnetic field strength are similar to the values in M82 (Blom et al. 1999). Indeed, the physical state of the two galaxies appears to be quite similar, as mentioned above. Goldshmidt and Rephaeli (1995) have used a specific wind model by Lerche and Schlickeiser (1980) to devise an Inverse Compton $`\gamma `$-ray emission model for this galaxy. More recent Chandra observations (Weaver et al. 2002) appear to show a heavily absorbed source of hard X-rays within the nuclear SB region, with an X-ray luminosity that exceeds the Eddington luminosity of a $`1.4M_{}`$ object within the nuclear SB region. The ionising source may be an intermediate-mass Black Hole, similar to a luminous compact X-ray source near the centre of M82 (Matsumoto et al. 2001). The absorbing optically thick molecular torus may collimate the SB-driven nuclear outflow. If this interpretation is correct, then these galaxies combine a starburst with AGN activity. As a consequence we must be careful not to over interprete the $`\gamma `$-ray results from these galaxies, because the hard X-ray source might itself be a source of $`\gamma `$-rays. Nevertheless, at the moment we will take the point of view that starbursts can be the dominant effect in many galaxies, like in the ones considered here, giving us the possibility to learn more about CR production in the Universe as a result of star formation alone.

##### Gamma-ray emission

Like for M82, EGRET has only obtained an upper limit for the $`\gamma `$-ray flux above a few 100 MeV (Blom et al. 1999). The OSSE instrument on CGRO, on the other hand, has marginally detected $`\gamma `$-ray emission below 1 MeV (Bhattacharya et al. 1994). Most recently, the CANGAROO collaboration (Itoh et al. 2002) has announced the detection of $`>0.5`$ TeV $`\gamma `$-rays with the CANGAROO II telescope at the high significance level of $`11\sigma `$. If confirmed this is an important result. The interpretation by Itoh et al. (2002) is in terms of Inverse Compton emission from a slightly extended source (for an estimated angular resolution of 0.23 degrees). The differential $`\gamma `$-ray spectrum was fitted with an expression $`E^{1.5}\mathrm{exp}(\sqrt{E}/a)`$, quite soft, in order to stay below the EGRET upper limit. This may be possible, even given the severe radiation losses of the electrons discussed before.

### 3. Clusters of galaxies

Galaxy clusters are the largest gravitationally bound objects in the Universe. In contrast to individual galaxies, the diffuse Intracluster (IC) gas dominates their visible mass by factors of a few compared to the stars in the constituent galaxies (e.g. Boehringer 2002). Rich clusters, our topic here, contain typically several hundred galaxies. The galaxy velocities are nevertheless so large that the main gravitational mass of clusters must be due to a Dark Matter. Due to Big Bang nucleosynthesis constraints this Dark Matter cannot be of a baryonic nature.

What makes galaxy clusters particularly interesting from a cosmological point of view, is the belief that they contain a representative sample of matter in the Universe, confining the primordial mixture during their formation. Typical temperatures of the hot IC gas amount to several keV, giving it as much or more kinetic and thermal energy as is contained in the galaxies. This gas is presumably representative of the larger part of the baryonic Dark Matter in the Universe.

Due to the enormous spatial dimensions of clusters and their turbulent state also the nonthermal particles should be confined over a Hubble time (Völk et al. 1996; see also Berezinsky et al. 1997). Ultrarelativistic nuclear particles do not suffer any significant losses during this time either. For this reason we have called the nuclear nonthermal component in clusters ”Cosmological Cosmic Rays”, because they have been stored in the clusters since the beginning of production with the onset of cosmic structure formation. The mass in energetic particles is of course negligible. On the other hand, their acceleration is inevitably associated with irreversible processes like shock dissipation. Whereas the thermal gas can cool radiatively, nuclear CRs do not. Observations of $`\pi ^0`$-decay $`\gamma `$-rays allow us therefore to study the entropy history of cluster formation in an unambiguous way.

From the chemical abundances of the IC gas, whose composition is about 40 percent Solar, a huge amount of star formation per galaxy must have occurred in evolved clusters like the Perseus and Coma clusters. This is in line with the dominance of giant elliptical and lenticular galaxies generally thought to be the result of interactions and mergers of (proto)spirals. The CR production in the large number of core collapse SNRs associated with these SBs must have been correspondingly large. The galactic winds should not only have adiabatically cooled and transported the particles into the IC medium but should also have strongly reaccelerated them in their termination shocks at the transition to the IC gas.

If one calculates the energy input into the IC gas from these processes, it follows that the main energisation of the IC medium, by a factor of more than 10, must have come from gravitational compression. As shown by Völk et al. (1996), this could in principle have come primarily through adiabatic $`pdV`$-work at contraction, after the main particle acceleration events as a result of the SBs, with little further entropy production. However this is an unlikely extreme, even though the most common shocks in the IC gas are not very strong, with typical Mach numbers of 4 – 5 (Miniati et al. 2000). Consequently, there is a comparably important possibility of (primary) particle acceleration in accretion shocks of subclusters onto the main cluster during ongoing cluster formation (e.g. Blasi 2002; Miniati 2002).

In addition to these ubiquitous and often large scale acceleration sites also the energy input from shock heating and acceleration by the large-scale jets from single FR II galaxies in clusters has been considered (e.g. Enßlin et al. 1997) and even argued to be dominant (e.g. Inoue and Sasaki 2001).

It will be important to distinguish the different processes. We believe that their joint action will lead to a multi-phase IC medium whose observational study is the key to a full understanding of the role of the various mechanisms that increase the entropy of the IC gas.

#### 3.1. Nonthermal radiation from clusters

Presumably nonthermal radiation from clusters has been observed as synchrotron emission at radio wavelengths (e.g. Ferretti and Giovannini 1996), and in the form of extreme UV (EUV) (e.g. Lieu et al. 1996; Berghöfer et al. 2000) and hard X-ray (Fusco-Femiano et al. 1999; Rephaeli et al. 1999) ’excesses’ (over a single-temperature plasma emission). No high energy $`\gamma `$-rays have been detected yet.

The IC magnetic field strengths in clusters are possibly very large, of the order of $`\mu `$G. They could originally be produced by dynamo action in starburst galaxies and subsequently be extended into the IC medium by galactic winds. The idea is old (Brecher and Burbidge 1972) and has been taken up independently by Völk and Atoyan (2000). At the present epoch the magnetic fields would be mostly cut off from the galaxies by magnetic reconnection and be sheared and possibly amplified by the IC flow field (Dolag et al. 1999).

The strength of the magnetic field is a decisive parameter for the interpretation of the radio, EUV and hard X-ray emissions. For ’low’ fields of the order of $`0.1\mu `$G, it is possible to interpret the hard X-ray excess – and also the EUV excess, observed in less than a handful clusters – as Inverse Compton emission by the radio synchrotron emitting electrons, assumed to be injected steadily into the IC medium for the last several billion years (Fusco-Femiano et al. 1999; Rephaeli et al. 1999). A break between the EUV and the hard X-ray regions could then be understood as a result of radiative losses.

For similar parameters and a mean gas density of $`n=10^3\mathrm{cm}^3`$, the extension of these spectra into the high energy $`\gamma `$-ray region is shown in Fig. 3, also including the electron Bremsstrahlung emission spectrum. There would be a good chance for GLAST to observe such emission.

However, this simple interpretation becomes untenable if the magnetic field is ’high’, i.e. of the order of $`1\mu `$G (Kim et al. 1990; Kronberg 1994), because then, for given radio intensity, the Inverse Compton emission is very small. For Coma it is then possible to understand the radio synchrotron emission as the result of a population of secondary electrons, produced by inelastic collisions of CR protons with IC gas nuclei and charged pion decay (Blasi and Colafrancesco 1999). Any Inverse Compton component responsible for the EUV excess must however be cut off before hard X-ray energies, in order not to overproduce radio synchrotron photons (Fig. 4). This is possible with a relic (from the early SB history) electron population that has been radiatively cooled to maximal Lorentz factors of the order of 100, while being subject to continuing gravitational compression (Sarazin and Lieu 1998). The hard X-ray excess remains unexplained in this context (Atoyan and Völk 2000).

##### Thermal origin of the hard X-ray ’excess’

This X-ray problem might have a very different solution than by nonthermal radiation! If we take into account that the SBs, forming the early type galaxies in the clusters, must have led to strong and fast galactic winds, with speeds of several 1000 km/sec, similar or in excess to those of present-day SB galaxies such as M82, then the dissipation in their termination shocks will result in thermal energies of the wind bubbles capable to produce thermal X-rays in the 10 keV range. Subsequent gravitational compression gives ample room for adiabatic heating to temperatures of the order of 100 keV. Such a ’very hot’ component of the IC gas would certainly survive until the present epoch, giving rise to a multi-phase IC medium in approximate pressure equilibrium. The relatively weak accretion and merger shocks can not produce such a component. Its formation would have to be accompanied by the acceleration of nuclear particles and electrons, so that one must expect radio synchrotron and thermal hard X-ray emissions to be correlated. Observations of this kind have been recently reported by Fusco-Femiano et al. (2003). Thus we suggest that the so-called ’excess’ hard X-ray emission is of thermal origin, as a relic of the early history of galaxy formation in the cluster.

#### 3.2. High energy $`\gamma `$-ray emission from clusters; observational perspectives

Assuming ’high’ magnetic fields in the cluster, and taking into account that the acceleration processes in the IC gas should produce a hadronic CR component as well, we can estimate the $`\gamma `$-ray emission. We shall do this for the case of the Coma cluster.

The relic electron component whose Inverse Compton emission is assumed to be responsible for the observed EUV ’excess’ is required to have an energy $`W_e5\times 10^{61}`$ erg presently. For a mean density $`N_{gas}10^3\mathrm{cm}^3`$ of the IC gas in Coma, this implies a large Bremsstrahlung flux up to about 200 MeV (Fig. 5). On the other hand, the Inverse Compton and Bremsstrahlung $`\gamma `$-ray fluxes produced by the recent secondary radio electrons are negligible. A total energy $`W_p3\times 10^{62}`$ ergs in energetic protons is required to produce these secondaries. This is less than 10 percent of the total thermal energy in the IC gas, and certainly about what one should expect for diffusive acceleration in the shocks heating the gas. These protons could indeed be ’old’, or ’cosmological’, since they do not cool or escape during a Hubble time. The hadronic $`\gamma `$-rays should dominate above the sharp Bremsstrahlung cutoff with a hard differential spectrum that should have a power law index of about 2. Although below the EGRET upper limit, the Bremsstrahlung and $`\pi ^0`$-decay $`\gamma `$-rays should be well within the reach of the future GLAST detector in the GeV range. At high energies, where no electronic emission can occur any more, ground-based arrays of imaging atmospheric Cherenkov telescopes like CANGAROO III, H.E.S.S., or VERITAS should be marginally able to detect objects like the Coma cluster in the TeV range, even though it is an extended source, with a radius in excess of 1 degree and with a flux level at 1 TeV of about 5 percent of the Crab Nebula.

### 4. References

Berezhko, E.G. et al., A&A, 395, 943 (2002)
Berezinsky, V.S. et al., ApJ, 487, 529 (1997)
Berghöfer, T. et al., ApJ, 535, 615 (2000)
Bhattacharya, D. et al., ApJ, 437, 173 (1994)
Blasi, P., astro-ph/0207361 (2002)
Balsi, P. and Colafrancesco, S., Astropart. Phys., 12, 169 (1999)
Blom, J.J. et al., ApJ, 516, 744 (1999)
Boehringer, H., Space Sci. Rev., 100, 49 (2002)
Brecher, K. and Burbidge, G.R., ApJ, 174, 253 (1972)
Breitschwerdt, D. et al., A&A, 245, 79 (1991)
Burbidge, G.R. and Brecher, K., Comments Astrophys. Space Phys., 3, 140 (1971)
Cappi, M. et al., A&A, 350, 777 (1999)
Carilli, C.L. et al., ApJ, 399, L59 (1992)
Chevalier, R.A. and Clegg, A.W., Nature, 317, 44 (1985)
Deiss, B.M. et al., A&A, 321, 55 (1997)
Dolag, K.M. et al., A&A, 348, 351 (1999)
Enßlin, Th. et al., ApJ, 477, 560 (1997)
Fusco-Femiano, R. et al., ApJ, 513, L21 (1999)
Fusco-Femiano, R. et al., A&A 398, 441 (2003)
Goldshmidt, O. and Rephaeli, Y., ApJ, 444, 113 (1995)
Inoue, S. and Sasaki, S., ApJ, 562, 618 (2001)
Itoh, C. et al., A&A, 396, L1 (2002)
Kim, K.-T. et al., ApJ, 355, 29 (1990)
Kronberg, P.P., Rep. Prog. Phys., 57, 325 (1994)
Lerche , I. and Schlickeiser, R., ApJ, 239, 1089
Lieu, R. et al, Science, 274, 1335 (1996)
Lisenfeld, U. et al., A&A, 314, 745 (1996)
Matsumoto, H. et al., ApJ, 547, L25 (2001)
Miniati, F. et al., ApJ, 542, 608 (2000)
Miniati, F., astro-ph/0206430
Ptuskin, V.L. et al., A&A, 321, 434 (1997)
Rephaeli, Y. et al., ApJ 511, L21 (1999)
Sarazin, C.L. and Lieu, R., ApJ, 494, 177 (1998)
Schaaf, R.. et al., A&A 336, 722 (1989)
Sreekumar, P. et al. ApJ, 464, 628
Strickland, D.K. et al., ApJ, 565, 689 (2002)
Völk, H.J., A&A, 218, 67 (1989)
Völk, H.J. et al., A&A, 213, L12 (1989)
Völk, H.J. et al., Space Sci. Rev., 75, 279 (1996)
Völk, H.J. et al., A&A, 396, 649 (2002)
Weaver, K.A. et al., ApJ, 576, L19 (2002)
Wunderlich, E. et al., A&AS, 69, 487 (1987)
Zirakashvili, V.L. et al., A&A, 311, 113 (1996)


text:
# The divergence of neighboring magnetic field lines and fast-particle diffusion in strong magnetohydrodynamic turbulence, with application to thermal conduction in galaxy clusters

## Abstract

Using direct numerical simulations, we calculate the rate of divergence of neighboring magnetic field lines in different types of strong magnetohydrodynamic turbulence. In the static-magnetic-field approximation, our results imply that tangled magnetic fields in galaxy clusters reduce the electron diffusion coefficient and thermal conductivity by a factor of $`510`$ relative to their values in a non-magnetized plasma.

The diffusion of fast particles in turbulent magnetized plasmas is important for fusion experiments, cosmic-ray propagation, and thermal conduction in galaxy-cluster plasmas fab94 ; bin81 ; ros89 ; tri89 ; chu93 ; tao95 ; pis96 ; cha98 ; mal01 ; jok73 ; ski74 ; nar01 ; gru02 . We consider particle diffusion in a static magnetic field, which is a reasonable first approximation for particles moving much faster than the $`𝐄\times 𝐁`$ velocity of field lines. The effects of turbulent bulk motions and field evolution have recently been considered by gru02 ; cha03 .

We consider magnetic fluctuations that possess an inertial range extending from an outer scale $`l_B`$ to a much smaller inner scale $`l_d`$ with the magnetic energy dominated by scales $`l_B`$. Except where specified, the discussion focuses on the case in which there is no mean magnetic field pervading the plasma. If a particle is tied to a single field line and travels a distance $`ll_B`$ along the static magnetic field, it takes $`l/l_B`$ random-walk steps of length $`l_B`$, resulting in a mean-square three-dimensional displacement of

$$(\mathrm{}x)^2=\alpha l_Bl,$$
(1)
where $`\alpha `$ is a constant of order unity. The values of $`\alpha `$ for the turbulence simulations used in this paper are listed in table 1. When there is a mean field $`𝐁_0`$ comparable to the rms field, $`(\mathrm{}x)^2`$ in equation (1) is interpreted as the mean-square displacement perpendicular to $`𝐁_0`$. If the particle’s motion along the field is diffusive with diffusion coefficient $`D_{}`$, then

$$l\sqrt{D_{}t},$$
(2)
and rec78 ; kro83

$$(\mathrm{}x)^2t^{1/2},$$
(3)
indicating subdiffusion \[$`Dlim_t\mathrm{}(\mathrm{}x)^2/6t0`$\].

In fact, a particle is not tied to a single field line. After traveling a short distance along the field, field gradients and scattering cause a particle to take a step comparable to its gyroradius $`\rho `$ across the magnetic field, from its initial field line, F1, to a new field line, F2. If the particle were to follow F2, it would separate from F1 because neighboring field lines tend to diverge. We call $`z_s`$ the distance the particle would have to follow F2 before separating from F1 by a distance $`l_B`$. (A particle typically separates from its initial field line by a distance $`l_B`$ after traveling a distance slightly less than $`z_s`$ since it drifts across the field continuously, but we ignore this effect in this paper.)

If a particle “takes a random step” of length $`mz_s`$ along the magnetic field, where $`m`$ is a constant of order a few, reverses direction, and then takes another random step of length $`mz_s`$ back along the field, it doesn’t return to its initial point. Part of the second step retraces part of the first step, but the remainder is uncorrelated from the first step. This loss of correlation leads to a Markovian random walk cha98 ; cha03 . When $`mz_sl_B`$, a single random step corresponds to a 3D displacement of

$$(\mathrm{}x)^2\alpha mz_sl_B.$$
(4)
When $`mz_s\lambda `$, where $`\lambda `$ is the Coulomb mean free path, a single step takes a time

$$\mathrm{}t\frac{m^2z_s^2}{D_{}}.$$
(5)
When $`mz_s`$ is only moderately greater than $`l_B`$ or $`\lambda `$, equations (4) and (5) remain approximately valid. During successive random steps, a particle will find itself in regions of differing magnetic shear, and thus $`z_s`$ will vary. The diffusion coefficient is given by $`D=(\mathrm{}x)^2/6\mathrm{}t`$ where $`\mathrm{}`$ is an average over a large number of steps cha43 . Ignoring factors of order unity,

$$DD_{}\frac{l_B}{L_\mathrm{S}}$$
(6)
as in rec78 ; cha98 ; cha03 , with

$$L_\mathrm{S}=\frac{z_s^2}{z_s}.$$
(7)
If there is a mean magnetic field comparable to the fluctuating field, equation (6) is recovered provided $`D`$ is replaced by $`D_{}`$, the coefficient of diffusion perpendicular to the mean field.

We treat field-line separation in a static field using magnetic-field data obtained from four direct numerical simulations of incompressible MHD turbulence. Key simulation parameters are given in table 1. The numerical method is described in mar01 . Each simulation uses Newtonian viscosity $`\nu `$ and resistivity $`\eta `$ and is run until a statistical steady state is reached. In simulations A1 and A2, the mean magnetic field is zero, and the initial magnetic field is dominated by large-scale fluctuations containing 10% of the maximum possible magnetic helicity. Turbulent fluctuations are sustained by non-helical random forcing of the velocity field, and the magnetic Prandtl number $`P_m\nu /\eta `$ equals 1. These simulations reach a statistical steady state with Kolmogorov-like kinetic and magnetic power spectra. In simulation B, $`P_m=75`$ and the field is amplified from an initially weak seed magnetic field by turbulent velocities sustained by non-helical random forcing. When the dynamo growth of the magnetic field saturates, the magnetic field is dominated by fluctuations on scales much smaller than the dominant velocity length scale, as in the simulations of marb . In simulation D, $`P_m=1`$ and the field is amplified from an initially weak seed magnetic field by turbulent velocities sustained by maximally helical random forcing. When the dynamo growth of the magnetic field saturates, the dominant magnetic field length scale is comparable to the dominant velocity length scale, as in the simulations of mara . “Simulation” A2rp is obtained by assigning each Fourier mode in simulation A2 a random phase without changing the mode’s amplitude. In simulations A1, A2, A2rp, and D we set $`\pi /l_B`$ equal to the maximum of $`kE_b(k)`$, where $`E_b(k)`$ is the power spectrum of the magnetic field \[the total magnetic energy is $`E_b(k)𝑑k`$\]. In simulation B we take $`\pi /l_B`$ to be the maximum of $`kE_v(k)`$, where $`E_v(k)`$ is the power spectrum of the velocity field. We set $`\pi /l_d`$ equal to the maximum of $`k^3E_b(k)`$.

We take a snapshot of the magnetic field in each simulation and introduce 2,000 pairs of field-line tracers whose initial separation $`𝐫_0`$ is perpendicular to the local field. We use linear interpolation to obtain the magnetic field between grid points and employ second-order Runge-Kutta to integrate field-lines. To improve the statistics in simulation A1, we use 2,000 field lines in each of seventeen snapshots of the field separated in time by an interval $`0.4l_B/u`$, where $`u`$ is the rms velocity. In simulations A2 and A2rp, we use 20,000 field lines in each of five snapshots of the magnetic field separated in time by an interval $`0.2l_B/u`$. We iteratively reduce the length step in the field-line integrations to achieve convergence.

The probability distribution function (PDF) of field-line separations in simulation A1 for $`r_0=2^9l_B`$ is plotted in figure 1 for different values of the distance $`z`$ a field-line pair is followed along the magnetic field. \[The probability is per unit $`\mathrm{log}_{10}(r/l_B)`$, so that the probability that $`r`$ lies in some interval is proportional to the area under the plotted curves.\] The PDF is highly non-Gaussian, and the tail of the distribution dominates the growth of $`r`$ when $`l_d<r<l_B`$. The maximum kurtosis $`r^4/r^2^2`$ is $`500`$ and occurs for $`z0.2`$. We find that decreasing $`r_0`$ increases the maximum kurtosis.

The growth of $`r`$ with $`z`$ in simulation A1 is shown in figure 2 for several values of $`r_0`$. There are three stages of growth similar to those described in previous theoretical treatments of field-line divergence nar01 ; ski74 : (1) an initial stage of exponential growth when $`r_0l_d`$, (2) approximate power-law growth with $`rz^a`$ for $`rl_B`$, and (3) $`rz^{1/2}`$, for $`r>l_B`$. However, some aspects of our results differ from previous studies. Within a single simulation, $`a`$ increases with decreasing $`r_0`$, which is probably related to the increasing prominence of the tail of the PDF as $`r_0`$ is decreased. In addition, stage 2 with $`rz^a`$ begins for $`r<l_d`$, perhaps because the field-line pairs with largest $`r`$, which dominate the growth of $`r`$, satisfy $`r>l_d`$ before $`r>l_d`$.

We seek to test the qualitative prediction of jok73 ; ski74 ; nar01  that $`z_s`$ and $`L_\mathrm{S}`$ asymptote to a value of order a few $`l_B`$ as $`r_0`$ is decreased towards $`l_d`$ in the large-$`l_B/l_d`$ limit. In figure 3, we plot $`z_s`$ for simulation A1 and simulation A2. The lower-resolution data of simulation A1 suggest the scaling $`z_s\mathrm{ln}(l_B/r_0)`$ for $`l_d<r_0<0.25l_B`$, in contradiction to the theoretical treatments. On the other hand, for simulation A2, the curve through the data is concave downward for $`l_d<r_0<l_B`$. Moreover, figure 3 shows that the simulation A1 data, and probably also the A2 data, have not converged to the high-Reynolds-number values of $`z_s`$ for $`l_B/16<r_0<l_B`$, values of $`r_0`$ that are within the inertial ranges ($`l_d`$ to $`l_B`$) of both simulations. Also, the slope $`dz_s/d(\mathrm{ln}(l_B/r_0))`$ for both $`r_0<l_d`$ and $`r_0<l_B/10`$ decreases significantly when $`l_B/l_d`$ is doubled. A comparison of the data for A1 and A2 thus suggests that in the large-$`l_B/l_d`$ limit $`z_s`$ asymptotes to a value of order several $`l_B`$ as $`r_0`$ is decreased towards $`l_d`$, as in jok73 ; ski74 ; nar01 . The same comments apply to the data for $`L_\mathrm{S}`$, which are plotted in figure 4.

We note that for $`r_0=l_B`$, $`z_s`$ and $`L_\mathrm{S}`$ are by definition 0. The numerical-simulation data points that appear to be plotted above $`l_B/r_0=1`$ actually correspond to $`r_0`$ just slightly smaller than $`l_B`$, indicating that $`z_s`$ and $`L_\mathrm{S}`$ are discontinuous at $`r_0=l_B`$ in the numerical simulations. The reason is that for $`r_0`$ just slightly less than $`l_B`$, some fraction of the field line pairs are initially converging and must be followed a significant distance before they start to diverge.

The thermal conductivity $`\kappa _T`$ in galaxy-cluster plasmas scales approximately like the diffusion coefficient of thermal electrons rec78 ; kro83 . For clusters, $`l_B/l_dl_B/\rho _i10^{13}`$, where $`\rho _i`$ is the proton gyroradius qua98 . In the large-$`l_B/l_d`$ limit, the numerical simulations and theoretical models indicate that $`L_\mathrm{S}`$ asymptotes to a value of order several $`l_B`$ as $`r_0`$ is decreased towards $`l_d`$, and $`L_\mathrm{S}`$ is not expected to increase appreciably as $`r_0`$ is further decreased from $`l_d=\rho _i`$ to $`\rho _e`$. Thus, $`L_\mathrm{S}(\rho _e)L_\mathrm{S}(l_d)`$. We take simulations A1 and A2 to be our best models of a galaxy-cluster magnetic field. We find that $`L_\mathrm{S}(l_d)11l_B`$ in A1 with $`l_B/l_d=23`$, and $`L_\mathrm{S}(l_d)1011l_B`$ in A2 with $`l_B/l_d=50`$. If $`l_B`$ is redefined so that $`2\pi /l_B`$ (instead of $`\pi /l_B`$) equals the maximum of $`kE_b(k)`$, then $`L_\mathrm{S}(l_d)7l_B`$ in A1 and $`L_\mathrm{S}(l_d)6.5l_B`$ in A2cha03 . Since the definition of $`l_B`$ contains an arbitrary constant of order unity, it is not clear which of the two definitions leads to the more accurate prediction of $`\kappa _T`$. Given this uncertainty, extrapolating these results to the large-$`l_B/l_d`$ limit suggests that $`L_\mathrm{S}(l_d)510l_B`$ in intracluster plasmas. Diffusion along the magnetic field can be suppressed by magnetic mirrors cha98 ; mal01  and wave pitch-angle scattering pis98 , but in clusters the Coulomb mean free path is sufficiently short that neither of these effects is very important. nar01 ; cha03  Thus, the parallel diffusion coefficient of thermal electrons is comparable to the thermal-electron diffusion coefficient in a non-magnetized plasma, $`D_0`$. In the static-magnetic-field approximation, equation (6) thus implies that $`D/D_00.10.2`$ for thermal electrons in intracluster plasmas, and that $`\kappa _T`$ is reduced by a factor of $`510`$ relative to the Spitzer thermal conductivity $`\kappa _\mathrm{S}`$ of a non-magnetized plasma. Heating of intracluster plasma from thermal conduction with $`\kappa _T0.10.2\kappa _\mathrm{S}`$ would be sufficient to balance radiative cooling in a some but not all clusters fab02 ; zak03 . More work is needed to determine the validity of the static-field approximation in clusters, to determine the factors of order unity that have been neglected in our phenomenological treatment of thermal conduction, and to clarify the effects of turbulent diffusion and turbulent resistivity on heat conduction in strong MHD turbulence.

We thank Steve Cowley and Alex Schekochihin for valuable discussions. This work was supported by NSF grant AST-0098086 and DOE grants DE-FG02-01ER54658 and DE-FC02-01ER54651 at the University of Iowa, with computing resources provided by the National Partnership for Advanced Computational Infrastructure and the National Center for Supercomputing Applications.


text:
# Level statistics of systems with infinitely many independent components based on the Berry-Robnik approach

## I Introduction

An important property of quantum-classical correspondence appears in the statistical property of energy levels of bounded quantum systems in the semiclassical limit. Universal behaviors are found in the statistics of unfolded energy levels at a given interval, which are the sequence of numbers uniquely determined by the energy levels using the mean level density obtained from the Thomas-Fermi rule. It is widely known that, for quantum systems whose classical counterparts are integrable (those systems will be referred to as classically integrable systems), the distribution of nearest-neighbor level spacing is characterized by the Poisson (exponential) distribution, while for quantum systems whose classical counterparts are strongly chaotic, the quantal level statistics are well characterized by the random matrix theory which gives level-spacing distribution obeying the Wigner distribution.

Level statistics for the integrable quantum systems has been theoretically studied by Berry-Tabor, Sinai, Molchanov, Minami, Bleher, Connors and Keating, and Marklof, and have been the subject of many numerical investigations. Still its mechanism is not well understood, the appearance of the Poisson distributions is now widely admitted as a universal phenomenon in generic integrable quantum systems.

As suggested, e.g., by Hannay (see the discussion of ), one possible explanation would be as follows: For an integrable system of $`f`$ degrees-of-freedom, almost every orbit is generically confined in each inherent torus, and the whole region in the phase space is densely covered by invariant tori as suggested by the Liouville-Arnold theorem. In other words, the phase space of the integrable system consists of infinitely many tori which have infinitesimal volumes in Liouville measure. Then, the energy level sequence of the whole system is a superposition of sub-sequences which are contributed from those regions. Therefore, if the mean level spacing of each independent subset is large, one would expect the Poisson distribution as a result of the law of small numbers. This scenario suggested by Hannay is based on the theory proposed by Berry and Robnik.

The Berry-Robnik theory relates the statistics of the energy level distribution to the phase-space geometry by assuming that the sequence of the energy spectrum is given by the superposition of statistically independent subspectra, which are contributed respectively from eigenfunctions localized onto the invariant regions in phase space. Formation of such independent subspectra is a consequence of the condensation of energy eigenfunctions on disjoint regions in the classical phase space and of the lack of mutual overlap between their eigenfunctions, and, thus, can only be expected in the semi-classical limit where the Planck constant tends to zero, $`\mathrm{}0`$. This mechanism is sometimes referred to as the principle of uniform semi-classical condensation of eigenstates. This principle states that the Wigner function of a semiclassical eigenstate is connected on a region in phase space explored by a typical trajectory of the classical dynamical system. In integrable systems, the phase space is folded into invariant tori, and the Wigner functions of the eigenstates tend to delta functions on these tori in the semiclassical limit. On the other hand, in a strongly chaotic system, almost all trajectories cover the energy shell uniformly, and hence the Wigner functions of eigenstates are expected to become a delta function on the energy shell as suggested by the quantum ergodicity theorem. Because of the suppression of the tunneling, each quantum eigenstate is folded into independent subsets in the semiclassical limit, and is expected to form independent spectral components. Indeed, formation of such independent components are checked numerically in a deep semi-classical regime.

In the Berry-Robnik approach, the overall level spacing distribution is derived along a line of mathematical framework by Mehta, as follows: Consider a system whose classical phase space is decomposed into $`N`$-disjoint regions. The Liouville measures of these regions are denoted by $`\rho _i(i=1,2,3,\mathrm{},N)`$ which satisfy $`_{i=1}^N\rho _i=1`$. Let $`E(S)`$ be the gap distribution function which stands for the probability that an interval $`(0,S)`$ contains no level. $`E(S)`$ is expressed by the level spacing distribution $`P(S)`$ as follows;

$$E(S)=_S^{\mathrm{}}𝑑\sigma _\sigma ^{\mathrm{}}P(x)𝑑x.$$
(1)
When the entire sequence of energy levels is a product of the statistically independent superposition of $`N`$ sub-sequences, $`E(S;N)`$ is decomposed into those of sub-sequences, $`E_i(S;\rho _i)`$,

$$E(S;N)=\underset{i=1}{\overset{N}{}}E_i(S;\rho _i).$$
(2)
In terms of the normalized level spacing distribution $`p_i(S;\rho _i)`$ of a sub-sequence, $`E_i(S;\rho _i)`$ is given by

$$E_i(S;\rho _i)=\rho _i_S^{\mathrm{}}𝑑\sigma _\sigma ^{\mathrm{}}p_i(x;\rho _i)𝑑x,$$
(3)
and $`p_i(S;\rho _i)`$ is assumed to satisfy

$$_0^{\mathrm{}}Sp_i(S;\rho _i)𝑑S=\frac{1}{\rho _i}.$$
(4)
Equations (2) and (4) relate the level statistics in the semiclassical limit with the phase-space geometry.

In most general cases, the level spacing distribution might be singular. In such a case, it is convenient to use its cumulative distribution function $`\mu _i`$;

$$\mu _i(S)=_0^Sp_i(x;\rho _i)𝑑x.$$
(5)
In addition to equations (2) and (4), we assume two conditions for the statistical weights:

* Assumption (i): The statistical weights of independent regions uniformly vanishes in the limit of infinitely many regions;

$$\underset{i}{\mathrm{max}}\rho _i0\text{as}N+\mathrm{}.$$
(6)
* Assumption (ii): The weighted mean of the cumulative distribution of energy spacing,

$$\mu (\rho ;N)=\underset{i=1}{\overset{N}{}}\rho _i\mu _i(\rho ),$$
(7)
converges in $`N+\mathrm{}`$ to $`\overline{\mu }(\rho )`$

$$\underset{N+\mathrm{}}{lim}\mu (\rho ;N)=\overline{\mu }(\rho ).$$
(8)
The limit is uniform on each closed interval: $`0\rho S`$.

In the Berry-Robnik theory, the statistical weights of individual components are related to the phase volumes of the corresponding invariant regions. This relation is satisfactory if the Thomas-Fermi rule for the individual phase space regions still holds, and thus supports eq.(4). Here we do not specify the validity of this problem and deal with the statistical weights as parameters.

Under assumptions (i) and (ii), eqs.(2) and (4) lead to the overall level spacing distribution whose gap distribution function is given by the following formula in the limit of $`N+\mathrm{}`$,

$$E_{\overline{\mu }}(S)=\mathrm{exp}\left[_0^S\left(1\overline{\mu }(\sigma )\right)𝑑\sigma \right],$$
(9)
where the convergence is in the sense of the weak limit. When the level spacing distributions of individual components are sparse enough, one may expect $`\overline{\mu }=0`$ and the level spacing distribution of the whole energy sequence reduces to the Poisson distribution,

$$E_{\overline{\mu }=0}(S)=\mathrm{exp}\left(S\right).$$
(10)
In general, one may expect $`\overline{\mu }0`$ which corresponds to a certain accumulation of the levels of individual components.

In the following sections, the above statement is proved and the limiting level spacing distributions are classified into three classes. One of them is the Poisson distribution as discussed in the original work by Berry and Robnik. The others are not Poissonian. We give examples of non-Poissonian limiting level spacing distributions in section III. In the concluding section, we discuss some relations between our results and other related works.

## II Limiting level spacing Distribution

### A Derivation of the limiting Gap distribution

In this section, starting from eqs.(2) and (4), and assumptions (i) and (ii) introduced in the previous section, we show that, in the limit of infinitely many components $`N+\mathrm{}`$, the gap distribution $`E(S;N)`$ converges to the distribution function (9) with $`\overline{\mu }`$. The convergence is shown as follows.

Following the procedure by Mehta(see appendix A.2 of Ref.), we rewrite $`E(S;N)`$ in terms of the cumulative level spacing distribution $`\mu _i(S)`$ of independent components:

$$E(S;N)=\underset{i=1}{\overset{N}{}}\left[\rho _i_S^+\mathrm{}𝑑\sigma \left(1\mu _i(\sigma )\right)\right].$$
(11)
Equation (1.4) and integration by parts lead to

$$_0^+\mathrm{}𝑑\sigma (1\mu _i(\sigma ))=\sigma (1\mu _i(\sigma ))|_0^+\mathrm{}+_0^+\mathrm{}\sigma p_i(\sigma )𝑑\sigma =\frac{1}{\rho _i},$$

where $`lim_{\sigma +\mathrm{}}\sigma (1\mu _i(\sigma ))=0`$ follows from the existance of average, and hence to

$$\rho _i_S^+\mathrm{}𝑑\sigma (1\mu _i(\sigma ))=1\rho _i_0^S𝑑\sigma (1\mu _i(\sigma )).$$
(12)
Since the convergence of $`_{i=1}^N\rho _i\mu _i(\sigma )\overline{\mu }(\sigma )`$ for $`N+\mathrm{}`$ is uniform on each interval $`\sigma [0,S]`$ by Assumption (ii), $`E(S;N)`$ has the following limit,

$`\mathrm{log}E(S;N)`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{N}{}}}\mathrm{log}\left[1\rho _i{\displaystyle _0^S}𝑑\sigma (1\mu _i(\sigma ))\right]`$ (13)
$`=`$ $`{\displaystyle \underset{i=1}{\overset{N}{}}}\left[\rho _i{\displaystyle _0^S}𝑑\sigma (1\mu _i(\sigma ))+O(\rho _i^2)\right]`$ (14)
$`=`$ $`{\displaystyle _0^S}𝑑\sigma \left[1\mu (\sigma ;N)\right]+{\displaystyle \underset{i}{\overset{N}{}}}O(\rho _i^2)`$ (16)
$`{\displaystyle _0^S}𝑑\sigma \left[1\overline{\mu }(\sigma )\right]\text{as}N+\mathrm{},`$
where we have used $`|\mu _i(\sigma )|1`$, $`\mathrm{log}(1+ϵ)=ϵ+O(ϵ^2)`$ in $`ϵ1`$, and the following property obtained from Assumption (i),

$$|\underset{i=1}{\overset{N}{}}O(\rho _i^2)|C\underset{i}{\mathrm{max}}\rho _i\underset{i=1}{\overset{N}{}}\rho _i=C\underset{i}{\mathrm{max}}\rho _i0\text{as}N+\mathrm{},$$
(17)
with $`C`$ a positive constant. Therefore, we have the desired result:

$$\underset{N+\mathrm{}}{lim}E(S;N)=E_{\overline{\mu }}(S)=\mathrm{exp}\left[_0^S\left(1\overline{\mu }(\sigma )\right)𝑑\sigma \right].$$
(18)
### B Weak convergence limit of the level spacing distribution

In this section, we show that in $`N+\mathrm{}`$ limit, the level spacing distribution $`P(S;N)`$ converges weakly to $`P_{\overline{\mu }}(S)`$,

$$P_{\overline{\mu }}(S)=\frac{d^2}{dS^2}E_{\overline{\mu }}(S).$$
(19)
According to the Helly’s theorem, the weak convergence of the level spacing distribution is defined by

$$\underset{N+\mathrm{}}{lim}_0^SP(x;N)𝑑x=_0^SP_{\overline{\mu }}(x)𝑑x.$$
(20)
Since each side of the above equation is rewritten as,

$$_0^SP(x;N)𝑑x=\left[\frac{d}{dx}E(x;N)\right]_0^S,_0^SP_{\overline{\mu }}(x)𝑑x=\left[\frac{d}{dx}E_{\overline{\mu }}(x)\right]_0^S,$$

the limit (20) is equivalent to

$$\underset{N+\mathrm{}}{lim}\frac{d}{dS}E(S;N)=\frac{d}{dS}E_{\overline{\mu }}(S).$$

The above equation is proved as follows: By using equation (11), we rewrite $`\frac{d}{dS}E(S;N)`$ in terms of the cumulative level spacing distribution function $`\mu _i(S)`$ of spectral components,

$$\frac{d}{dS}E(S;N)=E(S;N)\underset{i=1}{\overset{N}{}}\frac{1\mu _i(S)}{\rho _i_S^+\mathrm{}𝑑\sigma (1\mu _i(\sigma ))}.$$

In the limit $`N+\mathrm{}`$, one has $`E(S;N)E_{\overline{\mu }}(S)`$ as shown by equation (18) and,

$`{\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle \frac{1\mu _i(S)}{\rho _i_S^+\mathrm{}𝑑\sigma (1\mu _i(\sigma ))}}`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle \frac{\rho _i\rho _i\mu _i(S)}{1\rho _i_0^S𝑑\sigma (1\mu _i(\sigma ))}}`$ (21)
$`=`$ $`1{\displaystyle \underset{i=1}{\overset{N}{}}}\rho _i\mu _i(S)+{\displaystyle \underset{i=1}{\overset{N}{}}}O(\rho _i^2)`$ (23)
$`1\overline{\mu }(S)\text{as}N+\mathrm{},`$
where we have used equation (12), $`1/(1ϵ)=1+O(ϵ)`$ in $`ϵ1`$, and the limit (17). Therefore, we have the desired result:

$$\underset{N+\mathrm{}}{lim}\frac{d}{dS}E(S;N)=\left(1\overline{\mu }(S)\right)E_{\overline{\mu }}(S)=\frac{d}{dS}E_{\overline{\mu }}(S).$$
(24)
We remark that, when the limiting distribution function $`\overline{\mu }(S)`$ is differentiable, the asymptotic level spacing distribution is described as follows:

$$P_{\overline{\mu }}(S)=\left[(1\overline{\mu }(S))^2+\overline{\mu }^{}(S)\right]\mathrm{exp}\left[_0^S\left(1\overline{\mu }(\sigma )\right)𝑑\sigma \right].$$
(25)
### C Properties of the limiting level spacing distribution

Since $`\mu _i(S)`$ is monotonically increasing and $`0\mu _i(S)1`$, $`\overline{\mu }(S)`$ has the same properties. Then, $`1\overline{\mu }(S)0`$ for any $`S0`$ and one has

$$\frac{1}{S}_0^S𝑑\sigma (1\overline{\mu }(\sigma ))1\overline{\mu }(+\mathrm{})\text{as}S+\mathrm{}.$$
(26)
According to the above limit, the level spacing distribution is classified into the following three cases in the sense of weak limit:

* Case 1, $`\overline{\mu }(+\mathrm{})=0`$: The limiting level spacing distribution is the Poisson distribution. Note that this condition is equivalent to $`\overline{\mu }(S)=0`$ for $`{}_{}{}^{}S`$ because $`\overline{\mu }(S)`$ is monotonically increasing.
* Case 2, $`0<\overline{\mu }(+\mathrm{})<1`$: For large $`S`$ values, the limiting level spacing distribution is well approximated by the Poisson distribution, while, for small $`S`$ values, it may deviate from the Poisson distribution.
* Case 3, $`\overline{\mu }(+\mathrm{})=1`$: The limiting level spacing distribution deviates from the Poisson distribution for $`{}_{}{}^{}S`$, and decays as $`S+\mathrm{}`$ more slowly than does the Poisson distribution. This case will be referred to as a sub-Poisson distribution.

One has Case 1 if the individual cumulative distribution function $`\mu _i(S)`$ are bounded by a finite positive function $`g(S)`$,

$$\mu _i(S)\rho _i^\eta g(S),$$

for $`{}_{}{}^{}i`$ and $`\eta >0`$. Indeed, one has,

$$\mu (S;N)=\underset{i=1}{\overset{N}{}}\rho _i\mu _i(S)g(S)\underset{i}{\overset{N}{}}\rho _i\rho _i^\eta g(S)(\underset{i}{\mathrm{max}}\rho _i)^\eta 0\overline{\mu }(S),\text{as}N\mathrm{}.$$

More specifically, for example, one has Case 1 if the individual level spacing distributions are derived from scaled distribution functions $`f_i`$ as***Eq.(4) is described by rewriting $`x=\rho _iS`$ in the following way, $`1=S\rho _ip_i(S;\rho _i)𝑑S=xf_i(x)𝑑x`$, so that $`f_i(x)=\rho _ip_i(x;\rho _i)`$ is a dimensionless function. For instance, $`f_i(x)=\mathrm{exp}(x)`$ in the case of Poisson distribution $`p_i(S;\rho _i)=\rho _i\mathrm{exp}\left(\rho _iS\right)`$, and $`f(x)=\frac{\pi }{2}x\mathrm{exp}\left(\frac{\pi }{4}x^2\right)`$ in the case of the Wigner distribution $`p_i(S;\rho _i)=\frac{\pi }{2}\rho _i^2S\mathrm{exp}\left[\frac{\pi }{4}\rho _i^2S^2\right]`$.,

$$p_i(S;\rho _i)=\rho _if_i(\rho _iS),$$

where $`f_i`$ satisfy

$$_0^+\mathrm{}f_i(x)𝑑x=1,_0^+\mathrm{}xf_i(x)𝑑x=1,$$

and are uniformly bounded by a positive constant $`D`$: $`|f_i(S)|D`$ ($`1iN`$ and $`S0`$). Indeed, one then has

$`|\mu (S;N)|`$ $``$ $`\left|{\displaystyle \underset{i=1}{\overset{N}{}}}\rho _i{\displaystyle _0^S}p_i(x,\rho _i)𝑑x\right|`$ (27)
$``$ $`{\displaystyle \underset{i=1}{\overset{N}{}}}\rho _i^2{\displaystyle _0^S}\left|f_i\left(\rho _ix\right)\right|𝑑x`$ (28)
$``$ $`DS{\displaystyle \underset{i=1}{\overset{N}{}}}\rho _i^2DS\underset{i}{\mathrm{max}}\rho _i{\displaystyle \underset{i=1}{\overset{N}{}}}\rho _i0\overline{\mu }(S).`$ (29)
This includes the case studied by Berry and Robnik, where the gap distribution is a product of superposition of a single regular component characterized by the Poisson distribution and $`N`$ equivalent chaotic components characterized by the Wigner distribution, and the latter is expressed by the product of the scaled distributions as:

$$E^{BR}(S;N)=\mathrm{exp}\left(\rho _0S\right)\underset{i=1}{\overset{N}{}}E_i^{\text{WIGNER}}(S;\rho _i),$$
(30)
where the statistical weights are $`\rho _i=\frac{1\rho _0}{N}`$ and the individual level spacing distributions $`f_i`$ corresponding to the gap distributions $`E_i^{\text{WIGNER}}(S;\rho _i)=\text{erfc}\left(\frac{\sqrt{\pi }}{2}\rho _iS\right)`$ are given by the dimensionless Wigner distribution:

$$f_i(x)=\frac{\pi x}{2}\mathrm{exp}\left[\frac{\pi }{4}x^2\right].$$
(31)
Indeed, one has the Poisson distribution in $`N+\mathrm{}`$ limit:

$`E^{\text{BR}}(S;N)`$ $`=`$ $`\mathrm{exp}\left[\rho _0S+(N1)\mathrm{log}\text{erfc}\left({\displaystyle \frac{\sqrt{\pi }}{2}}{\displaystyle \frac{1\rho _0}{N1}}S\right)\right]e^S.`$ (32)
## III Example

As an example of the deviation from the Poisson distribution, we study rectangular billiard system whose energy levels are described by using positive integer numbers, $`m`$ and $`i`$ as follows,

$$ϵ_{m,i}=m^2+\alpha i^2,$$
(33)
where $`\alpha `$ is the aspect ratio of two sides of a billiard wall. For a given energy interval $`[ϵ,ϵ+\mathrm{\Delta }ϵ]`$, each energy level is classified into components according to the eigenvalues, $`i=1,2,3,\mathrm{},N`$;

$$N=\left[\sqrt{\frac{1+\gamma }{\alpha }ϵ\text{ }}\right],$$
(34)
where $`\gamma \mathrm{\Delta }ϵ/ϵ`$, and $`[x]`$ stands for the maximum integer which does not exceed $`x`$. The relative weight of each component is given by

$$\rho _i=\{\begin{array}{cc}\frac{4(1+\gamma )}{N\pi \gamma }\left(\sqrt{1\left(\frac{i}{N}\right)^2}\sqrt{\frac{1}{1+\gamma }\left(\frac{i}{N}\right)^2}\right)+O\left(\frac{1}{N^2}\right)\hfill & \hfill \text{if}i<\frac{N}{\sqrt{1+\gamma }},\\ \frac{4(1+\gamma )}{N\pi \gamma }\left(\sqrt{1\left(\frac{i}{N}\right)^2}\right)+O\left(\frac{1}{N^2}\right)\hfill & \hfill \text{if}\frac{N}{\sqrt{1+\gamma }}iN.\end{array}$$
(35)
As easily seen, $`\rho _i`$ satisfies the assumption (i);

$$\underset{i}{\mathrm{max}}\rho _i\frac{4}{N\pi }\sqrt{1+\frac{1}{\gamma }}+O\left(\frac{1}{N^2}\right)0\text{as}N+\mathrm{}.$$
(36)
Note that the limit of infinitely many components, $`N+\mathrm{}`$, corresponds to the high energy limit, $`ϵ+\mathrm{}`$(see Eq.(34)), which is equivalent to the semiclassical limit in physical systems. In this limit, the statistical weight of each sub-spectrium becomes sparse, since each element of $`\mu (S;N)`$, $`\rho _i\mu _i(S)`$, tends to zero: $`\rho _i\mu _i(S)\mathrm{max}_j\rho _j0`$.

Figures 1(a) and 1(b) show numerical results of the level-spacing distribution $`P(S)`$ for two values of $`\alpha `$, and figures 2(a) and 2(b) show the gap distribution function corresponding to figures 1(a) and 1(b), respectively. In case that $`\alpha `$ is far from rational, $`P(S)`$ and $`E(S)`$ are well approximated by the Poisson distribution while in case that $`\alpha `$ is close to a rational expressed as $`\alpha =p/q`$, where $`p`$ and $`q`$ are coprime positive integers, they deviate from the Poisson distribution.

In order to compare the non-Poisson distribution and the classification given in the previous section, we consider,

$$\stackrel{~}{\mu }(S;N)=1\frac{1_0^SP(x;N)𝑑x}{E(S;N)}.$$
(37)
When $`N+\mathrm{}`$, $`\stackrel{~}{\mu }(S;N)`$ approaches $`\overline{\mu }(S)`$, and this function distinguishes the three cases as follows: In Case 1 i.e., where the level spacing obeys the Poisson distribution, $`lim_{N+\mathrm{}}\stackrel{~}{\mu }(S;N)=0`$. In Case 2, $`lim_{N+\mathrm{}}\stackrel{~}{\mu }(S;N)c`$ as $`S+\mathrm{}`$ $`(0<c<1)`$, and in Case 3, where the sub-Poisson distribution is expected, $`lim_{N+\mathrm{}}\stackrel{~}{\mu }(S;N)1`$ as $`S+\mathrm{}`$.

Figure 3 shows $`\stackrel{~}{\mu }(S;N)`$ for different values of $`N`$. The dotted line $`\stackrel{~}{\mu }=0`$ exhibits the Poisson distribution. From this, one can think that $`\stackrel{~}{\mu }(S;N)`$ for $`N=61905,S10`$ well approximates $`lim_{N+\mathrm{}}\stackrel{~}{\mu }(S;N)`$.

Figure 4 shows $`\stackrel{~}{\mu }(S;N)`$ for the two values of $`\alpha `$ corresponding to figures 2(a)–2(b), respectively. In case that the numerical data is well characterized by the Poisson distribution (figures 1(a) and 2(a)), the corresponding function $`\stackrel{~}{\mu }(S;N)`$ agrees with $`0`$, while in case that deviates from the Poisson distribution (figures 1(b) and 2(b)), $`\stackrel{~}{\mu }(S;N)0`$ and $`\stackrel{~}{\mu }(S;N)c(0<c<1)`$ for $`S+\mathrm{}`$. Therefore, this result corresponds to the Case 2.

In this model, we have not yet observed the clear evidence of Case 3. Such a case is expected when there is stronger accumulation of the energy levels of individual components.

## IV Extended formalism of the Berry-Robnik distribution

In this section, we propose one possible extension of the Berry-Robnik distribution (30) for the level statistics of the nearly-integrable system with two degree-of-freedom. Since the classical phase space of the nearly-integrable system consists of regular and chaotic regions and the Liouville measures of the chaotic regions are larger than zero, $`\rho _i>0`$, this system does not support the assumption (i). However, the regular regions consist of infinitely many subsets, and our approach shown in section II is partially applicable to the spectral components corresponding to the regular regions.

Following to the assumption proposed by Berry and Robnik, the gap distribution functions in the nearly integrable system, which are contributed from the individual chaotic regions, are characterized by the Random Matrix Theory(RMT). Then one has,

$`\mathrm{log}E(S;N_1;N_2)`$ $`=`$ $`\mathrm{log}{\displaystyle \underset{i=1}{\overset{N_1}{}}}E_i^{\text{RMT}}(S;\rho _i){\displaystyle \underset{j=1}{\overset{N_2}{}}}E_j(S;\rho _j^{^{}})`$ (38)
$`=`$ $`{\displaystyle \underset{i=1}{\overset{N_1}{}}}\mathrm{log}E_i^{\text{RMT}}(S;\rho _i)+{\displaystyle \underset{j=1}{\overset{N_2}{}}}\mathrm{log}E_j(S;\rho _j^{^{}}),`$ (39)
where $`\rho _j^{^{}}=\rho _{j+N_1}`$, and $`E_j(S;\rho _j^{^{}})`$ denote the gap distribution functions corresponding to the subsets in the regular regions. As shown by eq.(2.6), $`E_j`$ has the following limit,

$$\underset{j=1}{\overset{N_2}{}}\mathrm{log}E_j(S;\rho _j^{^{}})\rho _0_0^S𝑑\sigma \left[1\overline{\mu }(\sigma )\right]\text{as}N_2+\mathrm{},$$

where $`\rho _0=_{j=1}^{N_2}\rho _j^{^{}}`$ and $`\overline{\mu }(\sigma )=lim_{N_2+\mathrm{}}_{j=1}^{N_2}\rho _j^{^{}}\mu _j(\sigma )`$. Accordingly, in the partial limit of $`N_1+\mathrm{},N_2+\mathrm{}`$, the original proposal for the gap distribution by Berry and Robnik is replaced by

$$\underset{N_2+\mathrm{}}{lim}E(S;N_1;N_2)=E_{\overline{\mu }}(S;N_1)=\mathrm{exp}\left[\rho _0_0^S(1\overline{\mu }(\sigma ))𝑑\sigma \right]\underset{i=1}{\overset{N_1}{}}E_i^{\text{RMT}}(S;\rho _i),$$
(40)
where $`\rho _0`$ denotes the total amount of the Liouville measures of the regular region in mixed phase space. The above distribution formula is classified into the following three cases; Case 1’, $`\overline{\mu }(+\mathrm{})=0`$: Berry-Robnik distribution, Case 2’, $`0<\overline{\mu }(+\mathrm{})<1`$: Berry-Robnik distribution for large $`S`$, but possibly not for small $`S`$, and Case 3’, $`\overline{\mu }(+\mathrm{})=1`$: A distribution function obtained by the superposition of spectral components obeying the sub-Poisson statistics and the Random matrix theory. From this classification, one can see that the new formula (40) admits deviations from the Berry-Robnik distribution when $`\overline{\mu }(+\mathrm{})0`$.

We remark that $`P(S;N_1,N_2)=\frac{d^2}{dS^2}E(S;N_1,N_2)`$ in the limit $`N_2+\mathrm{}`$ converges weakly to the limiting level spacing distribution: $`P_{\overline{\mu }}(S;N_1)=\frac{d^2}{dS^2}E_{\overline{\mu }}(S;N_1)`$, and when the limiting function $`\overline{\mu }(S)`$ is differentiable, the asymptotic level spacing distribution admits the following density:

$$P_{\overline{\mu }}(S;N_1)=\frac{d^2}{dS^2}\left[\mathrm{exp}\left(\rho _0_0^S(1\overline{\mu }(\sigma ))𝑑\sigma \right)\underset{i=1}{\overset{N_1}{}}E_i^{\text{RMT}}(S;\rho _i)\right]$$
(41)
The validity of the Berry-Robnik distribution has been checked for generic nearly-integrable systems by many numerical investigations. Among them, Prosen and Robnik found numerically for several systems that there is a high energy region in which the Berry-Robnik distribution formula (30) well approximates the level spacing distribution. This energy region is sometimes referred to as the Berry-Robnik regime. While they also found that the level spacing distribution in the low energy region deviates from the Berry-Robnik formula, and approximates the Brody distribution. This behavior was studied in terms of a fractional power dependence of the spacing distribution near the origin at $`S=0`$, which could be attributed to the localization properties of eigenstates on chaotic components. From the above classification, Case 1’:$`\overline{\mu }(+\mathrm{})=0`$ sould be satisfied in the Berry-Robnik regime. While Case 2’ and Case 3’ might propose another possibilities. When the spectral components corresponding to regular regions show strong accumulation, the gap distribution obeys the distribution formula (40) with $`0<\overline{\mu }(+\mathrm{})1`$, and the level spacing distribution shows deviations from the Berry-Robnik distribution. Therefore, this result might propose another possibility of the Berry-Robnik approach.

## V Conclusion and Discussion

In this paper, we investigated the gap distribution function of systems with infinitely many independent components and discussed the level-spacing statistics of classically integrable quantum systems. In the semiclassical limit, reflecting infinitely fine classical phase space structures, individual energy eigenfunctions are expected to be well localized in the phase space and contribute independently to the level statistics. Keeping this expectation in mind, we considered a situation in which the system consists of infinitely many components and each of them gives an infinitesimal contribution. And by applying the arguments of Mehta, and Berry and Robnik, the limiting level spacing distribution was obtained whose gap distribution function is described by a single monotonically increasing function $`\overline{\mu }(S)`$ of the level spacing $`S`$:

$$E_{\overline{\mu }}(S)=\mathrm{exp}\left[_0^S\left(1\overline{\mu }(\sigma )\right)𝑑\sigma \right]$$
(42)
The weak convergence limit of the level spacing distribution is classified into three cases; Case 1: Poissonian if $`\overline{\mu }(+\mathrm{})=0`$, Case 2: Poissonian for large $`S`$, but possibly not for small $`S`$ if $`0<\overline{\mu }(+\mathrm{})<1`$, and Case 3: sub-Poissonian if $`\overline{\mu }(+\mathrm{})=1`$. Thus, even when the energy levels of individual components are statistically independent, non-Poissonian level spacing distributions are possible.

In most general cases, the integral in equation (42) converges in $`S+\mathrm{}`$ and then $`lim_{S+\mathrm{}}E_{\overline{\mu }}(S)0`$, the limiting gap distribution $`E_{\overline{\mu }}(S)`$ does not work accurately. In such case, however, its differentiation (24) still work accurately in $`S+\mathrm{}`$ limit, and thus the above classification (Case 1–3) holds in general.

Note that the singular level spacing distribution can be taken into account in terms of non-smooth cumulative distributions. Such a singularity is expected when there is strong accumulation of the energy levels of individual components. For a certain class of systems, such accumulation is observable. One example is shown in section III where the results show clear evidence of Case 2. Another example is the two-dimensional harmonic oscillator whose level spacing distribution is non-smooth for arbitrary system parameter. The final example is studied by Shnirelman, Chirikov and Shepelyansky, and Frahm and Shepelyansky for a certain type of system which contains a quasi-degeneracy result from inherent symmetry(time reversibility). As is well known, the existence of quasi-degeneracy leads to the sharp Shnirelman peak at small spacings.

Finally, in section IV, we proposed one possible extension of the Berry-Robnik distribution for classically nearly-integrable quantum systems. This extension admitted deviations from the Berry-Robnik distribution when there is strong accumulation of the energy levels of spectral components. Such possibilities will be studied elsewhere.
Acknowledgments
The authors would like to thank Professor M. Robnik, Professor Y.Aizawa, Professor A.Shudo, and Dr G. Veble for their helpful advices and discussions. The authors also thank the Yukawa Institute for Theoretical Physics at Kyoto University. Discussions during the YITP workshop YITP-W-02-13 on ”Quantum chaos: Present status of theory and experiment” were useful to complete this work. This work is partly supported by a Grant-in-Aid for Scientific Research (C) from the Japan Society for the Promotion of Science.