text:
# UWThPh-2001-6 4-Neutrino mass schemes and the likelihood of (3+1)-mass spectra

## 1 Introduction

At present, there are three indications in favour of neutrino oscillations , namely the solar $`\nu _e`$ deficit , the atmospheric $`\underset{\mu }{\overset{()}{\nu }}`$  deficit and the result of the LSND experiment  hinting at $`\underset{\mu }{\overset{()}{\nu }}\underset{e}{\overset{()}{\nu }}`$ transitions. Whereas in the case of the first two indications several experiments agree on the existence of the effect, the third indication is found only by the LSND collaboration. Therefore, in many analyses the LSND result is left out. However, if all three indications in favour of neutrino oscillations are confirmed, for three mass-squared differences of different orders of magnitude ($`10^{10}\text{eV}^2<\mathrm{\Delta }m_{\mathrm{solar}}^2<10^7\text{eV}^2`$ or $`\mathrm{\Delta }m_{\mathrm{solar}}^210^5\text{eV}^2`$, $`\mathrm{\Delta }m_{\mathrm{atm}}^23\times 10^3`$ eV2, $`\mathrm{\Delta }m_{\mathrm{LSND}}^21`$ eV2) one needs a minimum of four neutrinos, three active ones and a sterile one . In that case a major revision of our picture of the lepton sector of the elementary particles would be necessary, with a mixing between the active and the sterile neutrinos; i.e.,

$$\nu _{\alpha L}=\underset{j=1}{\overset{4}{}}U_{\alpha j}\nu _{jL}\text{with}\alpha =e,\mu ,\tau ,s,$$
(1)
if we stick to the minimum of four neutrinos. In Eq. (1) the left-handed flavour fields are denoted by $`\nu _{\alpha L}`$ and the left-handed mass eigenfields by $`\nu _{jL}`$, and the $`4\times 4`$ neutrino mixing matrix $`U`$ is assumed to be unitary.

One of the most important issues in the context of 4-neutrino scenarios is the question of the 4-neutrino mass spectrum . There are two different spectral classes with very different properties: the first class contains four types and consists of spectra where three neutrino masses are clustered together, whereas the fourth mass is separated from the cluster by the mass gap needed to reproduce the LSND result;333This class contains the hierarchical mass spectrum. the second class has two types where two pairs of nearly degenerate masses are separated by the LSND gap. These two classes have been dubbed (3+1) and (2+2)-neutrino mass spectra, respectively . The main difference between these two classes is that, if a (2+2)-spectrum is realized in nature, the transition into the sterile neutrino is a solution of either the solar or the atmospheric neutrino problem, or the sterile neutrino has to take part in both, whereas with a (3+1)-spectrum it could be only slightly mixed with the active neutrinos and mainly provide a description of the LSND result.

It has been argued in the literature  that the (3+1)-spectra are strongly disfavoured by the data, whereas the (2+2)-spectra are the preferred ones, in agreement with all data showing evidence for neutrino oscillations and also with those where no such evidence has been found. Recently, in Refs.  this statement has been challenged because in the latest LSND analysis the allowed region in the $`\mathrm{\Delta }m^2`$$`A_{\mu ;e}`$ plane, where $`A_{\mu ;e}`$ is the LSND transition amplitude, has undergone a slight shift towards smaller mass-squared differences, which makes the (3+1)-spectra somewhat less disfavoured. Furthermore, in a 2-neutrino analysis of atmospheric neutrino oscillations the Super-Kamiokande data prefer $`\nu _\mu \nu _\tau `$ conversion over $`\nu _\mu \nu _s`$ . Moreover, there is some debate also in the solar neutrino problem whether the $`\nu _e\nu _s`$ transition is disfavoured in comparison with other solutions , though such a feature seems not to be borne out by a global analysis of the data . In any case, moving away from pure 2-neutrino considerations in the solar and atmospheric neutrino problems, transitions into active–sterile superpositions  give viable solutions to both problems within the (2+2)-spectral schemes, with features which will be tested in the future .

The arguments presented in Refs. , which disfavour the (3+1)-mass spectra, are based on exclusion curves from short-baseline (SBL) experiments, and solar and atmospheric neutrino data enter into this simplified analysis only through inequalities. The advantage of this approach is that its parameters are confined to the quantities444Note that in Ref.  the quantities $`c_\alpha =1d_\alpha `$ are used instead.

$$d_\alpha =|U_{\alpha 4}|^2(\alpha =e,\mu )$$
(2)
and the SBL or LSND mass-squared difference $`\mathrm{\Delta }m^2`$. For definiteness we assume that the mass separated by the LSND gap is $`m_4`$ and, therefore, $`\mathrm{\Delta }m^2=|m_4^2m_1^2|`$. It has turned out that the up-down asymmetry of atmospheric multi-GeV $`\mu `$-like events measured in the Super-Kamiokande experiment  is very suitable to constrain $`d_\mu `$ , whereas from the solar data it follows that $`d_e`$ must be small. The probability of SBL $`\nu _\mu \nu _e`$ transitions is given by the two-neutrino-like formula

$$P_{\nu _\mu \nu _e}=P_{\overline{\nu }_\mu \overline{\nu }_e}=A_{\mu ;e}\mathrm{sin}^2\frac{\mathrm{\Delta }m^2L}{4E},$$
(3)
where $`L`$ is the distance between source and detector and $`E`$ is the neutrino energy and

$$A_{\mu ;e}=4d_ed_\mu .$$
(4)
The LSND experiment gives an allowed region in the $`\mathrm{\Delta }m^2`$$`A_{\mu ;e}`$ plane.

However, the arguments of Refs.  are not based on a well-defined statistical procedure. Therefore, they remain on a semi-quantitative level and do not allow to assess a confidence level (CL) which quantifies the degree at which the (3+1)-spectra are excluded.

In this paper we make a step forward towards such an assessment. The main points to achieve our goal are the following:

* The aim is to arrive, for every SBL mass-squared difference $`\mathrm{\Delta }m^2`$, at a probability distribution solely in terms of $`d_e`$ and $`d_\mu `$; a suitable method for this purpose is given by the likelihood function in combination with the Bayesian approach (see, e.g., Refs. ).
* We make full use of the data of the SBL Bugey , CDHS  and KARMEN  experiments through the likelihood function.
* In the spirit of Refs. , all information pertaining to the atmospheric and solar mass-squared differences is included via inequalities. Within our probabilistic framework we are able to treat inequalities as prior probabilities or with a kind of a maximum likelihood method.
* In this way we treat the inequality following from the atmospheric up-down asymmetry and, similarly, we include also the result of the CHOOZ experiment ; the solar neutrino data allow a simpler treatment in the context of the Bugey data , as described in Refs. .
* Eventually, for every given SBL $`\mathrm{\Delta }m^2`$ and any CL $`\beta `$, we are able to calculate an upper bound $`A_\beta ^0(\mathrm{\Delta }m^2)`$ on the transition amplitude $`A_{\mu ;e}`$, and we can compare such bounds $`A_\beta ^0(\mathrm{\Delta }m^2)`$ with the 90% and 99% CL regions in the $`\mathrm{\Delta }m^2`$$`A_{\mu ;e}`$ plane found by the LSND experiment .

In the same framework we will also discuss the (2+2)-spectra.

The plan of the paper is as follows. In Section 2 we rederive the atmospheric up-down inequality  in a form which is suitable for our purpose. In Section 3 we introduce the likelihood function and the Bayesian approach, and describe how to incorporate inequalities; we apply the methods developed there to the CHOOZ result and the atmospheric up-down inequality. In Section 4 we explain how we calculate bounds on $`A_{\mu ;e}`$ and discuss each of the SBL experiments we use, together with their features which are important in this context. Details of technical nature are deferred to the appendix. Our main result, represented as a plot in the $`\mathrm{\Delta }m^2`$$`A_{\mu ;e}`$ plane, is also given in this section. In Section 5 we consider the (2+2)-spectra in the framework of our statistical approach and in Section 6 we draw our conclusions.

## 2 The atmospheric up-down asymmetry as a constraint on short-baseline neutrino oscillations

The most convincing evidence for $`\underset{\mu }{\overset{()}{\nu }}`$ disappearance in atmospheric neutrino experiments is given by the so-called up-down asymmetry

$$A_{\mathrm{ud}}=\frac{UD}{U+D}$$
(5)
measured by the Super-Kamiokande Collaboration , where $`U`$ and $`D`$ refer to the number of up-going and down-going $`\mu `$-like multi-GeV events, respectively. Quoting the number for fully contained events, after 1289 days of operation the result

$$A_{\mathrm{ud}}^{\mathrm{exp}}=0.327\pm 0.045\pm 0.004$$
(6)
was found . Adding statistical and systematic error in quadrature, one obtains

$$\mathrm{\Delta }A_{\mathrm{ud}}^{\mathrm{exp}}\sigma _A=0.045.$$
(7)
Because of the smallness of the systematic error this value is identical with the statistical error.

Let us now rederive the atmospheric up-down inequality . In the following we will not indicate antineutrinos but our arguments will hold for both, neutrinos and antineutrinos. With the assumption that downward-going atmospheric neutrinos do not oscillate with the frequency associated with $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ and that oscillations according to $`\mathrm{\Delta }m^2`$ are averaged out, we obtain

$$P_{\nu _\mu \nu _\mu }^D=d_\mu ^2+(1d_\mu )^2\text{and}P_{\nu _e\nu _\mu }^D=\frac{1}{2}A_{\mu ;e}.$$
(8)
Denoting the number of muon (electron) neutrinos and antineutrinos produced in the atmosphere by $`n_\mu `$ ($`n_e`$), it follows from Eq. (8) that

$$D=n_\mu [d_\mu ^2+(1d_\mu )^2]+\frac{1}{2}n_eA_{\mu ;e}.$$
(9)
For the upward-going neutrinos we have the inequalities

$$P_{\nu _\mu \nu _\mu }^Ud_\mu ^2\text{and}P_{\nu _e\nu _\mu }^U\frac{1}{4}A_{\mu ;e}.$$
(10)
We, therefore, have the inequality

$$Un_\mu d_\mu ^2+\frac{1}{4}n_eA_{\mu ;e}.$$
(11)
Note that for the (3+1)-spectra the amplitude $`A_{\mu ;e}`$ is given by Eq. (4). Since $`A_{\mathrm{ud}}`$ (5) is a monotonously decreasing function in $`U`$, using Eqs. (9) and (11) we obtain the so-called up-down inequality for $`\mu `$-like atmospheric events

$$A_{\mathrm{ud}}G(d_e,d_\mu )=\frac{(1d_\mu )^2+d_ed_\mu /r}{(1d_\mu )^2+2d_\mu ^2+3d_ed_\mu /r}.$$
(12)
In this equation we have defined $`r=n_\mu /n_e`$. The numerical value $`r2.8`$ can be read off from Fig. 3 in Ref.  of the Super-Kamiokande Collaboration. With similar arguments one can also find an upper bound

$$A_{\mathrm{ud}}H(d_e,d_\mu )=\frac{d_\mu (1d_\mu )d_ed_\mu /r}{1d_\mu (1d_\mu )+d_ed_\mu /r}.$$
(13)
Hence, the up-down asymmetry is confined to the interval

$$GA_{\mathrm{ud}}H.$$
(14)
Some remarks are at order. Our inequality (12) is a little different from the analogous inequality in Ref. , because in the present case we have not eliminated $`d_e`$; this is useful because our aim is to derive a probability distribution in $`d_e`$ and $`d_\mu `$ (see next section). The ratio $`r`$ has a slight dependence on the atmospheric zenith angle which is neglected here; however, since $`d_e`$ is confined to rather small values  by the result of the Bugey experiment , the terms containing $`r`$ are rather unimportant numerically. The assumptions for deriving Eqs. (12) and (13) are not exactly fulfilled: for down-going neutrinos with zenith angles around $`0^{}`$, oscillations according to small SBL mass-squared differences $`\mathrm{\Delta }m^2`$ around 0.2 eV2 are not completely averaged out; for down-going neutrinos with zenith angles close to $`90^{}`$, oscillations according to $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ do occur already, if the atmospheric mass-squared difference is large. We have checked that both effects do not change numerically the bound $`G`$ by more than a few percent in the worst case. Finally, matter effects, which are important for up-going neutrinos, do not affect the inequalities (10) and thus also not the inequality for $`U`$. Eq. (12) is, therefore, a well-established inequality which restricts mainly the allowed range of $`d_\mu `$. Due to Eq. (4), it will be used in the following to constrain the SBL amplitude $`A_{\mu ;e}`$.

## 3 The statistical treatment of inequalities

As mentioned in the introduction, a suitable method for the purpose of deriving a probability distribution in the variables $`d_e`$ and $`d_\mu `$ is given by the likelihood function combined with the Bayesian approach, which is defined as follows. Suppose one has a series of measurements $`x=(x_1,\mathrm{},x_n)`$ and $`r`$ parameters $`\theta =(\theta _1,\mathrm{},\theta _r)`$ to be estimated. Then the Bayesian approach allows to construct a “posterior” probability density in the parameter space via (see, e.g., Ref. )

$$p(\theta |x)=\frac{L(x|\theta )\pi (\theta )}{d^r\theta ^{}L(x|\theta ^{})\pi (\theta ^{})}.$$
(15)
In this expression, $`L(x|\theta )`$ is the likelihood function and $`\pi (\theta )`$ is the prior probability density associated with the parameters $`\theta `$, reflecting the state of knowledge of $`\theta `$ before the measurement.

### 3.1 Inequalities included as priors

Let us suppose now that we consider an observable $`Z`$ whose experimental value is $`z_{\mathrm{exp}}\pm \sigma _z`$ and its *true* value is $`z`$. We assume that the values of measurements of $`Z`$ are distributed according to a Gaussian distribution

$$L_Z=\frac{1}{\sqrt{2\pi }\sigma _z}\mathrm{exp}\left[\frac{1}{2}\left(\frac{z_{\mathrm{exp}}z}{\sigma _z}\right)^2\right]$$
(16)
around the true value. Suppose further that from some theoretical consideration we have the knowledge that the true value $`z`$ is bounded by

$$azb.$$
(17)
We conceive the true value $`z`$, which is otherwise unknown, as a parameter in our scenario and assign to it a prior probability density

$$\pi _Z(z)=\frac{1}{ba}\mathrm{\Theta }(bz)\mathrm{\Theta }(za),$$
(18)
where $`\mathrm{\Theta }`$ denotes the Heaviside function. Due to lack of further knowledge, we have assumed a flat prior probability density. Since we are not interested in a posterior probability density in the parameter $`z`$ we perform the integral

$$𝑑zL_Z(z)\pi _Z(z)=\mathrm{}_Z$$
(19)
with

$$\mathrm{}_Z=\frac{1}{2(ba)}\left\{\mathrm{erf}\left(\frac{bz_{\mathrm{exp}}}{\sqrt{2}\sigma _z}\right)+\mathrm{erf}\left(\frac{z_{\mathrm{exp}}a}{\sqrt{2}\sigma _z}\right)\right\},$$
(20)
where we have made use of the error function defined by

$$\mathrm{erf}(z)=\frac{2}{\sqrt{\pi }}_0^z𝑑te^{t^2}.$$
(21)
The function (20) represents then the relevant factor in the posterior probability density which takes into account the inequality (17).

To check if we have been lead to a meaningful expression (20), we want to discuss the behaviour of this function. Following from $`\mathrm{erf}(\mathrm{})=\mathrm{erf}(\mathrm{})=1`$, we find that for $`az_{\mathrm{exp}}b`$ we have $`\mathrm{}_Z1/(ba)`$ to a very good approximation, and for $`z_{\mathrm{exp}}b`$ or $`z_{\mathrm{exp}}a`$ we have $`\mathrm{}_A0`$. In this discussion, “much smaller” or “much larger” is defined in units of $`\sigma _z`$. If we assume that $`\sigma _z`$ becomes negligibly small, the function $`\mathrm{}_Z`$ approaches a step function with $`1/(ba)`$ being the height of the step. The edges of the step are smoothened out by a finite $`\sigma _z`$. Thus, we will have a maximal contribution to the posterior probability density for $`z_{\mathrm{exp}}`$ well inside the interval $`[a,b]`$, whereas for $`z_{\mathrm{exp}}`$ well outside the interval $`[a,b]`$ the function $`\mathrm{}_Z`$ is very close to zero. Thus we can be confident that the inequality (17) is reasonably well taken into account by our procedure.

### 3.2 Inequalities treated with a maximum likelihood method

Again we depart from Eqs. (16) and (17). Now we take into account our lack of knowledge about the true value of $`Z`$ by maximizing $`L_Z`$ as a function of $`z[a,b]`$. It is easy to check that one obtains

$`\underset{z[a,b]}{\mathrm{max}}L_ZL_Z^m=`$ (22)
$`{\displaystyle \frac{1}{\sqrt{2\pi }\sigma _z}}\mathrm{exp}\left\{{\displaystyle \frac{1}{2}}\left[\left({\displaystyle \frac{z_{\mathrm{exp}}b}{\sigma _z}}\right)^2\mathrm{\Theta }(z_{\mathrm{exp}}b)+\left({\displaystyle \frac{az_{\mathrm{exp}}}{\sigma _z}}\right)^2\mathrm{\Theta }(az_{\mathrm{exp}})\right]\right\}.`$
In Appendix A we will discuss how this method is related to the method in the previous section.

### 3.3 The treatment of the CHOOZ result and the atmospheric up-down inequality in our statistical approach

The CHOOZ experiment is a long-baseline $`\overline{\nu }_e`$ disappearance experiment. It measures the survival probability $`P_{\mathrm{CH}}`$, for which we derive the inequality

$$P_{\mathrm{CH}}=12d_e(1d_e)2|U_{e3}|^2(1d_e|U_{e3}|^2)(1\mathrm{cos}\varphi _{\mathrm{atm}})12d_e(1d_e).$$
(23)
We have used the abbreviation $`\varphi _{\mathrm{atm}}=\mathrm{\Delta }m_{\mathrm{atm}}^2L/2E`$. We conceive $`P_{\mathrm{CH}}`$ in Eq. (23) as the *true* survival probability, as opposed to the experimental value $`P_{\mathrm{CH}}^{\mathrm{exp}}=1.01`$ with the error $`\sigma _{\mathrm{CH}}=0.039`$ . Thus we have the range

$$0P_{\mathrm{CH}}12d_e(1d_e).$$
(24)
Making the substitutions

$$a0,b12d_e(1d_e),zP_{\mathrm{CH}},z_{\mathrm{exp}}P_{\mathrm{CH}}^{\mathrm{exp}},\sigma _z\sigma _{\mathrm{CH}},$$
(25)
the CHOOZ result can be included in the statistical analysis according to both methods described in the previous subsections.

For the atmospheric up-down inequality we substitute

$$aG,bH,zA_{\mathrm{ud}},z_{\mathrm{exp}}A_{\mathrm{ud}}^{\mathrm{exp}},\sigma _z=\sigma _A,$$
(26)
where $`G`$, $`H`$, $`A_{\mathrm{ud}}^{\mathrm{exp}}`$ and $`\sigma _A`$ are given by Eqs. (12), (13), (6) and (7), respectively. Notice that, since we have $`A_{\mathrm{ud}}^{\mathrm{exp}}0`$ and $`H0`$, the first term in the exponential of the expression (22) does not contribute here, and the first term in the expression (20) is 1 for all practical purposes.

The substitutions (25) and (26) allow us to define – according to Eqs. (16), (18), (20), (22) – the functions $`L_{\mathrm{CH}}`$, $`\pi _{\mathrm{CH}}`$, $`\mathrm{}_{\mathrm{CH}}`$, $`L_{\mathrm{CH}}^m`$ and $`L_{\mathrm{ud}}`$, $`\pi _{\mathrm{ud}}`$, $`\mathrm{}_{\mathrm{ud}}`$, $`L_{\mathrm{ud}}^m`$, referring to the CHOOZ experiment and the atmospheric up-down inequality, respectively. These functions will be used in the following discussions.

## 4 A bound on the LSND $`\underset{\mu }{\overset{()}{\nu }}\underset{e}{\overset{()}{\nu }}`$ transition amplitude

### 4.1 The statistical procedure

Let us now describe how to derive our desired probability distribution. We concentrate on the four parameters $`d_e`$, $`d_\mu `$, $`A_{\mathrm{ud}}`$ and $`P_{\mathrm{CH}}`$. We want to stress again that the latter two quantities are conceived as *true* values. Therefore, they are parameters within our procedure and we will treat them according to the methods described in Section 3, in order to arrive at a distribution solely in $`d_e`$ and $`d_\mu `$. The likelihood function is given by

$$L=L_{\mathrm{osc}}(d_e,d_\mu )\times L_{\mathrm{ud}}(A_{\mathrm{ud}})\times L_{\mathrm{CH}}(P_{\mathrm{CH}}),$$
(27)
where the first factor $`L_{\mathrm{osc}}`$ is the product of the likelihood functions of the Bugey, CDHS and KARMEN experiments. This likelihood function will be discussed in the next section and in Appendix B. For the treatment of parameters other than $`d_e`$ and $`d_\mu `$, which appear in the fitting procedure to the SBL experiments, see also Appendix B.

The physically allowed region $`_d`$ of $`d_e`$ and $`d_\mu `$ is described by the inequalities

$$_d:d_e0,d_\mu 0\text{and}d_e+d_\mu 1.$$
(28)
In order to incorporate it into our procedure, we define a function $`R(d_e,d_\mu )`$ such that $`R(d_e,d_\mu )=1`$ for $`(d_e,d_\mu )_d`$ and 0 otherwise. This function has thus the task of a prior which confines $`(d_e,d_\mu )`$ to the physically meaningful region. Adopting the maximum likelihood method of Section 3.2 in order to deal with the parameters $`A_{\mathrm{ud}}`$ and $`P_{\mathrm{CH}}`$ and combining this method with Eq. (15), we finally arrive at the desired probability distribution

$$p_m(d_e,d_\mu )=\frac{L_{\mathrm{osc}}(d_e,d_\mu )L_{\mathrm{ud}}^m(d_e,d_\mu )L_{\mathrm{CH}}^m(d_e)R(d_e,d_\mu )}{𝑑d_e^{}𝑑d_\mu ^{}L_{\mathrm{osc}}(d_e^{},d_\mu ^{})L_{\mathrm{ud}}^m(d_e^{},d_\mu ^{})L_{\mathrm{CH}}^m(d_e^{})R(d_e^{},d_\mu ^{})}.$$
(29)
The dependence of $`L_{\mathrm{ud}}^m(d_e,d_\mu )`$ and $`L_{\mathrm{CH}}^m(d_e)`$ on $`d_e`$ and $`d_\mu `$ comes in through the substitutions (26) and (25) of the boundaries $`a`$ and $`b`$.

Discussing now the Bayesian approach described in Section 3.1, we have the prior probability density

$$\pi (d_e,d_\mu ,A_{\mathrm{ud}},P_{\mathrm{CH}})=\frac{1}{2}R(d_e,d_\mu )\pi _{\mathrm{ud}}(d_e,d_\mu ,A_{\mathrm{ud}})\pi _{\mathrm{CH}}(d_e,P_{\mathrm{CH}}).$$
(30)
We stress that

$$𝑑A_{\mathrm{ud}}𝑑P_{\mathrm{CH}}\pi (d_e,d_\mu ,A_{\mathrm{ud}},P_{\mathrm{CH}})=\frac{1}{2}R(d_e,d_\mu )$$
(31)
is fulfilled. This equation tells us that, after integrating over the other variables, the prior probability density for $`d_e`$ and $`d_\mu `$ is uniform. Since these are the variables whose distribution we want to calculate, Eq. (31) assures us that with our choice of priors we have not introduced a bias in the distribution of $`d_e`$ and $`d_\mu `$.

According to the Bayesian approach in Section 3.1 we perform the integrations

$$𝑑A_{\mathrm{ud}}𝑑P_{\mathrm{CH}}L_{\mathrm{ud}}L_{\mathrm{CH}}\pi .$$
(32)
Finally, we obtain via this method the probability distribution

$$p_b(d_e,d_\mu )=\frac{L_{\mathrm{osc}}(d_e,d_\mu )\mathrm{}_{\mathrm{ud}}(d_e,d_\mu )\mathrm{}_{\mathrm{CH}}(d_e)R(d_e,d_\mu )}{𝑑d_e^{}𝑑d_\mu ^{}L_{\mathrm{osc}}(d_e^{},d_\mu ^{})\mathrm{}_{\mathrm{ud}}(d_e^{},d_\mu ^{})\mathrm{}_{\mathrm{CH}}(d_e^{})R(d_e^{},d_\mu ^{})}.$$
(33)
It is important to note that for every $`\mathrm{\Delta }m^2`$ we have such a distribution. The SBL mass-squared difference, which is hidden in $`L_{\mathrm{osc}}`$, is not on the same footing as $`d_e`$ and $`d_\mu `$ in our approach. For every given $`\mathrm{\Delta }m^2`$, we find restrictions on $`d_e`$ and $`d_\mu `$ from experiment555Note that in the analyses of the Bugey  and CDHS  experiments $`\mathrm{\Delta }m^2`$ is treated in the same way. which allow us finally to obtain a bound on the transition amplitude $`A_{\mu ;e}`$ for any given CL. Choosing a CL $`\beta `$, we find the corresponding bound on $`A_{\mu ;e}`$ by the prescription

$$\underset{4d_ed_\mu A_\beta ^0}{}𝑑d_e𝑑d_\mu p_j(d_e,d_\mu )=\beta (j=m,b).$$
(34)
For instance, if we choose $`\beta =0.99`$, with Eq. (34) we can calculate a number $`A_{0.99}^0`$ such that $`A_{\mu ;e}A_{0.99}^0`$ at 99% CL.

The two distributions $`p_m(d_e,d_\mu )`$ (29) and $`p_b(d_e,d_\mu )`$ (33) are, of course, different functions of the variables. The method to obtain bounds on $`A_{\mu ;e}`$ as given by Eq. (34) is the same for both distributions. In the following, we will comment on how much the bounds on $`A_{\mu ;e}`$ calculated with both distributions differ numerically.

### 4.2 A qualitative discussion of the SBL data

As mentioned previously, the SBL experiments are treated with the likelihood function

$$L_{\mathrm{osc}}(d_e,d_\mu )=L_{\mathrm{Bugey}}(d_e)\times L_{\mathrm{CDHS}}(d_\mu )\times L_{\mathrm{KARMEN}}(A_{\mu ;e}),$$
(35)
which enters Eq. (27). Here we want to discuss some features of $`L_{\mathrm{osc}}`$ and how these features influence the bound on $`A_{\mu ;e}`$. Let us start with the Bugey  and CDHS  experiments. Both are disappearance experiments and, therefore, the relevant SBL survival probabilities are given by

$$P_{\nu _\alpha \nu _\alpha }=P_{\overline{\nu }_\alpha \overline{\nu }_\alpha }=14d_\alpha (1d_\alpha )\mathrm{sin}^2\frac{\mathrm{\Delta }m^2L}{4E},$$
(36)
where $`\alpha =e`$ refers to the Bugey and $`\alpha =\mu `$ to the CDHS experiment. Both experiments have not found evidence for neutrino oscillations. In the case of the Bugey experiment the survival amplitude $`4d_e(1d_e)`$ is bounded by 0.1 or smaller values in the relevant range $`0.1\text{eV}^2\stackrel{<}{}\mathrm{\Delta }m^2\stackrel{<}{}\mathrm{\hspace{0.25em}10}\text{eV}^2`$ at 90% CL . Thus, $`d_e`$ is either very small or close to 1. Since for solar neutrinos the inequality

$$P_{\nu _e\nu _e}^{}d_e^2$$
(37)
holds , $`d_e`$ must be small666We want to stress that this is the only place where solar neutrino data enter our analysis. Moreover, our inference that $`d_e`$ must be small is independent of the actual solution to the solar neutrino problem. and in the fit of the Bugey data we make the approximation $`d_e(1d_e)d_e`$ in the survival probability. In principle, the inequality (37) should be included in our analysis with the help of one of the methods discussed in Sections 3.1 and 3.2, but in view of the smallness of the Bugey survival amplitude this is not necessary. In the case of the CDHS experiment, for $`\mathrm{\Delta }m^2\stackrel{>}{}\mathrm{\hspace{0.25em}0.3}`$ eV2, there is an analogous feature concerning $`d_\mu `$ ; there, the selection of the small $`d_\mu `$ is guaranteed by the up-down inequality (12) . For smaller $`\mathrm{\Delta }m^2`$ the CDHS restriction disappears and values of $`d_\mu `$ as large as 0.5 are allowed by the up-down inequality. Therefore, in $`L_{\mathrm{CDHS}}`$ we do not neglect $`d_\mu ^2`$ in the survival amplitude because this would not be justified.

In the KARMEN experiment, $`\overline{\nu }_\mu \overline{\nu }_e`$ transitions have not been observed and the KARMEN exclusion curve  cuts right through the region allowed by LSND. Therefore, it is important to take into account the KARMEN result in our analysis. For the details of our fit to the KARMEN data and the analogous details concerning the Bugey and CDHS experiments see Appendix B.

We include the CHOOZ result in our analysis for $`\mathrm{\Delta }m^20.05`$ eV2, where the bound on $`\mathrm{sin}^22\theta `$ in the CHOOZ plots is a straight line ; this guarantees that oscillations with the SBL $`\mathrm{\Delta }m^2`$ are averaged out and hence Eq. (23) is valid. The CHOOZ result has an effect on the $`A_{\mu ;e}`$ exclusion curve mainly for $`\mathrm{\Delta }m^2`$ where the Bugey bound on $`d_e`$ is not so strong, which is the case for $`\mathrm{\Delta }m^2`$ around 0.1 eV2 and $`\mathrm{\Delta }m^2\stackrel{>}{}\mathrm{\hspace{0.25em}5}`$ eV2.

### 4.3 The bound on $`A_{\mu ;e}`$

Fig. 1 represents the main result of this paper. In this figure we show the regions in the $`\mathrm{\Delta }m^2`$$`A_{\mu ;e}`$ plane allowed by LSND at 90% and 99% CL  and our exclusion curves for $`A_{\mu ;e}`$ with (3+1)-mass spectra at 90%, 95% and 99% CL. These exclusion curves have been calculated using the maximum likelihood approach of Section 3.2 for the treatment of the *true* values $`A_{\mathrm{ud}}`$ and $`P_{\mathrm{CH}}`$. This method gives – for most values of $`\mathrm{\Delta }m^2`$ – slightly weaker bounds than the Bayesian approach of Section 3.1; however, the numerical difference between the two methods is always smaller than 4%. We read off from Fig. 1 that the 99% CL LSND region has no overlap with the allowed region to the left of our 95% CL bound. This shows that the likelihood of (3+1)-mass spectra is not very high. Comparing the 99% CL exclusion curve with the 99% CL LSND region and confining ourselves to $`\mathrm{\Delta }m^2<10`$ eV2, we find overlaps at $`\mathrm{\Delta }m^26,\mathrm{\hspace{0.33em}1.7},\mathrm{\hspace{0.33em}0.9}`$ and approximately between 0.25 and 0.4 eV2. This agrees with the findings in Ref. , where 90% exclusion curves and the bound on $`A_{\mu ;e}`$ derived in Ref.  are compared with the 99% CL LSND region. The result of the inclusion of the CHOOZ result for $`\mathrm{\Delta }m^20.05`$ eV2 can be seen as a jump in our exclusion curves. From this jump the effect of the CHOOZ result can be read off for small mass-squared differences. Let us stress once more that in our analysis an exclusion curve for a given CL is the result of a well-defined statistical procedure and is obtained by including all available data other than the LSND data; such exclusion curves have a precise statistical meaning; an overlap of the 99% CL LSND region with the region to the left of an exclusion curve occurs only with the exclusion curve at 99% CL, but not at 95% CL.

There are claims in the literature that all three indications in favour of neutrino oscillations are compatible with *three* neutrinos . Such claims have been refuted by a combined fit to the data . In our analysis we arrive at the same conclusion, because our treatment of the (3+1)-spectra is also applicable in the case of three neutrinos, where one has (2+1)-mass spectra, since there is a SBL mass-squared difference $`\mathrm{\Delta }m^2`$ and another mass-squared difference much smaller than $`\mathrm{\Delta }m^2`$; assuming that $`m_3`$ is the mass separated from the other two by the LSND mass gap, we have to make the identification $`d_\alpha =|U_{\alpha 3}|^2`$ ($`\alpha =e,\mu `$). The values of $`\mathrm{\Delta }m^2`$ and $`A_{\mu ;e}`$, indicated in the papers of Ref.  as preferred solutions, lie in the region allowed by LSND at 99% and thus beyond our 95% CL exclusion curve in Fig. 1.

## 5 The (2+2)-neutrino mass spectra

If we want to treat the (2+2)-neutrino mass spectra with the method discussed in the previous two sections, the situation is more involved. In this case the $`\nu _\mu \nu _e`$ transition amplitude is given by

$$A_{\mu ;e}=4\left|\underset{j=1,2}{}U_{ej}U_{\mu j}^{}\right|^2=4\left|\underset{j=3,4}{}U_{ej}U_{\mu j}^{}\right|^2$$
(38)
and cannot be expressed by

$$d_\alpha =\underset{j=3,4}{}|U_{\alpha j}|^2(\alpha =e,\mu ),$$
(39)
as in the case of the (3+1)-spectra (see Eq. (4)). This suggests to perform an analysis with the five parameters $`d_e`$, $`d_\mu `$, $`A_{\mu ;e}`$, $`A_{\mathrm{ud}}`$ and $`P_{\mathrm{CH}}`$ in the case of the (2+2)-spectra.

Let us explore in more detail the relationship between the amplitude $`A_{\mu ;e}`$ and the parameters $`d_e`$ and $`d_\mu `$. From Eq. (38), using the Cauchy–Schwarz inequality, we readily obtain the inequality

$$A_{\mu ;e}4\mathrm{min}[d_ed_\mu ,(1d_e)(1d_\mu )].$$
(40)
This inequality implies that for every $`A_{\mu ;e}`$ the following region is allowed in the $`d_e`$$`d_\mu `$ plane:

$$(A_{\mu ;e}):\{\begin{array}{ccccc}\hfill \frac{1}{2}(1\sqrt{1A_{\mu ;e}})& & d_e& & \frac{1}{2}(1+\sqrt{1A_{\mu ;e}}),\hfill \\ \hfill A_{\mu ;e}/4d_e& & d_\mu & & 1A_{\mu ;e}/4(1d_e).\hfill \end{array}$$
(41)
For every $`A_{\mu ;e}`$ between 0 and 1, this is a region in the unit square, confined by two hyperbolas. Let us label the neutrino masses such that $`\mathrm{\Delta }m_{21}^2=\mathrm{\Delta }m_{\mathrm{atm}}^2`$ and $`\mathrm{\Delta }m_{43}^2=\mathrm{\Delta }m_{\mathrm{solar}}^2`$. Then, with the definitions (39), Eqs. (9) and (11) hold also for the (2+2)-mass spectra and we obtain the inequalities

$$A_{\mathrm{ud}}G^{}(d_e,d_\mu ,A_{\mu ;e})=\frac{(1d_\mu )^2+A_{\mu ;e}/4r}{(1d_\mu )^2+2d_\mu ^2+3A_{\mu ;e}/4r}$$
(42)
and

$$A_{\mathrm{ud}}H^{}(d_e,d_\mu ,A_{\mu ;e})=\frac{d_\mu (1d_\mu )A_{\mu ;e}/4r}{1d_\mu (1d_\mu )+A_{\mu ;e}/4r}$$
(43)
for the atmospheric up-down asymmetry. Furthermore, the survival probability for solar neutrinos is bounded by $`P_{\nu _e\nu _e}^{}\frac{1}{2}(1d_e)^2`$. Therefore, $`d_e`$ is close to one777If $`d_e`$ is small, then it follows that $`P_{\nu _e\nu _e}^{}\stackrel{>}{}\mathrm{\hspace{0.25em}0.5}`$, which is in disagreement with the result of the Homestake experiment (see first paper in Ref. ); moreover, for small $`d_e`$ we obtain the same bound on $`A_{\mu ;e}`$ as for the disfavoured (3+1)-mass spectra. and we make the approximation $`d_e(1d_e)1d_e`$ in the Bugey survival amplitude. Note that the CHOOZ inequality (23) also holds for the (2+2)-spectra.

In analogy to Eq. (34) in the (3+1)-case, our aim is to obtain a probability distribution in $`A_{\mu ;e}`$. Here we confine ourselves to the maximum likelihood method of Section 3.2 for the treatment of inequalities. In addition to the up-down and CHOOZ inequalities we have to take into account Eq. (41). According to the method of Section 3.2, we maximize with respect to $`A_{\mathrm{ud}}`$, $`P_{\mathrm{CH}}`$, $`d_e`$ and $`d_\mu `$. The maximization with respect to the first two variables leads to the functions $`L_{}^{}{}_{\mathrm{ud}}{}^{m}`$ and $`L_{\mathrm{CH}}^m`$, which are defined in Section 3.3. The prime on the first function indicates that now it is obtained from $`L_{\mathrm{ud}}^m`$ by the replacement $`GG^{}`$ and $`HH^{}`$. Generalizing the method of Section 3.2, we define

$$L^{}(A_{\mu ;e})=\underset{(d_e,d_\mu )(A_{\mu ;e})}{\mathrm{max}}\left[L_{\mathrm{osc}}^{}L_{}^{}{}_{\mathrm{ud}}{}^{m}L_{\mathrm{CH}}^m\right],$$
(44)
where $`L_{\mathrm{osc}}^{}`$ is given by Eq. (35), but in view of the appearance of $`A_{\mu ;e}`$ in $`L_{\mathrm{KARMEN}}`$ it is interpreted as a function of three variables. Finally, we obtain the probability distribution for $`A_{\mu ;e}`$ given by

$$p_m^{}(A_{\mu ;e})=\frac{L^{}(A_{\mu ;e})\mathrm{\Theta }(A_{\mu ;e})\mathrm{\Theta }(1A_{\mu ;e})}{𝑑A_{\mu ;e}^{}L^{}(A_{\mu ;e}^{})\mathrm{\Theta }(A_{\mu ;e}^{})\mathrm{\Theta }(1A_{\mu ;e}^{})}.$$
(45)
In Fig. 2 the 90% CL bound on $`A_{\mu ;e}`$ calculated with the distribution (45) is plotted. In this figure the 90% and 99% CL level regions of LSND are also depicted. The region to the left of the bound has an overlap area with the LSND region of 90% CL. This shows that the (2+2)-spectra are well compatible with all data. Also shown in Fig. 2 are the KARMEN and Bugey exclusion curves. As expected, for large $`\mathrm{\Delta }m^2`$ the solid line of our exclusion curve follows rather well the KARMEN exclusion curve, whereas for small $`\mathrm{\Delta }m^2`$ it tends to follow the Bugey curve. Due to the inclusion of the CHOOZ result above 0.05 eV2, the solid line is more restrictive there than the Bugey curve, whereas below 0.05 eV2 both curves are nearly identical. We have also investigated the exclusion curve of the (2+2)-spectra within the Bayesian approach to inequalities of Section 3.1. In this approach we have more freedom to make “reasonable” choices for priors than in the (3+1)-case, and, numerically, the $`A_{\mu ;e}`$ bounds tend to be a little more restrictive than the solid line in Fig. 2.

## 6 Conclusions

In this paper we have examined the (3+1) and (2+2)-classes of the 4-neutrino mass spectra within a rigorous statistical analysis. Since we do not have sufficient information concerning the final LSND data , we have chosen the approach to analyse all other available data and compare our result with the LSND result in a $`\mathrm{\Delta }m^2`$$`A_{\mu ;e}`$ plot. This approach suggests to extract that information from the solar and atmospheric data which is most relevant for SBL oscillations with respect to $`\mathrm{\Delta }m^2`$. This extraction is most appropriately done in the form of inequalities involving elements of the neutrino mixing matrix . Similarly, we have extracted the “SBL information” contained in the CHOOZ result. On the other hand, concerning the SBL experiments, we have fully included the data of the Bugey, CDHS and KARMEN experiments.

Since our aim has been to incorporate in the statistical analysis the atmospheric up-down inequality, Eq. (12) or (42), and the CHOOZ inequality (23), we have made use of the likelihood function in combination with the Bayesian approach to probability which allows us to derive probability distributions of the parameters which are to be estimated. In this context we have presented two possibilities for including inequalities involving parameters: one way is to treat them in the form of prior probability densities for which a “reasonable choice” has to be made; another way is a kind of maximum likelihood treatment. Numerically, we have compared both methods for exclusion curves in the case of the (3+1)-spectra and found that the difference is negligible.

With the method described in the paper we have obtained for every $`\mathrm{\Delta }m^2`$ a probability distribution for the SBL transition amplitude $`A_{\mu ;e}`$, from which we could derive bounds as a function of $`\mathrm{\Delta }m^2`$ on this amplitude for any CL. The results are shown in Figs. 1 and 2 for the (3+1) and (2+2)-neutrino mass spectra, respectively. In the latter case our 90% CL exclusion curve is close to the KARMEN exclusion curve down to $`\mathrm{\Delta }m^20.5`$ eV2; there it turns off and starts to come close to the exclusion curve given by the Bugey data. Thus, for the (2+2)-spectra our method reproduces more or less what one obtains by naively comparing the KARMEN and Bugey exclusion curves with the region allowed by LSND. Therefore, these spectra are well compatible with all the data. On the other hand, in the case of the (3+1)-spectra, our 95% CL bound has no overlap with the LSND region of 99% CL. We, therefore, come to the conclusion that this spectral class is rather unlikely, even with the recent change in the LSND region.

Thus we strengthen with the method presented here the claims made in Refs. . Should the LSND result be confirmed by the MiniBooNE Collaboration , then, in a 4-neutrino scheme, the sterile neutrino should make its appearance either in the solar or atmospheric neutrinos, or both.

## Acknowledgements

We thank C. Giunti and S.M. Bilenky for many useful discussions. We are very grateful to F. Dydak and J. Wotschack for valuable information about the CDHS experiment and to C. Walter for providing us with the latest results of the Super-Kamiokande up-down asymmetry. T.S. is supported by the Austrian Academy of Sciences.

## Appendix A The relationship between the Bayesian and the maximum likelihood methods for the treatment of inequalities

Here we want to elucidate the relationship between the methods for the inclusion of inequalities introduced in Sections 3.1 and 3.2. In the Bayesian approach to this problem we consider in general integrals of the type

$$I[\pi ]=_a^b𝑑zL_Z(z)\pi (z),$$
(46)
where $`L_Z`$ is given by Eq. (16). A prior $`\pi `$ is a positive and piecewise continuous function such that $`_a^b𝑑z\pi (z)=1`$. In Section 3.1 we have used the flat prior (18). If we denote the set of all prior probability densities on the interval $`[a,b]`$ by $`𝒫`$, then the following proposition holds.
Proposition: For the integral in Eq. (46) one has the upper limit

$$\underset{\pi 𝒫}{\mathrm{max}}I[\pi ]=L_Z^m,$$
(47)
where $`L_Z^m`$ is given by Eq. (22).
Proof: The proof is very straightforward. Let us first assume that $`z_{\mathrm{exp}}>b`$. Then we have

$$I[\pi ]<L_Z(b)_a^b𝑑z\pi (z)=L_Z(b).$$

Thus we have found an upper bound on $`I`$. If we can find a sequence of $`\pi _k`$ of prior probability densities such that

$$\underset{k\mathrm{}}{lim}I[\pi _k]=L_Z(b),$$

then the proposition is proven for the case $`z_{\mathrm{exp}}>b`$. Such a sequence is easy to give: any sequence $`\pi _k`$, for which $`\pi _k(z)=0`$ for $`zb1/k`$ holds, fulfills our purpose. For $`z_{\mathrm{exp}}<a`$ we proceed analogously. If $`z_{\mathrm{exp}}[a,b]`$, the upper bound on $`I`$ is $`L_Z(z_{\mathrm{exp}})`$. The sequence $`\pi _k`$ is constructed accordingly with the idea given before. Thus for the three cases $`z_{\mathrm{exp}}>b`$, $`z_{\mathrm{exp}}<a`$ and $`z_{\mathrm{exp}}[a,b]`$ we have the upper limits $`L_Z(b)`$, $`L_Z(a)`$ and $`L_Z(z_{\mathrm{exp}})`$, respectively. In summary, we just have obtained the expression (22). Q.E.D.

Thus the expression $`L_Z^m`$ is obtained from the Bayesian approach as the maximum over all possible priors. Sloppily speaking, the prior in $`I`$, which gives $`L_Z^m`$, is a delta function ($`\delta (zb)`$, $`\delta (za)`$ or $`\delta (zz_{\mathrm{exp}})`$); however, in order to choose the correct delta function one has to know the experimental value $`z_{\mathrm{exp}}`$ and, therefore, such a function does not deserve the name “prior” anymore. Though $`\mathrm{}_Z<L_Z^m`$ holds, it is not obvious which method gives the stronger constraint in an actual situation, because one has additional factors in the probability distribution one aims at, coming from additional data, and one has to normalize the combined distribution (see Eq. (15)). Moreover, the bounds $`a`$ and $`b`$ are in general functions of the parameters whose distribution we want to know. In the concrete situation discussed in Section 4 the maximum likelihood method represented by $`L_Z^m`$ gives a slightly weaker restriction on the transition amplitude $`A_{\mu ;e}`$ for most $`\mathrm{\Delta }m^2`$.

## Appendix B The analyses of the Bugey, CDHS and KARMEN experiments

In this appendix we describe how the data of the SBL experiments Bugey, CDHS and KARMEN is included in our analysis. We use as much information as can be recovered from the publications of the experimental groups to perform a fit to the data. The fact that we can reproduce the published 90% CL bounds in the case of 2-neutrino oscillations to a good accuracy inspires confidence in our analysis.

### B.1 Bugey

The Bugey experiment  searches for $`\overline{\nu }_e`$ disappearance at the three distances 15 m, 40 m and 95 m away from a nuclear reactor. The electron antineutrinos are detected through the reaction $`\overline{\nu }_e+pe^++n`$. As input data for our analysis we use Fig. 17 of Ref. , where the ratios of the observed events to the number of expected events in case of no oscillations are shown for the three positions in bins of positron kinetic energy in the range 1 MeV $`E_{e^+}`$ 6 MeV.

For the analysis we follow Eq. (9) of Ref.  and use the $`\chi ^2`$-function

$$\chi ^2=\underset{j}{}\left\{\underset{i=1}{\overset{N_j}{}}\frac{\left[(Aa_j+b(E_{ji}E_0))R_{ji}^{\mathrm{theo}}R_{ji}^{\mathrm{exp}}\right]^2}{\sigma _{ji}^2}+\frac{(a_j1)^2}{\sigma _a^2}\right\}+\frac{(A1)^2}{\sigma _A^2}+\frac{b^2}{\sigma _b^2}.$$
(48)
Here $`j=15,40,95`$ labels the three positions, $`i`$ the positron energy bins and $`N_{15}=N_{40}=25,N_{95}=10`$ are the numbers of bins at each position. For $`E_{ji}`$ we take the mean positron energy in bin $`ji`$. $`R_{ji}^{\mathrm{exp}}`$ is the ratio of measured to expected events in each bin with its statistical error $`\sigma _{ji}`$, both read off from Fig. 17 of Ref. . $`R_{ji}^{\mathrm{theo}}`$ is the theoretical prediction for this ratio in the case of oscillations, depending on the oscillation parameters, and we set $`R_{ji}^{\mathrm{theo}}=P_{\overline{\nu }_e\overline{\nu }_e}_{ji}`$, where $`P_{\overline{\nu }_e\overline{\nu }_e}`$ is given in Eq. (36). Various systematic uncertainties are taken into account by minimizing the $`\chi ^2`$, Eq. (48), with respect to the five parameters $`A`$, $`a_j`$ and $`b`$ for given oscillation parameters. In Ref.  the values $`\sigma _A=4.796\%,\sigma _a=1.414\%,\sigma _b=0.02`$ MeV-1 and $`E_0=1`$ MeV are given.

To perform the averaging of the survival probability we estimate the uncertainty in the flight length of the neutrinos because of the size of the production region and the detector to 3 m and we assume that the flux varies with the distance as $`L^2`$. In the relevant energy range antineutrino and positron energy are related by $`E_\nu =E_{e^+}+1.8`$ MeV to a very good approximation. For the purpose of averaging the survival probability over the energy range in one bin it is a good approximation to take neutrino flux, detection cross section and efficiency as constant with energy. The reason for this is that the energy bins are relatively small and only ratios of observed to expected events in each bin are considered. Furthermore, we assume a Gaussian resolution function for the positron energy measurement with variance 0.4 MeV.

The likelihood function in Eq. (35) which contains the information of the Bugey experiment is obtained from the $`\chi ^2`$ of Eq. (48) by

$$L_{\mathrm{Bugey}}(d_e)\mathrm{exp}\left(\frac{1}{2}\chi ^2\right).$$
(49)
The 90% CL bound on $`d_e`$ obtained from this likelihood function alone can be compared to the Bugey exclusion curve in the 2-neutrino case with the identification $`\mathrm{sin}^22\theta _{\mathrm{Bugey}}4d_e`$. The curve obtained by our analysis is shown in Fig. 2 and compares well with the originally published one .

In Eq. (36) contributions of oscillations because of $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ to the SBL disappearance amplitude are neglected. This approximation may not be exactly fulfilled in the Bugey experiment for small $`\mathrm{\Delta }m_{\mathrm{SBL}}^2`$ and large $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$. We have calculated the 90% CL bound from Bugey by taking into account also oscillations with $`\mathrm{\Delta }m_{\mathrm{atm}}^2=6\times 10^3`$ eV2, which is the 99% CL upper bound on $`\mathrm{\Delta }m_{\mathrm{atm}}^2`$ , and find that for $`\mathrm{\Delta }m_{\mathrm{SBL}}^2>0.04`$ eV2 the effect is smaller than 6%.

### B.2 CDHS

The CDHS experiment  searches for $`\nu _\mu `$ disappearance by comparing the number of events in the so-called back and front detectors at the distances $`L_{\mathrm{back}}=885`$ m and $`L_{\mathrm{front}}=130`$ m, respectively, from the neutrino source. The neutrinos are detected via muons produced in charged-current interactions. The data is given in form of the double ratios

$$R_{\mathrm{corr}}=\frac{N_{\mathrm{back}}/N_{\mathrm{front}}}{N_{\mathrm{back}}^{\mathrm{MC}}/N_{\mathrm{front}}^{\mathrm{MC}}},$$
(50)
where $`N_{\mathrm{back}}(N_{\mathrm{front}})`$ is the number of observed events in the back (front) detector and $`N_{\mathrm{back}}^{\mathrm{MC}}`$ and $`N_{\mathrm{front}}^{\mathrm{MC}}`$ are the corresponding quantities expected for no oscillations calculated by Monte Carlo.

The ratios (50) are given in 15 bins of “projected range in iron”. This is the distance traveled by the muon in the detector (consisting of iron) projected onto the detector axis, which has an angle of $`22^{}`$ relative to the neutrino beam axis. We calculate the range in iron $`r(E_\mu )`$ of a muon with energy $`E_\mu `$ by integrating Eq. (23.1) of Ref. . The muon energy intervals $`[E_\mu ^{(i1)},E_\mu ^{(i2)}]`$ corresponding to the intervals of projected range ($`r_{\mathrm{proj}}`$) for bin $`i`$ given in Table 1 of Ref.  are obtained by applying the relation $`r_{\mathrm{proj}}(E_\mu )r(E_\mu )\mathrm{cos}\vartheta \mathrm{cos}22^{}`$, where $`\vartheta 20^{}`$ is the average scattering angle of the muons .

We estimate the number of events in bin $`i`$ and at position $`p=`$back or front using

$$N_{ip}M_p_{E_\mu ^{(i1)}}^{E_\mu ^{(i2)}}𝑑E_\mu _{E_\mu }^{\mathrm{}}𝑑E_\nu _{L_p\mathrm{\Delta }L_p/2}^{L_p+\mathrm{\Delta }L_p/2}𝑑LL^2P_{\nu _\mu \nu _\mu }(L/E_\nu )\mathrm{\Phi }(E_\nu )\frac{d\sigma _{\mathrm{DIS}}(E_\nu ,E_\mu )}{dE_\mu }.$$
(51)
Here $`M_p`$ is the detector mass at position $`p`$. Taking into account the size of the decay tunnel (52 m) and the length relevant for detection in the back/front detector (72 m/22 m)  we have $`\mathrm{\Delta }L_{\mathrm{back}}=52+72`$ m and $`\mathrm{\Delta }L_{\mathrm{front}}=52+22`$ m. The disappearance probability $`P_{\nu _\mu \nu _\mu }(L/E_\nu )`$ depending on the oscillation parameters is given in Eq. (36) and $`\sigma _{\mathrm{DIS}}`$ is the cross section for the deep inelastic scattering process $`\nu _\mu +N\mu +X`$ (see, e.g., Ref. ). The neutrino flux $`\mathrm{\Phi }(E_\nu )`$ in the relevant neutrino energy range $`E_\nu \stackrel{>}{}\mathrm{\hspace{0.25em}0.7}`$ GeV is proportional to $`\mathrm{exp}(E_\nu /1\mathrm{GeV})`$ . Finally, we obtain for the theoretical prediction for the double ratios (50) in bin $`i`$

$$R_{\mathrm{theo}}^i=\frac{N_{i\mathrm{back}}}{N_{i\mathrm{front}}}\left(\frac{L_{\mathrm{back}}}{L_{\mathrm{front}}}\right)^2\frac{M_{\mathrm{front}}}{M_{\mathrm{back}}},$$
(52)
and we define the CDHS likelihood function by

$$L_{\mathrm{CDHS}}(d_\mu )\mathrm{exp}\left[\frac{1}{2}\underset{ij}{}(R_{\mathrm{corr}}^iR_{\mathrm{theo}}^i)(S^1)_{ij}(R_{\mathrm{corr}}^jR_{\mathrm{theo}}^j)\right].$$
(53)
Assuming total correlation between any two bins, the covariance matrix is given by $`S_{ij}=\delta _{ij}\sigma _i^2+\sigma _{\mathrm{syst}}^2`$. $`R_{\mathrm{corr}}^i`$ and its statistical error $`\sigma _i`$ are read off from Table 1 of Ref. . The overall systematic error in the ratio of event rates in the two detectors was estimated to $`\sigma _{\mathrm{syst}}=2.5`$% in Ref. .

The 90% CL bound on $`\mathrm{sin}^22\theta _{\mathrm{CDHS}}=4d_\mu (1d_\mu )`$ obtained by our analysis is shown in Fig. 3. It is very similar to the bound published by the CDHS collaboration . There are minor differences for small mass squared differences, which could have some effect for our bound on $`A_{\mu ;e}`$ in the region 0.2 eV$`{}_{}{}^{2}\stackrel{<}{}\mathrm{\Delta }m^2\stackrel{<}{}\mathrm{\hspace{0.25em}0.4}`$ eV2: our bound disappears at $`\mathrm{\Delta }m^20.3`$ eV2, whereas the CDHS bound extends down to approximately 0.24 eV2.

### B.3 KARMEN

The latest data of the $`\overline{\nu }_\mu \overline{\nu }_e`$ oscillation search in the KARMEN experiment is presented in Ref. . Analysing the data taken from February 1997 to March 2000 they find a total of 11 candidate events, in good agreement with the expected number of background events for no oscillations of $`12.3\pm 0.6`$. For our analysis we use the data resulting from the detection process $`\overline{\nu }_e+pe^++n`$. The positron spectrum $`S(E_{e^+})`$ expected for $`A_{\mu ;e}\mathrm{sin}^22\theta _{\mathrm{KARMEN}}=1`$ and $`\mathrm{\Delta }m^2=100`$ eV2 we take from Fig. 2(a) of Ref. . To estimate the number of events in a positron energy interval $`[E_1,E_2]`$ resulting from neutrino oscillations we follow Eq. (B1) of Ref. :

$$N^{\mathrm{osc}}=N_{E_1}^{E_2}𝑑E_{e^+}S(E_{e^+})_{L_1}^{L_2}𝑑LL^2P_{\overline{\nu }_\mu \overline{\nu }_e}(L/E_\nu ),$$
(54)
where the oscillation probability $`P_{\overline{\nu }_\mu \overline{\nu }_e}`$ is given in Eq. (3). For $`L_{1,2}`$ we take $`17.51.75`$ m  and antineutrino energy and positron kinetic energy are related by $`E_\nu =E_{e^+}+1.8`$ MeV. The normalization factor $`N`$ is fixed by requiring that for the total positron energy range $`E_1=16`$ MeV, $`E_2=52`$ MeV, full mixing ($`\mathrm{sin}^22\theta _{\mathrm{KARMEN}}=1`$) and $`\mathrm{\Delta }m^2100`$ eV2 the number of events resulting from oscillations is 2442 (see Table 1 of Ref. ).

From Fig. 2(b) of Ref.  we read off for each of the 9 positron energy bins the number of observed events $`N_i^{\mathrm{obs}}`$ and the number of background events expected for no oscillations $`B_i`$. Then we construct the likelihood function by using the Poisson distribution:

$$L_{\mathrm{KARMEN}}(A_{\mu ;e})=\underset{i=1}{\overset{9}{}}\frac{1}{N_i^{\mathrm{obs}}!}(N_i^{\mathrm{osc}}+B_i)^{N_i^{\mathrm{obs}}}e^{(N_i^{\mathrm{osc}}+B_i)}.$$
(55)
Here $`N_i^{\mathrm{osc}}`$ is calculated from Eq. (54) by choosing $`E_1`$ and $`E_2`$ according to the bin $`i`$. The 90% CL bound on $`A_{\mu ;e}`$ for a given $`\mathrm{\Delta }m^2`$ obtained from the probability distribution implied by Eq. (55) in the Bayesian approach is shown in Fig. 2. Our bound is very close to the one presented in Ref. .


format: text

text:
# Structure of wavefunctions in (1+2)-body random matrix ensembles

## I Introduction

In the last few years the study of quantum chaos in isolated finite interacting particle systems has shifted from spectral statistics to properties of wavefunctions and transition strengths (for example, electromagnetic and Gamow-Teller transition strengths in atomic nuclei, dipole strengths in atoms etc.). Working in this direction, several research groups have recognized recently that the two-body random matrix ensembles and their various extended versions form good models for understanding various aspects of chaos in interacting particle systems . In particular using the so called EGOE(1+2), embedded Gaussian orthogonal ensemble of (1+2)-body interactions defined by a mean-field one-body plus a chaos generating random two-body interaction, there are now several studies on the nature of occupancies of single particle states, strength functions (or local density of states), information entropy, transition strength sums and transition matrix elements of one-body transition operators, Fock-space localization etc., in the chaotic domain of interacting particle systems such as atoms , nuclei , quantum dots , quantum computers  and so on. Ref.  gives a recent overview of this subject. The focus in the present article is on two important measures of localization (in wavefunctions and transition strength distributions): (i) number of principal components NPC (or the inverse participation ratio IPR); (ii) localization length $`𝒍_H`$ as defined by the information entropy (Sinfo). It is well established that the NPC in wavefunctions characterizes various layers of chaos in interacting particle systems . In addition, for systems such as atomic nuclei, NPC for transition strengths is a measure of fluctuations in transition strength sums. Similarly the role of $`𝒍_H`$ in quantum chaos studies is well emphasized by Izrailev  and more significantly, using nuclear physics examples  it is well demonstrated that the wavefunction entropy Sinfo coincides with the thermodynamic entropy for many particle systems with two-body interactions of sufficient strength but only in the presence of a mean-field, i.e. in the chaotic domain but with a mean-field - therefore the significance of EGOE(1+2). Clearly deriving the predictions of EGOE(1+2) for NPC and $`𝒍_H`$ are of considerable importance. This problem was addressed in Refs. . In  results for NPC in wavefunctions, in the so called Breit-Wigner (BW) domain, are derived. On the other hand, in  results in the so called Gaussian domain (the BW and Gaussian domains are defined in Sect. II ahead) are derived for NPC and $`𝒍_H`$ in transition strength distributions with only the final results mentioned for wavefunctions. The purpose of the present paper is to give a detailed derivation of the results mentioned in  for NPC and $`𝒍_H`$ in wavefunctions and subject them to numerical tests. Now we will give a preview.

Section II gives some of the basic results for EGOE(1+2). In Section III, formulas for NPC and $`𝒍_H`$ in wavefunctions are derived by exploiting the Gaussian nature and the associated properties of strength functions in EGOE(1+2). Numerical tests of the theory are given in Section IV. Finally Section V gives concluding remarks.

## II Basic results for (1+2)-body random matrix ensembles

Given $`m`$ fermions in $`N`$ single particle states, assuming at the outset that the many particle spaces are direct product spaces of the single particle states, two-body random matrix ensembles (usually called TBRE) are generated by defining the hamiltonian $`H`$, which is $`2`$-body, to be a random matrix in the $`2`$-particle spaces and then propagating it to the $`\left(\begin{array}{c}N\\ m\end{array}\right)`$ dimensional $`m`$ \- particle spaces by using their geometry (direct product structure); often one considers a GOE representation in the $`2`$-particle spaces and then the TBRE is called EGOE(2); see  for more details. For a EGOE(2), with $`N>>m>>2`$, the normalized state density $`\rho (E)=\delta (HE)`$ takes Gaussian form and it is defined by its centroid $`ϵ=H`$ and variance $`\sigma ^2=(Hϵ)^2`$. In order to explicitly state that the state density is generated by the hamiltonian $`H`$, sometimes $`\rho (E)`$ is denoted as $`\rho ^H(E)`$ and similarly $`ϵ`$ as $`ϵ_H`$ and $`\sigma `$ as $`\sigma _H`$. Note that the averages $``$ are over the $`m`$-particle spaces ; in the nuclear physics examples, they are usually over the $`m`$-particle spaces with fixed angular momentum ($`J`$) and isospin ($`T`$) which are good quantum numbers. Just as the state density, given a transition operator $`𝒪`$, the normalized bivariate strength densities (matrix elements of $`𝒪`$ weighted by the state densities at the initial and final energies) $`\rho _{biv;𝒪}(E_i,E_f)=\left[𝒪^{}𝒪\right]^1𝒪^{}\delta (HE_f)𝒪\delta (HE_i)`$ take bivariate Gaussian form for EGOE(2) and it is defined by the centroids ($`ϵ_i`$, $`ϵ_f`$) and widths ($`\sigma _i`$, $`\sigma _f`$) of its two marginals and the bivariate correlation coefficient which is given by $`𝒪^{}[(Hϵ_f)/\sigma _f]𝒪[(Hϵ_i)/\sigma _i]/𝒪^{}𝒪`$. Thirdly, the level and strength fluctuations follow GOE. Moreover, with the Gaussian form for the state densities and bivariate Gaussian form for the strength densities, the strength sums $`E𝒪^{}𝒪E=_E^{}\left|E^{}𝒪E\right|^2`$ take the form of ratio of two Gaussians, $`E𝒪^{}𝒪E=𝒪^{}𝒪\rho _{𝒪^{}𝒪:𝒢}(E)/\rho _𝒢(E)`$ where $`\rho _{𝒪^{}𝒪:𝒢}(E)=𝒪^{}𝒪\delta (HE)`$ is defined by its centroid $`ϵ_{𝒪^{}𝒪}=𝒪^{}𝒪H/𝒪^{}𝒪`$ and variance $`\sigma _{𝒪^{}𝒪}^2=𝒪^{}𝒪H^2/𝒪^{}𝒪ϵ_{𝒪^{}𝒪}^2`$; $`𝒢`$ stands for Gaussian.

Hamiltonians for realistic interacting particle systems contain a mean-field part (one-body part $`h(1)`$) and a two-body residual interaction $`V(2)`$ mixing the configurations built out of the distribution of particles in the mean-field single particle states; $`h(1)`$ is defined by the single particle energies $`ϵ_i,i=1,2,\mathrm{},N`$ and $`V(2)`$ is defined by its two-particle matrix elements. Then it is more realistic to use EGOE(1+2), the embedded Gaussian orthogonal ensemble of random matrices of (1+2)-body hamiltonians where $`\{H\}=h(1)+\lambda \{V(2)\}`$; sometimes it is more convenient to use $`\alpha h(1)+\lambda \{V(2)\}`$. Here $`\{\}`$ denotes ensemble, $`\lambda `$ and $`\alpha `$ are free parameters and $`V(2)`$ in the two particle spaces is a GOE with unit matrix elements variance; note that in general $`h(1)`$ need not be fixed nor $`V(2)`$ a GOE (in this general case, the ensemble is called EE(1+2); see  for more details). At this stage it is important to stress that all the EGOE(2) results mentioned before are indeed applicable to EGOE(1+2) but only in the domain of chaos. Given $`(m,N)`$ and the average spacing $`\mathrm{\Delta }`$ (generated by $`h(1)`$) of the single particle states (without loss of generality one can put $`\mathrm{\Delta }=1`$) it is possible to find the critical $`\lambda `$ value $`\lambda _c`$ such that for $`\lambda \lambda _c`$ there is onset of chaos (GOE level fluctuations) in many ($`m>>1`$) particle spaces. In fact $`\lambda _c`$ is of the order of the spacing between $`m`$-particle mean-field basis states that are directly coupled by the two-body interaction; see the second and third reference in . For $`\lambda >\lambda _c`$, for example, it is well established that the transition strength sums in EGOE(1+2) follow the EGOE(2) forms; see Fig. 1c ahead. Refs.  give many numerical examples for this, drawn from EGOE(1+2) and more importantly for atomic nuclei in several parts of the periodic table (detailed discussion of the nuclear physics examples is given in the last reference of ). It should be mentioned that the Gausian forms of state and transition strength densities are used in  to derive simple formulas for NPC and $`𝒍_H`$ in transition strength distributions.

For deriving formulas for NPC and $`𝒍_H`$ in wavefunctions, most useful quantity is the strength function (or local density of states) $`F_k(E)`$. Given the mean-field basis states $`|k`$ with energies $`E_k=kHk`$, the eigenstates $`|E`$ can be expanded as $`|E=_kC_k^E|k`$. Then the strength function $`F_k(E)=\delta (HE)^k=_E^{}\left|C_k^E^{}\right|^2\delta (EE^{})`$ and therefore it gives information about the structure of the eigenfunctions. In order to proceed further, let us say that the $`E_k`$ energies are generated by a hamiltonian $`H_k`$ (the structure of $`H_k`$ is discussed ahead). With this, it is easy to identify $`F_k(E)`$ as a conditional density of the bivariate density $`\rho _{biv}(E,E_k)=\delta (HE)\delta (H_kE_k)`$. Taking degeneracies of $`E`$ and $`E_k`$ energies into account we have,

$$\begin{array}{ccc}\hfill \rho _{biv}(E,E_k)& =& \delta (HE)\delta (H_kE_k)\hfill \\ & =& (1/d)\underset{\alpha k,\beta E}{}\left|C_{k,\alpha }^{E,\beta }\right|^2\hfill \\ & =& (1/d)\overline{\left|C_k^E\right|^2}\left[d\rho ^H(E)\right]\left[d\rho ^{H_k}(E_k)\right]\hfill \\ \hfill & & \end{array}$$

$$\begin{array}{c}F_k(E)=\rho _{biv}(E,E_k)/\rho ^{H_k}(E_k)\hfill \\ \overline{\left|C_k^E\right|^2}=\rho _{biv}(E,E_k)/\left[d\rho ^H(E)\rho ^{H_k}(E_k)\right]\hfill \end{array}$$
(1)
In (1) $`d`$ stands for the dimensionality of the $`m`$ particle spaces and $`\overline{\left|C_k^E\right|^2}`$ is the average of $`\left|C_k^E\right|^2`$ over all the degenerate states. Let us now examine the structure of $`H_k`$ and $`\rho _{biv}(E,E_k)`$. Firstly it should be noted that the two-body interaction $`V(2)`$ can be decomposed into two parts $`V(2)=V^{[0]}+𝑽`$ so that $`h(1)+V^{[0]}`$ generates the $`E_k`$ energies (diagonal matrix elements of $`H`$ in the $`m`$-particle mean-field basis states). With $`m`$ particles in $`N`$ single particle states, there is a $`U(N)`$ group and with respect to this group $`V^{[0]}`$ contains a scalar part $`V^{[0],0}`$ (a function of $`m`$), an effective ($`m`$-dependent) one-body (Hartree-Fock like) part $`V^{[0],1}`$ and an irreducible 2-body part $`V^{[0],2}`$. The $`V^{[0],0}+V^{[0],1}`$ will add to $`h(1)`$ giving an effective one-body part of $`H`$; $`h(1)h(1)+V^{[0],0}+V^{[0],1}=𝒉`$. The important point now being that, with respect to a $`U(N)`$ norm, the size of $`V^{[0],2}`$ is usually very small compared to the size of $`𝒉`$ in the $`m`$-particle spaces. With this, $`H=𝒉+𝑽`$ and then the $`H_k`$ is nothing but $`𝒉`$. The piece $`𝑽=V(2)V^{[0]}`$ generates the widths and other shape parameters of $`F_k(E)`$. It should be added that with respect to the $`U(N)`$ norm $`𝒉`$ and $`𝑽`$ are orthogonal and therefore $`\sigma _H^2=\sigma _𝒉^2+\sigma _𝑽^2`$. Definition of $`V^{[0]}`$, a brief discussion of its $`U(N)`$ decomposition etc. are given in Appendix A. For EGOE(1+2), it is well known that the widths of $`F_k(E)`$ are in general constant; see  and Appendix A. The average variance of $`F_k(E)`$’s is given simply by

$$\overline{\sigma _k^2}=\sigma _𝑽^2=(d)^1\underset{\alpha \beta }{}\left|\alpha H\beta \right|^2$$

where $`\alpha `$ and $`\beta `$ are $`m`$-particle mean-field basis states indices. The results, (i) the norm of the $`V^{[0],2}`$ part is negligible and (ii) the widths of the strength functions are nearly constant (with little fluctuations) are well verified in a number of examples; see  and references in  for many nuclear physics examples. EGOE(1+2) discussions in the literature tacitly assume that $`𝒉`$ is $`h(1)`$ and $`𝑽`$ is $`V(2)`$ and the same is assumed from now on, i.e $`H=𝒉+\lambda 𝑽h(1)+\lambda V(2)`$. In addition to (i) and (ii), it is well verified in a number of numerical calculations that: (iii) $`F_k(E)`$’s exhibit a transition from BW to Gaussian form in the chaotic domain defined by $`\lambda >\lambda _{F_k}`$; usually $`\lambda _c<\lambda _{F_k}`$; see  for some analytical understanding of this result. The results (i), (ii) and (iii) clearly imply that the $`\rho _{biv}(E,E_k)`$ is a bivariate Gaussian and this result was first mentioned in . A numerical example for the BW to Gaussian transition in strength functions in EGOE(1+2) is shown in Fig. 1a. In this example $`\lambda _{F_k}0.2`$ and it is much larger than $`\lambda _c=0.06`$ obtained via the results for the Dyson-Mehta $`\overline{\mathrm{\Delta }}_3`$ level statistic shown in Fig. 1b. Thus there is onset of GOE fluctuations much before the $`F_k(E)`$’s start taking Gaussian form, i.e. $`\lambda _{F_k}>\lambda _c`$. Unlike the case with strength functions (also transition strength densities; see ), as mentioned before, strength sums start following the EGOE(2) form (i.e. $`\rho _{n_i:𝒢}(E)/\rho _𝒢(E)`$) from $`\lambda =\lambda _c`$. This is demonstrated in Fig. 1c where occupancies $`En_iE`$ as a function of $`E`$ are shown (they correspond to the strength sums generated by single state ($`i`$) destruction operators). As mentioned in the introduction, the nature of NPC (which is the inverse of IPR) in wavefunctions in the $`\lambda _c\lambda <\lambda _{F_k}`$ domain where $`F_k(E)`$ is of BW form (i.e. in the BW domain) was studied in  while the present article is concerned with the $`\lambda >\lambda _{F_k}`$ domain (i.e. the Gaussian domain) where $`F_k(E)`$ is of Gaussian form.

## III EGOE(1+2) results for NPC and $`𝒍_H`$ in wavefunctions

For EGOE(1+2), in the chaotic domain with $`\lambda >\lambda _{F_k}`$, one has from Sect. II the results: (i) $`E_k`$ are generated by $`H_k=h(1)`$, therefore the variance of $`\rho ^{H_k}(E_k)`$ is $`\sigma _𝒉^2`$; (ii) widths of the strength functions are constant and they are generated by $`V(2)`$, the average variance $`\overline{\sigma _k^2}=\sigma _𝑽^2`$; (iii) $`F_k(E)`$’s are Gaussian in form; (iv) $`F_k(E)`$ is a conditional density of the bivariate Gaussian $`\rho _{biv:𝒢}(E,E_k)`$. The correlation coefficient $`\zeta `$ of $`\rho _{biv:𝒢}(E,E_k)`$ is given by,

$$\zeta =\frac{(Hϵ_H)(H_kϵ_H)}{\sqrt{(Hϵ_H)^2(H_kϵ_H)^2}}=\sqrt{\left(1\frac{\overline{\sigma _k^2}}{\sigma _H^2}\right)}$$
(2)
Note that the centroids of the $`E`$ and $`E_k`$ energies are both given by $`ϵ_H=H`$. In (2) the second equality is obtained by using the orthogonality between $`h(1)`$ and $`V(2)`$ operators. It is immediately seen that the $`\zeta ^2`$ is nothing but the variance of $`E_k`$’s (the centroids of $`F_k(E)`$) normalized by the state density variance. The $`\rho _{biv:𝒢}(E,E_k)`$, which takes into account the fluctuations in the centroids of $`F_k(E)`$ and assumes that the variances are constant, is used to derive formulas for NPC and $`𝒍_H`$ in the wavefunctions (methods of taking into account variance fluctuations will be discussed ahead) $`\psi _E=|E`$ expanded in the mean-field basis defined by the states $`\varphi _k=|k`$. Let us first define NPC and $`𝒍_H`$,

$$\begin{array}{ccc}\hfill |E& =& \underset{k}{}C_k^E|k\hfill \\ \hfill & & \\ \hfill \text{(NPC)}_E& =& \left[\underset{k}{}\left|C_k^E\right|^4\right]^1,\hfill \\ \hfill 𝒍_H(E)& =& exp\left[(\text{S}^{info})_E\right]/(0.48d);\hfill \\ \hfill (\text{S}^{info})_E& =& \underset{k}{}\left|C_k^E\right|^2\mathrm{ln}\left|C_k^E\right|^2.\hfill \end{array}$$
(3)
In (3) $`0.48d`$ is the GOE value for $`\text{S}^{info}`$, thus $`𝒍_H=1`$ for GOE. Similarly NPC is $`d/3`$ for GOE.

In terms of the locally renarmalized amplitudes $`𝒞_k^E=C_k^E/\sqrt{\overline{\left|C_k^E\right|^2}}`$ where the bar denotes ensemble average with respect to EGOE(1+2), $`{\displaystyle \underset{k}{}}\left|C_k^E\right|^4={\displaystyle \underset{k}{}}\left|𝒞_k^E\right|^4\left(\overline{\left|C_k^E\right|^2}\right)^2`$. Then the ensemble averaged (NPC)E is obtained as follows,

$$\begin{array}{ccc}\hfill \overline{\underset{k}{}\left|C_k^E\right|^4}& \stackrel{EGOE(1+2)}{}& \underset{k}{}\overline{\left|𝒞_k^E\right|^4}\left(\overline{\left|C_k^E\right|^2}\right)^2\hfill \\ & =& 3\underset{k}{}\left(\overline{\left|C_k^E\right|^2}\right)^2\hfill \\ & =& \frac{(3/d)}{\left[\rho _𝒢^H(E)\right]^2}𝑑E_k\frac{\left[\rho _{biv:𝒢}(E,E_k)\right]^2}{\rho _𝒢^{H_k}(E_k)}=\frac{(3/d)}{\left[\rho _𝒢^H(E)\right]^2}𝑑E_k\rho _𝒢^{H_k}(E_k)\left[F_{k:𝒢}(E)\right]^2\hfill \\ \hfill & & \\ \hfill (\text{NPC})_E& =& (d/3)\sqrt{1\zeta ^4}exp\left\{\frac{\zeta ^2\widehat{E}^2}{1+\zeta ^2}\right\}\hfill \end{array}$$
(4)
The $`\widehat{E}`$ in (4) is the standardized $`E`$, i.e. it is zero centered and normalized to unit width, $`\widehat{E}=(Eϵ_H)/\sigma _H`$. In the first step in (4) the fact that EGOE exhibits average fluctuations separation (with little communication between the two) is used. For example, in the normal mode decomposition of the EGOE state denity, it is seen that the long wavelength parts generate the smoothed Gaussian density (with corrections) and the short wavelength parts the GOE fluctuations with damping of the intermediate ones (see  for detailed discussions on this important result). This allows one to carry out $`\left|𝒞_k^E\right|^4`$ ensemble average independent of the other smoothed (average) term. In the second line the Porter-Thomas form of local strength fluctuations is used and then $`\overline{\left|𝒞_k^E\right|^4}`$=3, a GOE result. In the third step the result in (1) is used. Then, the Gaussian forms, valid in the chaotic domain ($`\lambda >\lambda _{F_k}`$), of all the densities for EGOE(1+2) give the final formula (this result was quoted first in  without details). Before turning to the formula for the localization length $`𝒍_H`$, let us briefly discuss about the corrections to (4) due to the fluctuations in the variances of $`F_k(E)`$; the form with $`F_k(E)`$ shown explicitly, is written in (4) for this purpose and this form also allows one to understand the results in  as discussed ahead.

The correction to NPC due to $`\delta \sigma _k^2=\sigma _k^2\overline{\sigma _k^2}0`$ is obtained by using, for small $`\left|\delta \sigma _k^2\right|`$, the hermite polynomial expansion which gives , $`F_{k:𝒢}(E)F_{k:𝒢}(E)\left\{1+c_2(_k^21)\right\}`$ where $`c_2=\delta \sigma _k^2/2\overline{\sigma _k^2}`$ and $`_k=(EE_k)/\sqrt{\overline{\sigma _k^2}}`$. This corrected $`F_k(E)`$ is used in the integral form with $`F_k(E)`$ in (4). As NPC involves sum over all the $`|k`$ states, it is a valid assumption to treat $`(\delta \sigma _k^2)`$’s as random and therefore in $`[F_k(E)]^2`$ only the terms that are quadratic in $`(\delta \sigma _k^2)`$ will contribute (see ). Replacing $`\left[(\delta \sigma _k^2)/\overline{\sigma _k^2}\right]`$ by $`(\delta \sigma ^2)/\overline{\sigma _k^2}=\left[(d)^1\{_k(\delta \sigma _k^2)^2\}\right]^{1/2}/\overline{\sigma _k^2}`$ and substituting the corrected $`F_k(E)`$ for $`F_{k:𝒢}(E)`$ in (4), we get

$$\begin{array}{ccc}\hfill (\text{NPC})_E& =& \frac{(3/d)}{\left[\rho _𝒢^H(E)\right]^2}𝑑E_k\rho _𝒢^{H_k}(E_k)\left[F_{k:𝒢}(E)\right]^2\left(1+\frac{(\delta \sigma ^2)}{2\overline{\sigma _k^2}}(_k^21)\right)^2\hfill \\ & =& \frac{d}{3}\sqrt{1\zeta ^4}exp\left\{\frac{\zeta ^2\widehat{E}^2}{1+\zeta ^2}\right\}\left\{1+\frac{1}{4}\left[\frac{(\delta \sigma ^2)}{\sigma _H^2}\right]^2X(E)\right\}^1;\hfill \end{array}$$
(5)
$$X(E)=\frac{1}{(1+\zeta ^2)^4}\left[\widehat{E}^42\frac{(1+\zeta ^2)(12\zeta ^2)}{1\zeta ^2}\widehat{E}^2+\left(\frac{1+\zeta ^2}{1\zeta ^2}\right)^2(1+2\zeta ^4)\right]$$

The $`\delta \sigma ^2`$ corection term in (5) is valid only when the fluctuations in the variances of $`F_k(E)`$’s are small (this is in general always true). For small $`\zeta `$ values, the formula for NPC in (5) reduces to the expression given recently, for EGOE(2), by Kaplan and Papenbrock ; they use ideas related to ’scar theory’. In the EGOE(1+2) hamiltonian $`H=h(1)+\lambda V(2)`$, with $`\lambda \mathrm{}`$ one obtains EGOE(2) and then it is clear from the definition in (2) that in this limit $`\zeta 0`$. More precisely, with $`N>>m>>1`$, $`\zeta ^2\left(\begin{array}{c}N\\ 2\end{array}\right)^1`$ and $`\left[(\delta \sigma ^2)/\sigma _H^2\right]^2\left[\left(\begin{array}{c}m\\ 2\end{array}\right)\left(\begin{array}{c}N\\ 2\end{array}\right)\right]^1`$ for $`\{H\}=\{V(2)\}`$; see appendix A. Therefore for finite $`N`$, the correlation coefficient and the variance corrections are small but non zero and in the large $`N`$ limit they are zero giving the GOE result as pointed out in . As we add the mean-field part to the EGOE(2), $`\zeta `$ increases and at the same time the variance correction decreases; see Appendix A. Thus the formula (5) with the ($`\delta \sigma ^2)`$ term is important only for small $`\zeta `$. Eq. (4) is accurate for reasonably large $`\zeta `$ (say for $`\zeta 0.3`$) as in the examples discussed in . All these results are well tested by the numerical examples in Sect. IV.

Proceeding exactly as in (4), formula for the localization length $`𝒍_H`$ as a function of the excitation energy is derived. Using the definition (3), writing $`\left|C_k^E\right|^2`$ in terms of $`\left|𝒞_k^E\right|^2`$ and $`\overline{\left|C_k^E\right|^2}`$, using the GOE results $`\overline{\left|𝒞_k^E\right|^2}`$=1 and $`\overline{\left|𝒞_k^E\right|^2\mathrm{ln}(\left|𝒞_k^E\right|^2)}=\mathrm{ln}0.48`$, applying the last equality in (1) and replacing all the densities by their corresponding Gaussian forms, converting the sum in (3) into an integral and finally carrying out the integration, the expression for $`𝒍_H`$ in wavefunctions is obtained,

$$\begin{array}{ccc}\hfill 𝒍_H(E)& \stackrel{EGOE(1+2)}{}& 𝑑E_k\frac{\rho _{biv:𝒢}(E,E_k)}{\rho _𝒢^H(E)}\mathrm{ln}\left\{\frac{\rho _{biv:𝒢}(E,E_k)}{\rho _𝒢^{H_k}(E_k)\rho _𝒢^H(E)}\right\}\hfill \\ & & \\ & =& \sqrt{1\zeta ^2}exp\left(\frac{\zeta ^2}{2}\right)exp\left(\frac{\zeta ^2\widehat{E}^2}{2}\right)\hfill \end{array}$$
(6)
The result in (6) was reported in  without details. By rewriting the integral in (6) in terms of $`F_k(E)`$ and making small $`(\delta \sigma ^2)`$ expansion just as in the case of NPC, the formula incorporating corrections due to fluctuations (with respect to $`k`$) in the variances of $`F_k(E)`$ is derived following the arguments that led to (5). Neglecting higher order terms in $`\left[(\delta \sigma ^2)/\sigma _H^2\right]^2`$, the final result is,

$$𝒍_H(E)=\sqrt{1\zeta ^2}exp\left(\frac{\zeta ^2}{2}\right)exp\left(\frac{\zeta ^2\widehat{E}^2}{2}\right)\left(1\frac{1}{8}\left[\frac{(\delta \sigma ^2)}{\sigma _H^2}\right]^2Y(E)\right);$$
(7)
$$Y(E)=\frac{1}{\left(1\zeta ^2\right)^2}\left\{\left(1\zeta ^2\right)^2\left(\widehat{E}^21\right)^2+4\zeta ^2\left(1\zeta ^2\right)\widehat{E}^2+2\zeta ^4\right\}$$

## IV Numerical tests

NPC and $`𝒍_H`$ are calulated for a EGOE(1+2) with 6 particles in 12 single particle states and the results are compared with (4-7) in Fig. 2. In the numerical calculations, the single particle energies $`(i+1/i),i=1,2,\mathrm{},12`$ define $`h(1)`$ (as in  and Fig. 1), in the two-particle space $`V(2)`$ is a GOE (calculations use 25 members) with unit matrix elements variance and the hamiltonian ensemble is $`\{H_{\alpha ,\lambda }\}=\alpha h(1)+\lambda \{V(2)\}`$. The value of $`\lambda =0.2`$ is fixed so that, for $`\alpha 1`$ the level fluctuations are of GOE; i.e. one is in the chaotic domain (see  and Fig. 1). Results for $`\alpha `$ = 0, 0.5 and 1 in Figs. 2a,b clearly demonstrate that the EGOE(1+2) formulas based on the bivariate Gaussian form for $`\rho _{biv}(E,E_k)`$ are excellent. In these examples $`\zeta `$ values are 0.16, 0.59 and 0.82 respectively. The $`(\delta \sigma ^2)`$ correction is seen to be important only for the case with $`\alpha =0`$. In fact the $`\left[(\delta \sigma ^2)/\sigma _H^2\right]^2`$ values for the three cases considered are $`0.121\times 10^1`$, $`0.545\times 10^2`$ and $`0.134\times 10^2`$. Thus, for realistic fermion models that are represented by EGOE(1+2) (with $`\lambda >\lambda _{F_k}`$), the correction due to variance fluctuations is expected to be significant only in the situation $`\zeta `$ is small. Extension of EGOE(2) with explicit inclusion of spin degrees of freedom (each single particle level is taken to be doubly degenerate with $`s_z=\pm \frac{1}{2}`$; see the third reference in ) was considered and for a system of 6 fermions in 7 levels (i.e. $`m=6`$, $`N=7\times 2`$) with total $`S_z=0`$, giving $`d=1225`$, NPC was calculated as a function of the excitation energy in ; we call this model EGOE(2)-S. In this example, as given in , $`\zeta =0.3`$ and $`\left[(\delta \sigma ^2)/\sigma _H^2\right]^2=0.052`$. Thus, here the corrections due to variance fluctuations are non-negligible (the situation in this case is similar to the $`\alpha =0`$ case in Fig. 2) and applying (5) gives excellent description, as shown in Fig. 3a, of the results for NPC reported in  for the EGOE(2)-S model. Returning to Fig. 2, it should be mentioned that there are differences between the numerical results and the predictions based on (4,6) even for the cases with $`\zeta `$ = 0.59 and 0.82. These may be due to the departures of $`\rho _{biv}(E,E_k)`$ from the bivariate Gaussian form. An important observation from (4,6) is, at the spectrum center NPC=$`(d/3)\sqrt{1\zeta ^4}`$ and $`𝒍_H=\sqrt{1\zeta ^2}exp(\zeta ^2/2)`$. Therefore for $`\zeta ^2`$ close to 0.8 or large, there will be large deviations from GOE even at the spectrum center for a system described by EGOE(1+2). This is clearly seen in the $`\alpha =1`$ case in Fig. 2; here $`\zeta =0.82`$. Finally, it should be mentioned that the EGOE(1+2) calculations for the $`N=14`$, $`m=7`$ system (the case considered in Fig. 1) are also carried out and the results are seen to be essentially same as in Figs. 2a,b.

Let us now turn to the nuclear shell model which is a realistic interacting fermion model. There are shell model results for the $`(2s1d)`$ shell (here after called $`sd`$ shell) nuclei 28Si  and 22Na (see  and the second reference in ) for NPC and $`𝒍_H`$ in wavefunctions. For 28Si the 839 dimensional $`J=0,T=0`$ space (with six protons and six neutrons in the $`sd`$ shell) and the 3243 dimensional $`J=2,T=0`$ space are considered. Similarly, for 22Na the 307 dimensional $`J=2,T=0`$ space (with three protons and three neutrons in the $`sd`$ shell) is considered. The results for these nuclei are analyzed using (4,6) as briefly discussed in . In all the $`sd`$ shell examples, $`\zeta 0.60.7`$ and therefore the situation is similar to the $`\alpha =0.5`$ case in Fig. 2. Thus, in these examples the departures from GOE at the spectrum center are no more than 10% but away from the center, there are large departures. The shell model NPC and $`𝒍_H`$ for $`sd`$ shell nuclei are seen to be well described by the EGOE forms in (4,6). For further confirming this, NPC is evaluated for 24Mg in the 325 dimensional $`J=0,T=0`$ space (with 4 protons and 4 neutrons in the $`sd`$ shell) and the results are shown in Fig. 3b; here $`\zeta =0.68`$. It can be concluded that the deviations of the $`sd`$ shell model results from GOE clearly imply that the shell model hamiltonians are well represented by EGOE(1+2) (with $`\lambda >\lambda _{F_k}`$) but not by GOE. It is also seen that the corrections due to $`(\delta \sigma ^2)`$ are small for $`(sd)`$ shell nuclei (note that here $`\zeta `$ is large); in the 24Mg example, $`\left[(\delta \sigma ^2)/\sigma _H^2\right]^2=0.024`$. In order to further substantiate the EGOE description of the structure of nuclear shell model wavefunctions, we have analyzed using (6) the $`𝒍_H(E)`$ vs $`E`$ results reported recently in  for $`2p1f`$ shell (hereafter called $`pf`$ shell) nuclei 50Ca and 46Sc. In the case of 50Ca the 2051 dimensional $`J=6,T=5`$ space (with 10 protons in the $`pf`$ shell) and in 46Sc the 2042 dimensional $`J=1,T=2`$ space (with one proton and 5 neutrons in the $`pf`$ shell) are considered and a modern large shell model code was used for obtaining the $`𝒍_H`$ values. The shell model results for $`𝒍_H`$ in Fig. 9 of , via (2), determine $`\zeta `$ to be 0.96 and 0.92 respectively for the 50Ca and 46Sc examples; results for 46Sc are shown in Fig. 3c. From the definition (2) but employing averages over $`mT`$ spaces (instead of $`mJT`$ spaces), we obtain the $`\zeta `$ values 0.91 and 0.89 respectively. It should be pointed out that given the single particle energies and the two-body matrix elements of the shell model hamiltonians, it is easy to calculate $`\zeta `$ in fixed $`mT`$ spaces using trace propagation methods (by extending (A.3) and (A.4)) . The $`pf`$ shell examples are similar to the $`\alpha =1`$ case in Fig. 2 and therefore, as expected, one sees large departures from GOE even at the spectrum center. Finally, it is seen from the shell model examples in Fig. 3 and the EGOE examples in Fig.2 that further corrections to the results in (4-7) need to be worked out but this is not attempted in this paper. Similarly, study of the nature of fluctuations in NPC and $`𝒍_H`$ is postponed for future.

## V Conclusions

Wavefunction structure given by the EGOE(1+2) random matrix ensemble $`\{H\}=h(1)+\lambda \{V(2)\}`$ is studied by deriving compact formulas for NPC and $`𝒍_H`$. They are based on: (i) the Gaussian form for strength functions $`F_k(E)`$’s and the bivariate Gaussian form for $`\rho _{biv}(E,E_k)`$ (with $`F_k(E)`$ being a conditional density of $`\rho _{biv}(E,E_k)`$) which are valid in the chaotic domain defined by $`\lambda >\lambda _{F_k}`$; (ii) there is average-fluctuations separation (with little communication between the two) in energy levels and strengths with local strength fluctuations following the Porter-Thomas law; (iii) there is a significant unitary group decomposition of the hamiltonian. With EGOE(1+2), the NPC and $`𝒍_H`$ take Gaussian forms as a function of the excitation energy and they are defined by the bivariate correlation coefficient $`\zeta `$ which measures the variance of the distribution of centroids of $`F_k(E)`$’s relative to the state density variance. Theory for incorporating corrections due to fluctuations in the variances (with $`k`$) of $`F_k(E)`$ is also given. For small $`\zeta `$, the present formulation gives back the results for pure EGOE(2) (i.e. in the $`\lambda \mathrm{}`$ limit of EGOE(1+2)) as derived in  recently. The formulas derived for NPC and $`𝒍_H`$ are subjected to numerical EGOE(1+2) tests with $`\zeta `$ changing from 0.1 to 0.8. These and the analysis of the results for a EGOE(2)-S example and some nuclear shell model examples, clearly point out that isolated finite realistic interacting particle systems, in the chaotic domain ($`\lambda \lambda _{F_k}`$), will have wavefunction structure as given by EGOE(1+2). Finally, in the theory given by (4,6), NPC and $`𝒍_H`$ depend on just one parameter and this appears to be an aspect of ’geometric chaos’ (see  for a recent discussion on ’geometric chaos’).

###### Acknowledgements.

This work has been partially supported by DST(India).  Apendix A  Let us consider a system of $`m`$ fermions in $`N`$ single particle states with a (1+2)-body hamiltonian $`H=h(1)+V(2)`$ where $`h(1)`$ is specified by the single particle energies $`ϵ_i`$ (with $`i`$ denoting the $`i^{\text{th}}`$ single particle state, $`h(1)=_iϵ_in_i`$ where $`n_i`$ are number operators) and $`V(2)`$ by the two-body matrix elements $`V_{ijkl}=klV(2)ij`$. The two-body interaction can be seperated into $`V(2)=V^{[0]}+𝑽`$ where $`V^{[0]}`$ is given by
$$V^{[0]}=\underset{i<j}{}V_{ijij}n_in_j(A.1)$$

The $`h(1)+V^{[0]}`$ generates the $`F_k(E)`$ centroids $`E_k`$. With $`N`$ single particle states, there is a $`U(N)`$ group generated by the $`N^2`$ operators $`a_i^{}a_j`$ where $`a_i^{}`$ and $`a_j`$ are one particle creation and destruction operators respectively. With respect to this $`U(N)`$ group, $`V^{[0]}`$ decomposes into $`\nu =0,1,2`$ parts and their explicit structure is (for a given $`m`$),

$$\begin{array}{c}V^{[0],0}=\left(\begin{array}{c}m\\ 2\end{array}\right)\left(\begin{array}{c}N\\ 2\end{array}\right)^1V^0;V^0=\underset{i<j}{}V_{ijij}\hfill \\ V^{[0],1}=\frac{m1}{N2}\underset{i}{}\zeta _in_i;\zeta _i=\underset{ji}{}\left(V_{ijij}V^0\right)\hfill \\ V^{[0],2}=V^{[0]}V^{[0],0}V^{[0],1}(A.2)\hfill \end{array}$$

Similarly the $`h(1)`$ operator will have $`\nu =0,1`$ parts; $`h^0=mϵ^0`$ where $`ϵ^0=(N)^1_iϵ_i`$ and $`h^1(1)={\displaystyle \underset{i}{}}ϵ_i^1n_i`$ where $`ϵ_i^1=ϵ_iϵ^0`$. Finally it is to be noted that $`𝑽`$ behaves as a $`\nu =2`$ operator.

The $`U(N)`$ norm (in the $`m`$-particle spaces) of an operator $`𝒪`$ is defined by $`𝒪_m=\sqrt{(𝒪𝒪^m)^{}(𝒪𝒪^m)^m}`$. An important theorem is that the $`\nu =0,1,2`$ parts of $`H`$ are orthogonal with respect to this $`U(N)`$ norm. For a $`\nu =1`$ operator $`𝒪(1)=_ie_in_i`$, the norm square is simply given by

$$𝒪(1)_m^2=\frac{m(Nm)}{N(N1)}\underset{i}{}e_i^2(A.3)$$

Similarly for a $`\nu =2`$ operator $`𝒪(2)`$,

$$𝒪(2)_m^2=\frac{m(m1)(Nm)(Nm1)}{2(N2)(N3)}𝒪^{}(2)𝒪(2)^2(A.4)$$

Using (A.3) and (A.4) one can calculate the norms of $`h^1+V^{[0],1}`$ and $`V^{[0],2}`$ and in general the later is very small compared to the former. Then $`h(1)+V^{[0]}𝒉=_i\xi _in_i`$ where $`\xi _i=ϵ_i^1+\frac{m1}{N2}\zeta _i`$ (note that at the end we add the spectrum centroid generating part $`h^0+V^{[0],0}`$ to $`𝒉`$). Thus, neglecting the $`V^{[0],2}`$ part, the centroids of $`F_k(E)`$’s are generated by $`𝒉`$ and the variances by $`𝑽`$. As $`𝒉`$ and $`𝑽`$ are orthogonal, $`\sigma _H^2=\sigma _𝒉^2+\sigma _𝑽^2`$. These variances, in $`m`$-particle spaces, follow easily from (A.3,A.4). See  for further details.

Let us consider a EGOE(1+2) hamiltonian $`H=\alpha h(1)+\lambda V(2)`$ with unit spacing between the $`ϵ_i`$’s and the $`V_{ijkl}`$ taken as zero centered Gaussian variables with unit variance. In the $`N>>m>>1`$ situation one can study the behaviour of $`\zeta ^2`$ and $`(\delta \sigma ^2)`$ as follows. The correlation coefficient $`\zeta ^2=\sigma _𝒉^2/\sigma _H^2`$ and, neglecting the contributions of $`V(2)`$ to $`\sigma _𝒉`$, one gets $`\sigma _𝒉^2(mN^2/12)\alpha ^2`$. Similarly $`\sigma _𝑽^2\left(\begin{array}{c}m\\ 2\end{array}\right)\left(\begin{array}{c}N\\ 2\end{array}\right)\lambda ^2`$. Here (A.3) and (A.4) are used. Therefore $`\zeta ^2=[(1+3m(\lambda /\alpha )^2]^1`$ and this expression gives 0.51 and 0.76 for the $`\alpha =0.5`$ and 1 cases in Fig. 2. They compare well with the exact numbers given in Fig. 2. However this estimate fails in the situation $`\alpha 0`$. For $`\alpha =0`$ the $`𝒉`$ has to be replaced by $`V^{[0]}`$ and then the $`E_k`$ energies are a sum of $`\left(\begin{array}{c}m\\ 2\end{array}\right)`$ zero centered Gaussian variables each with variance $`\lambda ^2`$. This together with the $`\sigma _𝑽^2`$ expression, gives $`\zeta ^2\left(\begin{array}{c}N\\ 2\end{array}\right)^1`$ for $`\alpha 0`$ as pointed out in . The number quoted for the $`\alpha =0`$ case in Fig. 2 is close to this estimate. Finally an estimate for $`[(\delta \sigma ^2)/\sigma _H^2]^2`$ is obtained from (A.4) by noting that $`\sigma _𝑽^2`$ is a sum of $`K\left(\begin{array}{c}m\\ 2\end{array}\right)\left(\begin{array}{c}N\\ 2\end{array}\right)`$ $`\chi ^2`$-variables and therefore $`[(\delta \sigma ^2)/\sigma _𝑽^2]^2=2/K`$ as given first in . Then, $`\sigma _𝑽^2=(1\zeta ^2)\sigma _H^2`$ gives the final result $`[(\delta \sigma ^2)/\sigma _H^2]^22(1\zeta ^2)/\left(\begin{array}{c}m\\ 2\end{array}\right)\left(\begin{array}{c}N\\ 2\end{array}\right)`$.

FIGURE CAPTIONS

Fig. 1 Strength functions $`F_k(E)`$, Dyson-Mehta $`\overline{\mathrm{\Delta }}_3`$ statistic for level fluctuations and occupancies $`En_iE`$ for EGOE(1+2) for various values of the interaction strength $`\lambda `$ in $`\{H\}=h(1)+\lambda \{V(2)\}`$ for a system of 7 fermions in 14 single particle states (due to computational constraints, here only one member is considered just as in ); the matrix dimension is 3432. The single particle energies used in the calculations are $`ϵ_i=(i+1/i),i=1,2,\mathrm{},14`$ just as in . (a) the histograms are EGOE(1+2) results for the strength functions, continuous curves are BW fit and the dotted curves are Gaussian for $`\lambda 0.1`$ and Edgeworth corrected Gaussian  for $`\lambda >0.1`$. In constructing the strength functions, $`|C_k^E|^2`$ are summed over the basis states $`|k`$ in the energy window $`\widehat{E}_k\pm \mathrm{\Delta }`$ and then the ensemble averaged $`F_{\widehat{E}_k}(\widehat{E})`$ vs $`\widehat{E}`$ is constructed as a histogram; the value of $`\mathrm{\Delta }`$ is chosen to be 0.025 for $`\lambda 0.1`$ and beyond this $`\mathrm{\Delta }=0.1`$. Here $`\widehat{E}_k=(E_kϵ_H)/\sigma _H`$ and in the figure $`\widehat{E}_k=0`$. Note that for $`\lambda _{F_k}0.2`$ there is BW to Gaussian transition. (b) The $`\overline{\mathrm{\Delta }}_3(L)`$ statistic for overlapping intervals of length $`L40`$ are compared with Poisson and GOE values. For $`\lambda 0.06`$, there is Poisson to GOE transition in the $`\overline{\mathrm{\Delta }}_3`$ statistic. (c) The wavy curves are numerical EGOE(1+2) results for occupancies and the smoothed curves with $`\lambda 0.06`$ correspond to the results of EGOE(2) theory (ratio of Gaussians). Note that for $`\lambda <0.06`$ there are wide fluctuations in occupancies and the smoothed forms here are meaningless. All the results are shown for the lowest 6 single particle states. Results similar to those in the figure, for the $`N=12,m=6`$ case are reported in .

Fig. 2 (a) Number of principal components NPC and (b) the localization length $`𝒍_H`$ in wavefunctions for a system of 6 interacting particles in 12 single particle states (matrix dimension is 924). Here, for conveniance, the EGOE(1+2) hamiltonian is changed to $`\{H_{(\alpha ,\lambda )}\}=\alpha h(1)+\lambda \{V(2)\}`$. Numerical EGOE(1+2) results correspond to filled circles. The continuous curves correspond to the theory (4) for NPC and (6) for $`𝒍_H`$. For the case with $`\alpha =0`$, the dashed curves correspond to the theory (5) for NPC and (7) for $`𝒍_H`$. For other cases the correction due to variance fluctuations is negligible and hence only the results of (4,6) are shown in the figure. Note that NPC=$`d/3`$ and $`𝒍_H=1`$ for GOE. See text for further details.

Fig. 3 (a) NPC for the EGOE(2)-S model described in the text compared with the results given by (4,5). The filled circles are for the numerical EGOE(2)-S calculations reported in . The continuous and dashed curves represent (4) and (5) respectively. (b) NPC for the $`(sd)`$ shell nucleus 24Mg compared with (4). The shell model calculations are same as in . In this example (4) and (5) give almost identical results and hence the curve corresponding to (5) is not shown in the figure. (c) Shell model results for 46Sc for $`𝒍_H`$ reported in  compared with the theoretical curve given by (6) with $`\zeta =0.92`$.