format: text

text:
# The exotic double-charge exchange 𝜇⁻→𝑒⁺ conversion in nuclei

## 1 Introduction

Modern gauge theories (grand unified theories, supersymmetry, etc.), which go beyond the standard model (SM), predict a great number of processes which violate the lepton and/or lepton family (flavor) quantum numbers -. One can distinguish three categories of such processes: the purely leptonic ($`\mu e+\gamma `$, $`\mu 3e`$, etc.), the semileptonic hadron decays ($`K_L^0\mu ^\pm e^\pm `$, etc.), and the semileptonic which take place in the presence of nuclei. Among the most interesting such processes are those which take place in a muonic atom -. One exotic such process is the muon-to-positron conversion,

$$\mu ^{}+(A,Z)e^++(A,Z2),$$
(1)
which violates the muonic ($`L_\mu `$), electronic ($`L_e`$) and total lepton ($`L`$) quantum numbers -. Another exotic process is the muon-to-electron conversion,

$$\mu ^{}+(A,Z)e^{}+(A,Z),$$
(2)
which violates only the lepton-family quantum numbers $`L_\mu `$ and $`L_e`$ .

Both of these processes can experimentally be studied simultaneously, since both of them have the same intrinsic background and the same initial state (a muon at rest in the innermost 1S orbit of a muonic atom). In the present work we will focus our attention on reaction (1).

In recent years, continuous experimental efforts have been devoted for the measurement of the branching ratio $`R_{\mu e^+}`$ defined as the ratio of the $`(\mu ^{},e^+`$) conversion rate divided by the total rate of the ordinary muon capture reaction , i.e.

$$R_{\mu e^+}=\mathrm{\Gamma }(\mu ^{}e^+)/\mathrm{\Gamma }(\mu ^{}\nu _\mu ).$$
(3)
Up to now only upper limits have been set and the best limit found is for the 48Ti nucleus at TRIUMF and PSI  yielding the values

$$R_{\mu e^+}4.6\times 10^{12}\text{[15]}$$

and

$$R_{\mu e^+}4.4\times 10^{12}\text{[13]},$$

respectively. This limit could be further improved by future experiments, at PSI (SINDRUM II experiment), which aims to push the sensitivity of the branching ratio $`R_{\mu e^+}`$ to $`10^{14}`$, and at Brookhaven (MECO experiment)  with expected sensitivity about three to four orders of magnitude below the existing experimental limits .

Traditionally the exotic $`\mu ^{}e^\pm `$ processes were searched by employing medium heavy (like 48Ti and 63Cu)  or very heavy (like 208Pb and 197Au)  targets. For technical reasons the MECO target has been chosen to be the light nucleus 27Al . The MECO experiment, which is planned to start soon at the Alternating Gradient Synchrotron (AGS), is going to employ a new very intense $`\mu ^{}`$ beam and a new detector . The basic feature of this experiment is the use of a pulsed $`\mu ^{}`$ beam to significantly reduce the prompt background from $`\pi ^{}`$ and $`e^{}`$ contaminations.

The best upper limit for the $`\mu ^{}e^{}`$ conversion branching ratio $`R_{\mu e^{}}`$ set up to the present has been extracted at PSI (SINDRUM II experiment) as

$$R_{\mu e^{}}^{Au}5.0\times 10^{13},\text{ for }^{197}\text{Au target }\text{[16]}.$$
(4)
For the 48Ti target the determined best limit is

$$R_{\mu e^{}}^{Ti}6.1\times 10^{13}\text{[13]},$$
(5)
while for the 208Pb target the extracted limit is

$$R_{\mu e^{}}^{Pb}4.6\times 10^{11}\text{[14]}.$$
(6)
Processes (1) and (2) are very good examples of the interplay between nuclear and particle physics in the area of physics beyond the standard model. Moreover, the ($`\mu ^{},e^+`$) conversion has many similarities with the neutrinoless double $`\beta `$-decay ($`0\nu \beta \beta `$),

$$(A,Z)e^{}+e^{}+(A,Z+2),$$
(7)
and especially with its sister electron to positron conversion , which violate the lepton-flavor ($`L_e`$) and total lepton ($`L`$) quantum numbers. Both reactions (1) and (7) involve a change of charge by two units and thus they cannot occur in the same nucleon. Both of them are forbidden, if lepton number is absolutely conserved. One can show that, if either of these processes is observed, the neutrinos must be massive Majorana particles. In spite of the many similarities, however, these double charge exchange processes have some significant differences which can be briefly summarized as follows:

(i) Due to the nuclear masses involved, neutrinoless double beta decay can occur only in specific nuclear systems for which single beta decay is absolutely forbidden, due to energy conservation, or greatly hindered due to angular momentum mismatch. These systems, with the possible exception of $`{}_{}{}^{48}Ca`$, have complicated nuclear structure. Furthermore, the neutrinoless double beta decay can lead mainly to the ground state and, in some cases, to few low-lying excited states of the final nucleus. Such constraints are not imposed on process (1), due to the rest energy of the disappearing muon.

(ii) From the corresponding experiments, in conjunction with appropriate nuclear matrix elements as input, one may extract lepton violating parameters, which depend on flavor. Thus, in the framework of the neutrino mixing models, the amplitudes of neutrinoless $`\beta \beta `$ decay and ($`\mu ^{},e^+`$), for leptonic currents of the same chirality, are proportional to some combination of the neutrino masses, which are different from each other. The same is true in the case of the mass independent lepton violating parameters $`\eta `$ and $`\lambda `$ appearing in the amplitude, if the leptonic currents are of opposite chirality. One does not know a priori which flavor combination is favored.

(iii) The long wavelength approximation does not hold in the case of ($`\mu ^{},e^+`$) conversion, since the momentum of the outgoing $`e^+`$ is quite high. Thus, the effective two-body operator responsible for the ($`\mu ^{},e^+`$) conversion is strongly energy dependent and more complicated than the corresponding one for the $`0\nu \beta \beta `$ decay. On the other hand, one can in this case choose a target, consistent with the standard experimental requirements, so that the nuclear structure required is the simplest possible.

(iv) Neutrinoless double beta decay has the experimental advantage that there are no other fast competent channels. The only other channel is the two-neutrino double-beta decay, which, however, is also of second order in weak interaction and kinematically suppressed.

In any case we view the two processes, ($`\mu ^{},e^+`$) conversion and $`0\nu \beta \beta `$ decay, as providing useful complementary information. As is well known in the case of neutrinos one needs the following information: a) the mixing matrix, b) the three eigenmasses and c) the two independent Majorana phases. These phases are present even in a CP conserving theory since they can take values $`0`$ and $`\pi `$. Study of CP violation in the leptonic sector is beyond the goals of the present paper. The interested reader is referred to the literature, e.g. see .

The neutrino oscillation data can, in principle, determine the mixing matrix and the two independent mass-squared differences (e.g. $`\delta m_{21}^2=m_2^2m_1^2`$, $`\delta m_{31}^2=m_3^2m_1^2`$). This has been done in a number of papers, see e.g. the recent work by Haug et al.  and references therein. They cannot determine, however, the scale of the masses, e.g. the lowest eigenvalue $`m_1`$ and the two relative Majorana phases. Once $`m_1`$ is known, one can find $`m_2=[\delta m_{21}^2+m_1^2]^{1/2}`$ and $`m_3=[\delta m_{31}^2+m_1^2]^{1/2}`$. The masses can be chosen positive by a proper redefinition of the neutrino fields, i.e. by absorbing the sign into the Majorana phases. If the mixing matrix is known, the mass $`m_1`$ can be determined from triton decay as follows:

$$(m_\nu )_{ee}m_\nu =|\underset{j=1}{\overset{3}{}}U_{ej}^{}U_{ej}m_j^2|^{1/2},$$
(8)
since the Majorana phases do not appear (this experiment cannot differentiate between Dirac and Majorana neutrinos). In the above expressions $`U`$ is the mixing matrix.

The two lepton violating processes provide two independent linear combinations of the masses and the Majorana phases. In fact

$$m_\nu _{ee}m_\nu =|\underset{j=1}{\overset{3}{}}U_{ej}U_{ej}e^{i\lambda _j}m_j|$$
(9)
and

$$m_\nu _{\mu e+}=|\underset{j=1}{\overset{3}{}}U_{\mu j}^{}U_{ej}^{}e^{i\lambda _j}m_j|.$$
(10)
It is clear from the above expressions that only the two relative CP phases are measurable.

So, the three experiments together can specify all parameters not settled by the neutrino oscillation experiments. So, if possible, all three of them should be pursued.

Strictly speaking $`(\mu ^{},e^+)`$ and neutrinoless double beta decay should be treated as two-step processes by explicitly constructing the states of the intermediate $`(A,Z\pm 1)`$ systems. It has been found , however, that, for neutrinoless double beta decay, since the energy denominators are dominated by the momentum of the virtual neutrino, closure approximation with some average energy denominator works very well. We expect this approximation to also hold in the case of $`(\mu ^{},e^+)`$ conversion to sufficient accuracy. We will, therefore, replace the intermediate nuclear energies by some suitable average. By summing over the partial rates of all allowed final states of the nucleus $`(A,Z2)`$ we obtain the total rate. This will then be compared to that obtained by invoking closure  with some appropriate mean energy $`E_f`$ of the final states.

So far, theoretically the ($`\mu ^{},e^+`$) process has been investigated  on the exclusive reactions

$$^{40}Ca+\mu ^{}e^++^{40}Ar(gs)$$
(11)
$$^{58}Ni+\mu ^{}e^++^{58}Fe(gs).$$
(12)
In these studies the partial $`g.s.g.s.`$ transition rate was calculated by performing microscopic calculations of these nuclear matrix elements. On the other hand, the total transition strength to all final states (inclusive process) was estimated along the lines of closure approximation and ignoring $`4`$-body terms .

In the present article we apply the shell-model techniques in the investigation of the ($`\mu ^{},e^+`$) conversion reaction

$$^{27}Al+\mu ^{}e^++^{27}Na.$$
(13)
To this end, we construct all the needed wave functions and calculate the rate not only to the ground state but also to the excited states of the final nucleus. As a first step we limited ourselves to the mass mechanism for light neutrinos. So, only the multipolarities associated with the Fermi and Gamow-Teller type operators were relevant. For practical reasons we had to limit ourselves to states lying below some excitation energy ($`25`$ MeV) of the final nucleus 27Na. The distribution of the strengths for the most important multipolarities versus the excitation energy of the daughter nucleus in the reaction (13) is also studied.

The paper is organized as follows. In sect. 2, an extensive presentation of the relevant expressions occurring in the formal description of the $`\mu ^{}e^+`$ transition operators is given. In sect. 3, we deal with the expressions of the branching ratios. In sect. 4, we discuss the evaluation of the inclusive $`\mu ^{}e^+`$ matrix elements by means of explicit construction of the needed nuclear wave functions in the framework of the $`sd`$ shell model. In sect. 5, the results obtained for the transition matrix elements in the case of $`{}_{}{}^{27}Al(\mu ^{},e^+)^{27}Na`$ are presented and discussed. Also the spreading of the contributions due to the occurrence of various multipoles are described. Our conclusions are summarized in sect. 6.

## 2 Brief theoretical formulation of the $`\mu ^{}e^+`$ conversion operators

### 2.1 Effective $`\mu ^{}e^+`$ conversion Lagrangian in Gauge models

¿From the particle physics point of view, processes like $`\mu ^{}e^\pm `$ conversions are forbidden in the SM by total-lepton and/or lepton-flavor (muonic and electronic) quantum number conservation. They have been recognized long time ago as important probes for studying the lepton and lepton-flavor changing charged-current interactions -, which go beyond the standard model.

There are several possible elementary particle mechanisms which can mediate the lepton-violating process (1). The mechanisms which have been studied theoretically are:

(i) Those mediated by a massive Majorana neutrino. In this case we have two possibilities.

1) The chiralities of the two leptonic currents are the same. Then, the amplitude in the case of light neutrinos is proportional to some average neutrino mass or to some average of the inverse of the neutrino mass, if it is heavy.

2) The chiralities of the leptonic currents are opposite. Then, the amplitude is not explicitly dependent on the neutrino mass. It vanishes, however, if the neutrinos are not Majorana particles. This mechanism is significant only in the case of light neutrinos.

(ii) Those accompanied by massless or light physical Higgs particles (majorons).

(iii) Those involving more exotic intermediate Higgs particles.

(iv) Those mediated by intermediate supersymmetric (SUSY) particles.

From a nuclear physics point of view one has to be a bit careful when the intermediate particles are very heavy. The elementary amplitude is constructed at the quark level, but the calculation is performed at the nucleon level. If, in going from the quark to the nucleon level, the nucleons are treated as point-like particles, the nuclear matrix elements are suppressed due to the presence of short range correlations. To avoid this suppression, a cure has been proposed  which treats the nucleons as composite particles described by a suitable form factor. A different approach is to consider mechanisms, which involve particles other than nucleons in the nuclear soup. Such are, e.g., mechanisms whereby the processes (1) and (7) are mediated by the decay of the doubly-charged $`\mathrm{\Delta }^{++}`$ $`(3/2,3/2)`$ resonance  present in the nuclear medium, i.e.

$$\mu ^{}+\mathrm{\Delta }^{++}n+e^+$$
(14)
This process, however, to leading order does not contribute to $`0^+0^+`$ $`\beta \beta `$ decays, but it may contribute to ($`\mu ^{},e^+`$). The other process is induced by pions in flight between the two nucleons  according to the following modes:

(a) The 1-pion mode represented by the reactions:

$$\mu ^{}+pn+\pi ^{}+e^+,\pi ^{}+pn$$
(15)
where the protons and neutrons are bound in the nucleus.

(b) The 2-pion mode represented by the reactions

$$pn+\pi ^+,\pi ^+\mu ^{}\pi ^{}+e^+,\pi ^{}+pn$$
(16)
### 2.2 The transition operators at nuclear level

The gauge models mentioned in the previous subsection give rise to a plethora of effective transition operators. Their essential isospin, spin and radial structure is given as follows.

1) The isospin structure is quite simple, i.e. of the form $`\tau _{}(i)\tau _{}(j)`$ where i and j label the nucleons participating in the process.

2) The spin structure is given in terms of the operators

$$W_{S1}(ij)=1(Fermi)$$
(17)
$$W_{S2}(ij)=\sigma _i\sigma _j(GamowTeller)$$
(18)
$$W_{S3}(ij)=3(\sigma _i\widehat{r}_{ij})(\sigma _j\widehat{r}_{ij})\sigma _i\sigma _j(Tensor)$$
(19)
$$W_{A1}(ij)=ı\sigma _i\times \sigma _j$$
(20)
$$W_{A2}(ij)=\sigma _i\sigma _j$$
(21)
The operator $`\widehat{r}_{ij}`$ is determined below.

3) The orbital part can be expressed in terms of the quantities:

(a) The momentum of the emitted positron ($`p_e`$) obtained from the kinematics of the reaction (1). One finds that

$$p_e|𝐩_e|=m_\mu ϵ_b+QE_x$$
(22)
where $`Q=M(A,Z)M(A,Z2)`$ is the atomic mass difference between the initial, $`(A,Z)`$ and final, $`(A,Z2)`$ nucleus, $`ϵ_b`$ is the binding energy of the muon at the muonic atom ($`ϵ_b0.5`$ MeV for 27Al), $`E_x`$ is the excitation energy ($`E_x=E_fE_{gs}`$) of the final nucleus and $`m_\mu `$ is the muon mass ($`m_\mu =105.6`$ MeV).

(b) The relative ($`𝐫_{ij}`$) and center of mass ($`𝐑_{ij}`$) coordinates which are written as

$$𝐫_{ij}=𝐫_i𝐫_j,\widehat{r}_{ij}=\frac{𝐫_{ij}}{|𝐫_{ij}|},r_{ij}=|𝐫_{ij}|,$$

and

$$𝐑_{ij}=\frac{1}{2}(𝐫_i+𝐫_j),\widehat{R}_{ij}=\frac{𝐑_{ij}}{|𝐑_{ij}|},R_{ij}=|𝐑_{ij}|,$$

The radial part of the operator contains:

(i) the spherical Bessel functions $`j_l(p_er_{ij}/2)`$ and $`j_{}(p_eR_{ij})`$ resulting from the decomposition of the outgoing positron and

(ii) a function $`f(r)`$ of the relative coordinate ($`r=r_{ij}`$) given by:

$$f(r)=\frac{R_0}{r}F(r)\mathrm{\Psi }_{cor}(r),$$
(23)
where the constant $`R_0`$ represents the nuclear radius ($`R_0=1.2A^{1/3}`$). The function $`\mathrm{\Psi }_{cor}(r)`$ is some reasonable two-nucleon correlation function, e.g. of the type

$$\mathrm{\Psi }_{cor}(r)=1e^{ar^2}(1br^2)$$
(24)
with $`a=1.1fm^2`$ and $`b=0.68fm^2`$.

We mention that, strictly speaking in the above radial function the muon wave function must also be included. In fact, if muon is considered as relativistic particle an additional lepton-spin dependence appears in the transition operator. In the case of the light nucleus 27Al studied in the present work, however, the muon is in the $`1s`$ atomic orbit. Its wave function varies very slowly inside the nucleus and thus it can be replaced by its average value. Anyway this average value drops out if the same approximation is assumed for the $`\mu ^{}`$-capture, i.e. for the denominator of the branching ratio $`R_{\mu e^+}`$.

As we have already indicated, the radial function $`F(r)`$ depends on the specific mechanism assumed for the $`\mu ^{}e^+`$ conversion process to occur. The following cases are of interest.

(i) In the case of light Majorana neutrinos, when the leptonic currents are left-handed, $`F(r)`$ takes the form

$$F(r)=\frac{2}{\pi }_0^{\mathrm{}}\frac{sinx}{x\alpha +ıϵ}𝑑x+\frac{2}{\pi }_0^{\mathrm{}}\frac{sinx}{x+\delta }_e𝑑x$$
(25)
The quantities $`\delta _e`$ and $`\alpha `$ are given in terms of the nuclear masses and the average excitation energy of the intermediate states, $`E_{xn}`$, as

$$\delta _e=[E_{xn}+M(A,Z1)M(A,Z)+p_e]r$$

$$\alpha =[m_\mu +M(A,Z)M(A,Z1)E_{xn}]r$$

Note that $`\delta _e`$ depends on the positron momentum. The first term of $`F(r)`$ in Eq. (25) can be written as

$$\frac{2}{\pi }_0^{\mathrm{}}\frac{sinx}{x\alpha +ıϵ}𝑑x=\frac{2}{\pi }P_0^{\mathrm{}}\frac{sinx}{x\alpha }𝑑xı2sin\alpha $$
(26)
As can be seen, the amplitude has now an imaginary part, a fact that was missed in earlier calculations . The principal value integral can be written in an equivalent form as follows:

$$\frac{2}{\pi }P_0^{\mathrm{}}\frac{sinx}{x\alpha }𝑑x=2cos\alpha 1+\frac{2}{\pi }\alpha _0^{\mathrm{}}\frac{sinx}{x(x+\alpha )}𝑑x$$
(27)
The latter expression is more convenient for numerical integration techniques. By combining Eqs. (25)-(27) we can write $`F(r)`$ with separate real and imaginary parts as

$$F(r)=2cos\alpha 1+\frac{2}{\pi }_0^{\mathrm{}}sinx\left[\frac{\alpha }{x(x+\alpha )}+\frac{1}{x+\delta }_e\right]𝑑xı2sin\alpha .$$
(28)
It is worth remarking that, in the case of $`0\nu \beta \beta `$ decay $`\alpha 0`$ therefore $`F(r)=1`$. This simplifies quite well the calculations in the $`0\nu \beta \beta `$ decay process.

(ii) In the case of light Majorana neutrinos, when the leptonic currents are of opposite chirality we have $`F(r)F^{^{}}(r)=r(d/dr)F(r)`$. The same situation occurs in the context of R-parity violating supersymmetric interactions mediated by light Majorana neutrinos in addition to other SUSY particles.

(iii) For heavy intermediate particles, e.g. heavy Majorana neutrinos, we will examine two modes:

1) Only nucleons are present in the nucleus. Then the function $`F(r)`$ reads

$$F(r)=\frac{1}{48}\frac{m_A^2}{m_em_p}x_A(x_A^2+3x_A+3)e^{x_A},x_A=m_Ar$$
(29)
with $`m_e`$, $`m_p`$ the masses of electron and proton respectively. It should be mentioned that the above expression was obtained for neutrino masses much heavier than the proton mass provided that the nucleon is not point like. The particular expression holds, if the nucleon is assumed to have a finite size adequately described  by a dipole shape form factor with characteristic mass $`m_A`$ taking the value $`m_A=0.85GeV/c^2`$

2) The process is mediated by pions in flight between the two interacting nucleons. Then one distinguishes two possibilities :

(a) The 1-pion mode Eq. (15). In this mode $`F(r)`$ is replaced by $`F_{1\pi }^i`$, $`i=GT`$ (Gamow-Teller), $`T`$ (Tensor) where

$$F_{1\pi }^{GT}(x)=\alpha _{1\pi }e^x,F_{1\pi }^T(x)=\alpha _{1\pi }(x^2+3x+3)e^x/x^2$$
(30)
with $`xx_\pi =m_\pi r`$ ($`m_\pi `$ denotes the pion mass) and $`\alpha _{1\pi }=1.4\times 10^2`$. In this case the radial functions are the same with those entering the neutrinoless double beta decay.

(b) The 2-pion mode Eq. (16). Now the radial functions are very different from those appearing in $`0\nu `$ $`\beta \beta `$ decay, since the momentum carried away by the outgoing lepton is not negligible. It is, however, obtained from those entering the neutrinoless double beta decay, via the substitution:

$$F(r)j_l(x_e/2)_0^1j_l((\xi 1/2)x_e)F_{2\pi }^i([\xi (1\xi )x_e^2+x_\pi ^2]^{1/2})𝑑\xi $$
(31)
($`x_e=m_er`$) where $`F_{2\pi }^i,i=GT,T`$ are given by

$$F_{2\pi }^{GT}(x)=\alpha _{2\pi }(x2)e^x,F_{2\pi }^T(x)=\alpha _{2\pi }(x+1)e^x$$
(32)
with $`\alpha _{2\pi }=2.0\times 10^2`$.

### 2.3 Irreducible Tensor Operators

In this section we are going to exhibit the structure of the various irreducible tensor operators relevant to our calculations. We characterize them by the set of quantum numbers $`l,,\lambda ,\mu ,\mathrm{\Lambda },L,S,J`$, some of which may be redundant in some special cases. Thus $`l`$ and $``$ refer to the multipolarity of the outgoing lepton in the relative and center of mass systems. $`\lambda `$ is the orbital rank of the operator in the relative coordinates (the corresponding rank in the CM system is $`\mu `$). $`L`$ is the total orbital rank, $`S`$ is the spin rank and $`J`$ the total angular momentum rank. Finally $`\mathrm{\Lambda }`$ is the rank of the spherical harmonic describing the momentum of the outgoing lepton. The latter couples to zero or 1 with the J-rank, i.e. ($`\mathrm{\Lambda },J)k`$, $`k=0,1`$. Some details on how these operators are combined to give the nuclear matrix elements will be discussed in the Appendix.

We will begin with operators appearing when the chiralities of the two leptonic currents involved are the same. This covers the case of minimal left-handed extensions of the standard model.

One encounters Fermi-type operators, $`\mathrm{\Omega }_F`$, of the form

$`\mathrm{\Omega }_F=`$ $`{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$ (33)
$`\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^\mu (\widehat{R}_{ij})\right]^J,`$
$`\lambda =l,\mu =,S=0,J=L`$
($`\mathrm{\Lambda }`$ is redundant). The Gamow-Teller operators, $`\mathrm{\Omega }_{GT}`$, are similarly written as

$`\mathrm{\Omega }_{GT}=`$ $`{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$ (34)
$`\left[\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^\mu (\widehat{R}_{ij})\right]^L(\sqrt{3})\left[\sigma _i\sigma _j\right]^0\right]^J,`$
$`\lambda =l,\mu =,S=0,J=L`$
($`\mathrm{\Lambda }`$ is redundant). Note that

$$\sigma _i\sigma _j=\sqrt{3}\left[\sigma _i\sigma _j\right]_0^0$$
(35)
The first spin antisymmetric operator is

$`\mathrm{\Omega }_{A1}=`$ $`{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$ (36)
$`\left[\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^\mu (\widehat{R}_{ij})\right]^L(\sqrt{2})\left[\sigma _i\sigma _j\right]^1\right]^J,`$
$`\lambda =l,\mu =,S=1,J=L,|L\pm 1|`$
($`\mathrm{\Lambda }`$ is redundant). Note that

$$ı\sigma _i\times \sigma _j=(\sqrt{2})\left[\sigma _i\sigma _j\right]^1$$
(37)
The second spin antisymmetric operator is

$`\mathrm{\Omega }_{A2}=`$ $`{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$ (38)
$`\left[\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^\mu (\widehat{R}_{ij})\right]^L\left(\sigma _i\sigma _j\right)\right]^J,`$
$`\lambda =l,\mu =,S=1,J=L,|L\pm 1|`$
($`\mathrm{\Lambda }`$ is redundant). Note that each operator must be overall symmetric with respect to interchange of the particle indices. So, in those cases in which the spin operator is of rank unity, $`l`$ must be odd. In the special case of $`0^+0^+`$ neutrinoless double beta decay, only the Fermi and Gamow-Teller operators occur.

We are now going to consider the case in which the theory contains both $`R`$ (Right) and $`L`$ (Left) currents. If both currents are right-handed the above results hold, but the relevant neutrinos are heavy. We thus need only consider the case in which we have $`LR`$ interference in the leptonic sector. This may be important in the case of light neutrinos. As we have already mentioned, this also occurs in the context of R-parity violating supersymmetric interactions, which, in addition to other SUSY particles, involve intermediate light Majorana neutrinos. The amplitude now is proportional to the 4-momentum of the intermediate neutrino. The time component has a structure similar to the above with a different energy dependence. The corresponding operators are indicated by putting a “prime” over the corresponding ones for the mass term. This point will not be further pursued, since their contribution is suppressed. Its space component, after the Fourier transform, gives an amplitude proportional to the gradient of the Fourier transform of the previous case. We thus get the above operators, to be denoted by $`\mathrm{\Omega }_F^{^{}}`$ $`\mathrm{\Omega }_{GT}^{^{}}`$, and $`\mathrm{\Omega }_{A2}^{^{}}`$ (associated with the term linear in the spin), with $`F(r)`$ replaced by $`F^{^{}}(r)`$. Now the overall operator has an overall rank of a vector, obtained by the coupling of the two operator ranks $`J`$ and $`\mathrm{\Lambda }`$. In this case, in addition to operators of the above form, we encounter an operator of spin rank two, which is of the form:

$`\mathrm{\Omega }_T^{^{}}=`$ $`{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f^{^{}}(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$ (39)
$`\left[\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^\mu (\widehat{R}_{ij})\right]^L\left[\sigma _i\sigma _j\right]^2\right]^J,`$
$`\lambda =|l\pm 1|,\mu =,S=2,J=L,|L\pm 1|,|L\pm 2|`$
($`\mathrm{\Lambda }`$ is redundant).

As it has already been mentioned, in the case of heavy intermediate particles one may have to consider pions in flight between nucleons. Then one encounters only Gamow-Teller and tensor operators except that now the radial part is different (see Eq. (30)-(32)).

In the special case of $`0^+0^+`$ neutrinoless double beta decay mediated by light neutrinos one can invoke the long wavelength approximation. Thus, to leading order one finds (up to normalization constants and possibly factors of $`p_e`$) the familiar operators:

$$\mathrm{\Omega }_F=\underset{ij}{}\tau _{}(i)\tau _{}(j)f(r_{ij})(Fermi)$$
(40)
$$\mathrm{\Omega }_{GT}=\underset{ij}{}\tau _{}(i)\tau _{}(j)f(r_{ij})\sigma _i\sigma _j(GamowTeller)$$
(41)
$$\mathrm{\Omega }_{A2}^{^{}}=\underset{ij}{}\tau _{}(i)\tau _{}(j)f^{^{}}(r_{ij})(\sigma _i\sigma _j).(ı\widehat{r}_{ij}\times \widehat{R}_{ij})$$
(42)
$$\mathrm{\Omega }_T^{}=\underset{ij}{}\tau _{}(i)\tau _{}(j)f^{^{}}(r_{ij})\left[3(\sigma _i\widehat{r}_{ij})(\sigma _j\widehat{r}_{ij})\sigma _i\sigma _j\right](Tensor)$$
(43)
## 3 Branching ratio

The branching ratio $`R_{\mu e^+}`$ of the ($`\mu ^{},e^+`$) reaction defined in Eq. (3) contains the LFV-parameters of the specific gauge model assumed. These parameters, are entered in $`R_{\mu e^+}`$ via a single lepton-violating parameter $`n_{eff}`$. Under some reasonable assumptions these parameters can be separated from the nuclear physics aspects of the problem. As has been pointed out , the branching ratio $`R_{\mu e^+}`$ takes the form

$$R_{\mu e^+}=\rho |\eta _{eff}|^2\frac{1}{A^{2/3}Zf_{PR}(A,Z)}\left(\frac{(p_e)_{max}}{m_\mu }\right)^2\underset{f}{}\left(\frac{p_e}{(p_e)_{max}}\right)^2|_{if}|^2$$
(44)
with $`\eta _{eff}=m_\nu _{\mu e+}/m_e`$ The parameters $`\rho `$ and $`\eta _{eff}`$, which depend on the gauge model adopted and the mechanism prevailing, are expected to be very small due to the fact that $`\mu ^{}e^+`$ conversion is a lepton violating second order weak process . In this definition, the total muon capture rate has been written in terms of the well-known Primakoff function $`f_{PR}(A,Z)`$  which takes into account the effect of the nucleon-nucleon correlations on the total muon capture rate. $`|_{if}|^2`$ denotes the square of the partial transition nuclear matrix element between an initial $`|i`$ and a final $`|f`$ state. This can be written as

$$|_{if}|^2=\frac{1}{2J_i+1}\underset{M_fM_i}{}|J_fM_f|\mathrm{\Omega }|J_iM_i|^2$$
(45)
In our case $`|i=|g.s.`$, i.e. the ground state of the initial nucleus. The summation in Eq. (44) runs over all states of the final nucleus lying up to $`25`$ MeV. ($`|f|J_\rho ^\pi `$. We consider the final nuclear states as having well-defined spin (J) and parity ($`\pi `$). The index $`\rho `$ counts the multipole states).

For the Fermi and Gamow-Teller contributions, which are expected to be the most important contributions, the square of the matrix element $`|_{if}|^2`$ is written as

$$|_{if}|^2=\frac{1}{2J_i+1}|\left(\frac{f_V}{f_A}\right)^2f||\mathrm{\Omega }_F||i(gs)f||\mathrm{\Omega }_{GT}||i(gs)|^2$$
(46)
where $`f_V`$ and $`f_A`$ are the usual vector and axial vector coupling constants $`(f_A/f_V=1.25)`$. By combining Eqs. (44) and (46) we see that, for the evaluation of the branching ratio $`R_{\mu e^+}`$, we have to calculate the reduced matrix elements $`f\mathrm{\Omega }_Fi`$ and $`f\mathrm{\Omega }_{GT}i`$ for $`|i=|g.s.`$ and $`|f`$ any accessible state of the final nucleus. In the present work these states have been constructed (for 27Al and 27Na systems) in the framework of the shell model as is described in the next section.

## 4 The shell model nuclear wave functions

The evaluation of the reduced matrix elements $`f\mathrm{\Omega }_Fi`$ and $`f\mathrm{\Omega }_{GT}i`$ for the $`\mu ^{}e^+`$ conversion requires reliable nuclear wave functions. These were obtained in the framework of the $`1s0d`$ shell model using the realistic effective interaction of Wildenthal . As we have already mentioned, this interaction has been tested over many years. It is known to accurately reproduce many nuclear observables for s-d shell nuclei. The Wildenthal two-body matrix elements as well as the single particle energies are determined by least square fits to experimental data in the region of the periodic table with $`A=1739`$.

The eigenstates of the daughter nucleus 27Na were evaluated in the isospin representation. For each spin $`J_f`$ with $`T=5/2`$ the first states reaching up to $`E_x=25`$ MeV in excitation energy were calculated. On the other hand, for 27Al we evaluated the ground state $`(5/2)_{gs}^+`$ with $`T=1/2`$. We also constructed all the excited (positive parity) states up to 5 MeV in order to check our predictions against experiment. In Fig. 1 we present the calculated and measured  low-energy spectrum of 27Al up to 5 MeV. As can be seen from this figure, the spectrum of 27Al is well reproduced. In the case of the unstable 27Na isotope the comparison between theory and experiment can not be accomplished due to lack of experimental data.

For the special case of the reaction (13) studied in the present work, since $`M(^{27}Al)M(^{27}Na)=10.6`$ MeV, the momentum transfer at which our matrix elements must be computed, is given by

$$p_e=(p_e)_{max}E_x/c,(p_e)_{max}=94.5(MeV/c)$$
(47)
## 5 Results and discussion

As we have mentioned, the primary purpose of the present work is the calculation of the total $`\mu ^{}e^+`$ reaction rate by summing over partial transition strengths. As a first step we restricted ourselves in the light neutrino mechanism. As a result we only dealt with the Fermi type and the Gamow-Teller type operators discussed in Sect. 2. In other words we evaluated the reduced matrix elements

$$M_F=f\mathrm{\Omega }_Fi(gs)$$
(48)
and

$$M_{GT}=f\mathrm{\Omega }_{GT}i(gs)$$
(49)
for the transitions between the initial $`|i=(5/2)_{gs}^+`$ and all the final $`|f|J_\rho ^\pi `$ states up to 25 MeV. We evaluated the contributions arising from both the real part and imaginary part of the radial function F(r) Eq. (25) which occur in the case of light Majorana neutrinos.

In the calculation of Fermi $`|M_F|^2`$ and Gamow-Teller $`|M_{GT}|^2`$ strengths we found that the contribution from the multipolarities $`L=2`$ and $`L=4`$ is, in general, quite small compared to that of $`L=0`$. This becomes obvious by glancing at Table 1 where the total Fermi and Gamow-Teller strengths with respect to multipolarities $`L`$ are listed. As can be seen from Table 1 the Gamow-Teller strength is almost 9 times greater than the Fermi one which means that it dominates the total strength.

In order to compare the branching ratio originating from the $`g.s.g.s.`$ transition with that associated with the transition to all final $`(5/2)^+`$ states, we define, for convenience, the ratio

$$\lambda \frac{R_{gs}}{R}=\frac{|_{gsgs}|^2}{_fS_{gsf}},$$
(50)
where

$$S_{gsf}=\left(1\frac{E_x}{(p_e)_{max}}\right)^2|_{gsf}|^2,$$
(51)
The quantity $`\lambda `$ gives the portion of the $`g.s.g.s.`$ contribution relative to the total rate here computed by the sum over all partial transitions included in our model space. For the $`g.s.g.s.`$ transition $`p_e=94.5`$ MeV/c according to Eq. (47). Since $`m_ec<<p_e`$ we can consider the approximation $`p_eE_e/c`$, which is equivalent to neglecting the electron mass ($`m_e`$) in the kinematics of the reaction (13).

For the explicit contribution of the excited states to the branching ratio we present in Table 2 the sum $`_fS_{gsf}`$ for each set of excited states of given $`|J_f`$. The second column of Table 2 corresponds to the contribution from the real part of Eq. (25) while the third column represents the total contribution coming from the real and imaginary parts of Eq. (25). As can be seen from Table 2, the ground state transition exhausts a large portion (41%) of the total $`\mu ^{}e^+`$ reaction rate. As can also be seen from Table 2, the main contribution to the total rate comes from the $`(5/2)^+`$ multipole states which contribute about 98% of the total rate. The rest, $`2\%`$, originates mainly from the $`(3/2)^+`$, $`(7/2)^+`$ and $`(9/2)^+`$ states. Another interesting conclusion stemming out of the results of Table 2 is the fact that, the contribution coming from the imaginary component of the radial function F(r) Eq. (25) dominates the total branching ratio. In fact the contribution of the real part is about 20 times smaller than the corresponding imaginary one.

In order to have an insight about the origin of this difference we studied the behavior as a function of the relative coordinate r of the following three quantities:

(i) $`I_{nn^{}}(r)=R_{n0}R_{n^{}0}r^2bIm\{F(r)\}`$,

(ii) $`S_{nn^{}}(r)=R_{n0}R_{n^{}0}r^2bRe\{F(r)\}`$,

where $`n,n^{^{}}`$ indicate the the nodes in the relative coordinate of the two nucleon wave function. $`Im\{F(r)\}`$ and $`Re\{F(r)\}`$ are the imaginary and real parts, respectively, of $`F(r)`$ \[see Eq. (28)\]. Since there is no interference between the real and imaginary parts we find it convenient to define $`I_{nn^{}}`$ with the (-) sign. For calculations within the $`0d1s`$ shell we make the plausible assumption that the main contribution comes from zero angular momentum states in the relative motion.

(iii) $`D_{nn^{}}(r)=R_{n0}R_{n^{}0}r^2b`$, which expresses the nuclear densities and corresponds to putting $`F(r)1`$ (this is the case of $`0\nu \beta \beta `$-decay).

$`R_{nl}`$ denotes the radial part of the harmonic oscillator wave functions entering the matrix elements of Eqs. (48) and (49). Thus, $`D_{nn^{}}`$, $`I_{nn^{}}`$ and $`S_{nn^{}}`$ are dimensionless quantities of which the squares determine the magnitudes of the two-body matrix elements involved in our evaluations of the partial transition rates. The results obtained for $`D_{nn^{}}(r)`$, $`I_{nn^{}}(r)`$ and $`S_{nn^{}}(r)`$ with $`n,n^{}=0,1`$, are plotted in Fig. 2 from where we conclude the following:

The variation of $`I_{nn^{}}(r)`$ which contains the imaginary part of $`F(r)`$ is in all cases in phase with $`D_{nn^{}}(r)`$. The peaks of these quantities appear at about the same value of r for each case. On the contrary, the variation of $`S_{nn^{}}(r)`$ which contain the real part of $`F(r)`$, is not always in phase with $`D_{nn^{}}(r)`$ and $`I_{nn^{}}(r)`$. Also the peaks of $`S_{nn^{}}(r)`$ do not appear at the same positions of $`r`$ as those of $`D_{nn^{}}(r)`$ and $`I_{nn^{}}(r)`$. This picture appears in all other cases resulting by using the various $`R_{nl}`$ wave functions entering our calculations, even though the absolute strengths of $`D_{nn^{}}(r)`$, $`I_{nn^{}}(r)`$ and $`S_{nn^{}}(r)`$ are much smaller compared to those of Fig. 2.

Obviously, the area bounded by the r-axis and $`I_{nn^{}}`$ or $`S_{nn^{}}`$ gives a measure of the contribution of the imaginary or real part of the $`\mu ^{}e^{}`$ operators, respectively, A rough estimation of this area gives for the ratio of imaginary to real a value of about 4-5 which justifies the factor of about 15-25 between imaginary and real part contributions to the partial sum evaluation of the total $`\mu ^{}e^{}`$ rate. Of course, this difference of the imaginary part need not apply when the model is expanded outside our space, i.e. to include $`1\mathrm{}\omega `$ and $`2\mathrm{}\omega `$ excitations.

At this point we found interesting to study the spreading of the contributions through the excitation spectrum of the final nucleus. To this end, in Fig. 3 we plot the distribution of $`S_{gsf}`$ for all $`(5/2)^+`$ states. Similar pictures are obtained for the other multipole states $`(1/2)^+`$, $`(3/2)^+`$, $`(7/2)^+`$, $`(9/2)^+`$, $`(11/2)^+`$, $`(13/2)^+`$, even though their contribution is quite smaller than that of $`(5/2)^+`$. The common feature of these distributions is the fact that for each multipolarity the main contribution comes from low-lying states and that excitations lying above $``$10 MeV contribute negligible quantities. In some cases, as e.g. $`(7/2)^+`$ the contribution comes from a narrow window of the excitation spectrum of the final nucleus.

At this stage we should mention that the total rates can also be obtained in the context of the closure approximation. We remind that according to the closure approximation the contribution of each individual state is effectively taken into account by assuming a mean excitation energy $`\overline{E}=E_fE_{gs}`$, and using the completeness relation $`_f|ff|=1`$. Therefore

$$\underset{f}{}|f|\mathrm{\Omega }|i|^2=i|\mathrm{\Omega }^+\mathrm{\Omega }|i$$

The matrix element $`i|\mathrm{\Omega }^+\mathrm{\Omega }|i`$ can be written as a sum of two pieces: a two-body term and a four-body one, that is

$$i|\mathrm{\Omega }^+\mathrm{\Omega }|i=i|(\mathrm{\Omega }^+\mathrm{\Omega })_{2b}|i+i|(\mathrm{\Omega }^+\mathrm{\Omega })_{4b}|i$$
(52)
The total rate evaluation thus involves only the $`g.s.`$ of the initial nucleus. The 4-body piece was neglected in the earlier calculations. Thus, following the work of Ref. , the relevant two-body operator for the process studied in the present work takes the form:

$$\mathrm{\Omega }^+\mathrm{\Omega }=\underset{i\mathrm{\#}j}{}\tau _{}(i)\tau _{}(j)\tau _+(i)\tau _+(j)[\frac{f_V^2}{f_A^2}\sigma (i).\sigma (j)][\frac{f_V^2}{f_A^2}\sigma (i).\sigma (j)](\frac{R_0}{r_{ij}})^2$$
(53)
The above equation can be written as follows:

$$\mathrm{\Omega }^+\mathrm{\Omega }=\underset{i\mathrm{\#}j}{}\tau _{}(i)\tau _+(i)\tau _{}(j)\tau _+(j)[9+(\frac{f_V^2}{f_A^2})^22(\frac{f_V^2}{f_A^2}+1)\sigma (i).\sigma (j)](\frac{R_0}{r_{ij}})^2$$
(54)
To proceed further one makes the simplifying assumption that the matrix element is dominated by the spin singlet states. Then the isospin operator counts the number of proton pairs in the initial nucleus. The radial part is estimated by assuming a uniform nuclear two body density, i.e.

$$(\frac{R_0}{r_{ij}})^2=\frac{1}{(2R_0)^3}_{r_c}^{2R_0}(\frac{R_0}{r})^2r^2𝑑r$$
(55)
where $`r_c`$ is the hard core radius, assuming for simplicity that the short range correlations are described by a simple step function. In any case since $`r_cR_0`$, the short range correlations can be neglected and to a good approximation the value of the above integral is $`1/4`$. We thus find that

$$i|\mathrm{\Omega }^+\mathrm{\Omega }|i=Z(Z1)(\frac{f_V^2}{f_A^2}+3)^2\frac{1}{4}$$
(56)
The above matrix element must be multiplied by the kinematical factor $`[p_e/(p_e)_{max}]^20.8`$, taking an average excitation energy of about $`20MeV`$. We thus find:

$$i|\mathrm{\Omega }^+\mathrm{\Omega }|i0.8\times 13\times 12\times (3.7)^2\times 0.25=430$$
(57)
It is clear, therfore, that the contribution of the excited states found by our state-by-state calculation in the present work is much smaller than the predictions found previously by employing closure approximation as outlined above.

The disagreement appeared between closure approximation and the present state-by-state calculation can be attributed to the following reasons:

i) The closure approximation takes into account not only the contribution of $`0\mathrm{}\omega `$ space but also includes excitations $`En\mathrm{}\omega ,n1`$, as well as the continuum spectrum. A possible extension of the $`sd`$ model space to include two-particle two-hole states with the above harmonic oscillator excitations is practically impossible.

ii) In closure approximation the second term in Eq. (52) which includes the four-body forces and which is very complicated, was not taking into account in the previous calculations. Of course, the obvious question raised is: how important the contribution of four-body terms is ?

iii) From the ordinary $`\mu `$-capture it is known that the results of the simple closure approximation are quite sensitive to the assumed mean excitation energy $`\overline{E}`$. Since the spectrum of the final nucleus reached by the operators of $`\mu ^{}e^+`$ conversion has not yet been observed, $`\overline{E}`$ may be naively estimated from the spectrum of the final nucleus.

Anyway putting our nuclear matrix elements and all the other nuclear physics input into Eq. (44) we get

$$R_{\mu e^+}=2.4\rho |\eta _{eff}|^2.$$
(58)
Combining this with the present experimental limit

$$R_{\mu e^+}4.4\times 10^{12}\text{[13]},$$

we get

$$\rho |\eta _{eff}|^21.8\times 10^{12}$$
(59)
The quantity $`\rho |\eta _{eff}|^2`$ depends on the particle model adopted as well as the prevailing mechanism for this process and it will not be further discussed in this work.

## 6 Summary and Conclusions

In the present work we have investigated the exotic double-charge exchange neutrinoless muon-to-positron conversion in the presence of nuclei, $`(A,Z)\mu ^{}e^+(A,Z2)`$. The appropriate operators have been constructed considering a variety of gauge model predictions, not only light intermediate Majorana neutrinos. We have chosen to study the experimentally interesting nucleus 27Al, since it is going to be used as a stopping target in the Brookhaven experiment. This nucleus is expected to be described rather well within the $`1s0d`$ shell model, since one can use the well-tested Wildenthal $`1s0d`$ interaction.

As a first step we restricted to the calculation of the rates in the case of light intermediate neutrinos.

¿From our results on the reaction $`{}_{}{}^{27}Al`$$`(\mu ^{},e^+)`$$`{}_{}{}^{27}Na`$, we can conclude the following:

(i) The contribution coming from the Gamow–Teller component of the $`\mu ^{}e^+`$ operator dominates the total rate matrix elements.

(ii) The contribution of $`L=0`$ multipolarity dominates the total rate.

(iii) The total strength, resulting by summing over partial transition matrix elements included in our model space, is much smaller than that found previously by using closure approximation.

(iv) The contribution arising from the imaginary component of the radial part in the $`\mu ^{}e^+`$ conversion operator is much larger than the one of the real component.

(v) The portion of the $`g.s.g.s.`$ contribution (which is proportional to the matrix element $`|_{gsgs}|^2`$) into the total rate is 41%. This is good news, since, eventually, the experiments will focus on the ground state.

(vi) By putting our nuclear physics input into Eq. (44) we obtain the limit $`\rho |\eta _{eff}|^21.8\times 10^{12}`$. The specific limits on $`\rho `$ and $`\eta _{eff}`$ depend on the particle model assumed and the prevailing mechanism for the $`\mu ^{}e^{}`$ conversion process.

Appendix

According to the gauge models mentioned in Sect. 2 the transition operator $`\mathrm{\Omega }`$ can be given in terms of the following components :

$$\mathrm{\Omega }_{Sa}=\underset{ij}{}\tau _{}(i)\tau _{}(j)e^{ı𝐩_e𝐫_i}f(r_{ij})W_{Sa}(ij),a=1,2$$
(60)
$$\mathrm{\Omega }_{Aa}=\underset{ij}{}\tau _{}(i)\tau _{}(j)e^{ı𝐩_e𝐫_i}f(r_{ij})W_{Aa}(ij),a=1,2$$
(61)
$$\mathrm{\Omega }_T^{^{}}=\underset{ij}{}\tau _{}(i)\tau _{}(j)e^{ı𝐩_e𝐫_i}f(r_{ij})W_{S3}(ij)$$
(62)
$$\mathrm{\Omega }_{A_2}^{^{}}=\underset{ij}{}\tau _{}(i)\tau _{}(j)e^{ı𝐩_e𝐫_i}f(r_{ij})W_{A3}(ij)$$
(63)
where

$$W_{S1}(ij)=1$$
(64)
$$W_{S2}(ij)=\sigma _i\sigma _j=\sqrt{3}\left[\sigma _i\sigma _j\right]_0^0$$
(65)
$$W_{S3}(ij)=3(\sigma _i\widehat{r}_{ij})(\sigma _j\widehat{r}_{ij})\sigma _i\sigma _j=\sqrt{6}\left[\sqrt{4\pi }Y^2(\widehat{r}_{ij})\left[\sigma _i\sigma _j\right]^2\right]_0^0$$
(66)
$$W_{A1}(ij)=ı\sigma _i\times \sigma _j=(\sqrt{2})\left[\sigma _i\sigma _j\right]^1$$
(67)
$$W_{A2}(ij)=\sigma _i\sigma _j$$
(68)
$$W_{A3}(ij)=(\sigma _i\sigma _j)(ı\widehat{r}_{ij}\times \widehat{R}_{ij})=\sqrt{\frac{2}{3}}\left[\left[\sqrt{4\pi }Y^1(\widehat{r}_{ij})\sqrt{4\pi }Y^1(\widehat{R}_{ij})\right]^1(\sigma _i\sigma _j)\right]_0^0$$
(69)
By Applying the usual multipole decomposition procedure the operators $`\mathrm{\Omega }_{Sa}`$, $`\mathrm{\Omega }_{Aa}`$, $`\mathrm{\Omega }_T^{^{}}`$ and $`\mathrm{\Omega }_{A_2}^{^{}}`$ read

$$\mathrm{\Omega }_{Sa}=\underset{L\mathrm{\Lambda }J}{}\sqrt{4\pi }Y^\mathrm{\Lambda }(\widehat{p}_e)O_{Sa}^{(L,S)J}\delta _{\mathrm{\Lambda }L}\delta _{LJ}$$
(70)
$`\mathrm{\Omega }_{Aa}={\displaystyle \underset{LJ\mathrm{\Lambda }}{}}\left[\sqrt{4\pi }Y^\mathrm{\Lambda }(\widehat{p}_e)O_{Aa}^{(L,S)J}\right]^1\delta _{\mathrm{\Lambda }L}`$ (71)
$$\mathrm{\Omega }_T^{^{}}=\underset{L\mathrm{\Lambda }J}{}\sqrt{4\pi }Y^\mathrm{\Lambda }(\widehat{p}_e)O_T^{(L,S)J}\delta _{\mathrm{\Lambda }J}$$
(72)
$$\mathrm{\Omega }_{A2}^{^{}}=\underset{L\mathrm{\Lambda }J}{}\sqrt{4\pi }Y^\mathrm{\Lambda }(\widehat{p}_e)O_{A2}^{{}_{}{}^{}(L,S)J}\delta _{\mathrm{\Lambda }J}$$
(73)
The operators $`O_{Sa}^{(L,S)J}`$, $`O_{Aa}^{(L,S)J}`$, $`O_T^{(L,S)J}`$ and $`O_{A2}^{{}_{}{}^{}(L,S)J}`$ are given by the following equations

$`O_{S1}^{(L,S)J}\mathrm{\Omega }_F={\displaystyle \underset{l}{}}A_{lL}{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$
$`\left[\sqrt{4\pi }Y^l(\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})\right]^J,S=0,J=L`$ (74)
$`O_{S2}^{(L,S)J}\mathrm{\Omega }_{GT}={\displaystyle \underset{l}{}}A_{lL}{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$
$`\left[\left[\sqrt{4\pi }Y^l(\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})\right]^L(\sqrt{3})\left[\sigma _i\sigma _j\right]^0\right]^J,S=0,J=L`$ (75)
$`O_{A1}^{(L,S)J}={\displaystyle \underset{l}{}}B_{lL}{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$
$`\left[\left[\sqrt{4\pi }Y^l(\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})\right]^L(\sqrt{2})\left[\sigma _i\sigma _j\right]^1\right]^J,S=1,J=L,|L\pm 1|`$ (76)
$`O_{A2}^{(L,S)J}={\displaystyle \underset{l}{}}B_{lL}{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$
$`\left[\left[\sqrt{4\pi }Y^l(\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})\right]^L\left(\sigma _i\sigma _j\right)\right]^J,S=1,J=L,|L\pm 1|`$ (77)
$`O_T^{(L,S)J}={\displaystyle \underset{l\lambda }{}}E_{l\lambda L}^\mathrm{\Lambda }{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f^{^{}}(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$
$`\left[\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})\right]^L\left[\sigma _i\sigma _j\right]^2\right]^J,\lambda =l,|l\pm 2|,S=2,\mathrm{\Lambda }=J`$ (78)
$`O_{A_2}^{{}_{}{}^{}(L,S)J}={\displaystyle \underset{l}{}}{\displaystyle \underset{k\lambda }{}}{\displaystyle \underset{\mu }{}}\mathrm{\Theta }_{k\lambda \mu L}^{\mathrm{\Lambda }l}{\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij})j_l({\displaystyle \frac{p_er_{ij}}{2}})j_{}(p_eR_{ij})`$
$`\left[\left[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^\mu (\widehat{R}_{ij})\right]^L\left(\sigma _i\sigma _j\right)\right]^J,S=1,J=\mathrm{\Lambda }`$ (79)
where

$`E_{l\lambda L}^\mathrm{\Lambda }=D_{l\lambda L}^\mathrm{\Lambda }A_{l\mathrm{\Lambda }}`$ (80)
$`D_{l\lambda L}^\mathrm{\Lambda }=(1)^{l++\mathrm{\Lambda }}\sqrt{30\widehat{L}\widehat{\mathrm{\Lambda }}\widehat{l}}\left\{\begin{array}{ccc}2& l& \lambda \\ & L& \mathrm{\Lambda }\end{array}\right\}\text{ }\left\{\begin{array}{ccc}\mathrm{\Lambda }& 2& L\\ 2& \mathrm{\Lambda }& 0\end{array}\right\}\text{ }20l0|\lambda 0`$ (85)
$`\mathrm{\Theta }_{k\lambda \mu L}^{\mathrm{\Lambda }l}=A_{l\mathrm{\Lambda }}I_{k\lambda \mu L}^{\mathrm{\Lambda }l}`$ (86)
$`I_{k\lambda \mu L}^{\mathrm{\Lambda }l}=3\sqrt{2}(1)^{l+\lambda +1}\widehat{k}\sqrt{\widehat{l}\widehat{}\widehat{L}\widehat{\mathrm{\Lambda }}}10l0|\lambda 0010|\mu 0\left\{\begin{array}{ccc}\mathrm{\Lambda }& 1& k\\ 1& L& 1\end{array}\right\}`$ (89)
$`\left\{\begin{array}{ccc}1& l& \lambda \\ & k& \mathrm{\Lambda }\end{array}\right\}`$ $`\left\{\begin{array}{ccc}\lambda & & k\\ 1& L& \mu \end{array}\right\}`$ $`\left\{\begin{array}{ccc}\mathrm{\Lambda }& 1& L\\ 1& \mathrm{\Lambda }& 0\end{array}\right\}`$ (96)
with

$$A_{l\mathrm{\Lambda }}=\sqrt{\frac{\widehat{l}\widehat{}}{\widehat{\mathrm{\Lambda }}}}l00|\mathrm{\Lambda }0(1+(1)^l)ı^{l+}$$
(97)
and

$$B_{l\mathrm{\Lambda }}=\sqrt{\frac{\widehat{l}\widehat{}\widehat{J}}{3\widehat{\mathrm{\Lambda }}}}(1)^{J+1}l00|\mathrm{\Lambda }0(1(1)^l)ı^{l+}$$
(98)
($`\widehat{a}2a+1`$).

We are now going to discuss the operators entering the leptonic R-L interference and in some SUSY mechanisms. Now the relevant operators may have a time component, which is small and, in any case, except for their radial part is the same with the $`F`$ and $`GT`$ discussed above. They also have a space component, which is proportional to $`\widehat{r}_{ij}`$. They are vectors, which yield a scalar, when combined with the leptonic current. They are of the form

$$\sigma _i(\sigma _j\widehat{r}_{ij})+(\sigma _i\widehat{r}_{ij})\sigma _j(\sigma _i\sigma _j)\widehat{r}_{ij}\text{ and }i(\sigma _i\sigma _j)\times \widehat{r}_{ij}$$

The above operators can also be written as

$$\sigma _i(\sigma _j\widehat{r}_{ij})+(\sigma _i\widehat{r}_{ij})\sigma _j(\sigma _i\sigma _j)\widehat{r}_{ij}=\omega ^{^{}}(0)+\omega ^{^{}}(2)$$
(99)
and

$$i(\sigma _i\sigma _j)\times \widehat{r}_{ij}=\omega ^{^{}}(1)$$
(100)
where

$$\omega ^{^{}}(k)=\alpha (k)[\sqrt{4\pi }Y^1(\widehat{r}_{ij})T^k(spin)]^1,k=0,1,2$$
(101)
$$\alpha (0)=\frac{1}{3\sqrt{3}},T^0(spin)=\sigma _i\sigma _j$$
(102)
$$\alpha (2)=\frac{2\sqrt{5}}{3},T^2(spin)=[\sigma _i\sigma _j]^2$$
(103)
$$\alpha (1)=\sqrt{\frac{2}{3}},T^1(spin)=\sigma _i\sigma _j$$
(104)
The above operators are accompanied by the lepton outgoing waves

$$O_{ij}=exp(i𝐩_e𝐑_{ij})[exp(i𝐩_e𝐫_{ij})/2+(1)^{k+1}exp(i𝐩_e𝐫_{ij})/2]$$
(105)
The phase of the second term guarantees that the combined operator is overall symmetric under the exchange of the particles $`i`$ and $`j`$. The last operator can be brought into the form

$$O_{ij}=\underset{l\lambda \mathrm{\Lambda }}{}[\beta (lk\mathrm{\Lambda })j_l(p_er_{ij}/2)j_{}(p_eR_{ij})[[\sqrt{4\pi }Y^l(\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})]^\mathrm{\Lambda }\sqrt{4\pi }Y^\mathrm{\Lambda }(\widehat{p}_e)]^0$$
(106)
where

$$\beta (lk\mathrm{\Lambda })=[1+(1)^{l+k+1}]i^{l+}(1)^{l+}[(2l+1)(2+1]^{1/2}l00|\mathrm{\Lambda }0$$
(107)
Combining the above factors we obtain

$`\mathrm{\Omega }^{^{}}(k)={\displaystyle \underset{LJ\mathrm{\Lambda }}{}}\left[\sqrt{4\pi }Y^\mathrm{\Lambda }(\widehat{p}_e)O^{(L,S)J}(k)\right]^1`$ (108)
where

$`O^{(L,S)J}(k)`$ $`=`$ $`\alpha (k){\displaystyle \underset{i<j}{}}\tau _{}(i)\tau _{}(j)f(r_{ij}){\displaystyle \underset{l\lambda }{}}\beta (lk\mathrm{\Lambda })\gamma (l,\lambda ,,L,k,J,\mathrm{\Lambda })j_l(p_er_{ij}/2)`$ (109)
$`j_{}(p_eR_{ij})[[\sqrt{4\pi }Y^\lambda (\widehat{r}_{ij})\sqrt{4\pi }Y^{}(\widehat{R}_{ij})]^LT^k(spin)]^J`$
with

$`\gamma (l,\lambda ,,L,k,J,\mathrm{\Lambda })`$ $`=`$ $`[(3(2L+1)(2l+1)(2J+1)]^{1/2}`$ (115)
$`10l0|\lambda 0\left\{\begin{array}{ccc}1& k& 1\\ J& \mathrm{\Lambda }& L\end{array}\right\}\text{ }\left\{\begin{array}{ccc}1& l& \lambda \\ & L& \mathrm{\Lambda }\end{array}\right\}`$
FIGURE CAPTIONS

FIGURE 1 The calculated (right) and measured (left)  energy spectrum for the lowest positive parity states of 27Al.

FIGURE 2. Variation of the quantities $`D_{nn^{}}`$, $`I_{nn^{}}`$ and $`S_{nn^{}}`$ as functions of the relative coordinate r (see Sect. 5) assuming that the relative zero angular momentum states dominate. Since the relative motion is indicated in the plots one sees some contribution beyond the nuclear radius. In this figure only the most prominent cases for $`n,n^{}=0,1`$ are shown. The behavior of $`D_{nn^{}}(r)`$, $`I_{nn^{}}(r)`$ and $`S_{nn^{}}(r)`$ in the other cases needed for our calculations is similar.

FIGURE 3. Distribution of the transition strength $`S_{gsf}`$ for all $`J^\pi =(5/2)^+`$ states. Both the real and imaginary parts of Eq. (25) are included.


text:
# Probability density function of turbulent velocity fluctuation

## I INTRODUCTION

While velocity differences in turbulence $`u(x+\delta x)u(x)`$ have attracted much interest, velocity fluctuations $`u(x)`$ themselves are also fundamental to describing the turbulence. Usually it is observed that the probability density function (PDF) of the velocity fluctuations is close to Gaussian B53 ; VC69 . However, some experiments yield a sub-Gaussian PDF NWLMF97 ; SD98 ; J98 , which has a less pronounced tail than a Gaussian PDF. The reason remains controversial because there have been no systematical studies.

The observed Gaussianity had been explained by applying the cental limit theorem to the Fourier transformation of the velocity fluctuations B53 . This theorem ensures that a sum of many independent random variables has a Gaussian PDF KS77 . However, owing to the steep energy spectrum of turbulence, the Fourier transforms have considerably different magnitudes. In this case, the central limit theorem is not applicable J98 . We first reconsider the condition for the velocity fluctuations to have a Gaussian PDF (Sec. II).

For studies of velocity fluctuations, an experimental approach is preferable to the popular direct numerical simulation. Although the direct numerical simulation is useful in studying small-scale motions of turbulence, they are not of our interest because the velocity fluctuations are dominated by energy-containing large-scale motions. Those in the direct numerical simulation suffer from its artificial boundary condition, initial condition, and forcing.

We obtain experimental data of velocity fluctuations in a wind tunnel downstream of a turbulence-generating grid (Sec. III). The grid turbulence is not homogeneous in the mean-wind direction and thus allows us to study developing, fully developed, and decayed states of turbulence by increasing the distance between the anemometer and the grid. We find that the PDF of velocity fluctuations changes accordingly from a sub-Gaussian to Gaussian, and to hyper-Gaussian (Sec. IV). Here a hyper-Gaussian PDF has a more pronounced tail than a Gaussian PDF. We discuss the observed behavior by using the velocity fluctuations themselves, the velocity differences, and the Fourier transforms.

## II THEORY FOR GAUSSIANITY

Suppose that velocity fluctuations $`u(x)`$ are measured in a turbulent flow repeatedly over the range $`0x<L`$. The length $`L`$ is much greater than the correlation length $`l_u=u(x+\delta x)u(x)d\delta x/u^2`$. Here $``$ denotes an ensemble average. The measurements serve as realizations of the turbulence. Each of them is expanded into a Fourier series as

$$u(x)=\underset{n=1}{\overset{\mathrm{}}{}}a_n(u)\mathrm{cos}\left(\frac{2\pi nx}{L}\right)+b_n(u)\mathrm{sin}\left(\frac{2\pi nx}{L}\right).$$
(1)
The Fourier transforms $`a_n(u)`$ and $`b_n(u)`$ for $`n1`$ have Gaussian PDFs over the realizations BP01  (see also Ref. FB95  for a more mathematical explanation). This is because, for example, the Fourier transform $`a_n(u)`$ is obtained as

$`a_n(u)`$ $`=`$ $`{\displaystyle \frac{2}{L}}{\displaystyle _0^L}u(x)\mathrm{cos}\left({\displaystyle \frac{2\pi nx}{L}}\right)dx`$ (2)
$`=`$ $`{\displaystyle \frac{2}{L}}\left({\displaystyle _0^{L/m}}\mathrm{}dx+{\displaystyle _{L/m}^{2L/m}}\mathrm{}dx+\mathrm{}+{\displaystyle _{(m1)L/m}^L}\mathrm{}dx\right),`$
with $`1mn`$. The segment size $`L/m`$ is set to be large enough so that the correlation $`u(x+\delta x)u(x)`$ has converged to zero at $`\delta x=L/m`$. The integrations $`_0^{L/m}\mathrm{}dx`$, $`_{L/m}^{2L/m}\mathrm{}dx`$, …, and $`_{(m1)L/m}^L\mathrm{}dx`$ are regarded as independent random variables of the same magnitude. Then $`a_n(u)`$ has a Gaussian PDF as a consequence of the central limit theorem. The variance $`a_n^2(u)`$ is equal to the energy spectrum $`E_n`$ 111 Note that the norm of the basis functions is not unity. If the norm were unity, the variance of the transform would be equal to $`E_n/2`$..

The above discussion is not applicable to the Fourier transforms $`a_n(u)`$ and $`b_n(u)`$ for $`n1`$, i.e., those for large wavelengths. Nevertheless, as far as the wavelength is finite, it is possible to show Gaussianity of the corresponding Fourier transform by increasing the data length $`L`$ and hence the $`n`$ value. In the limit $`L\mathrm{}`$, the transforms $`a_n(u)`$ and $`b_n(u)`$ for $`n1`$ become zero and do not contribute to the velocity fluctuations $`u(x)`$. We are able to assume safely that all the Fourier transforms have Gaussian PDFs. Since this fact is independent of detailed dynamics, it is universal.

If and only if all the Fourier transforms $`a_n(u)`$ and $`b_n(u)`$ are statistically independent of each other, the Gaussianity of them leads to Gaussianity of the velocity fluctuations $`u(x)`$. To demonstrate this, we use the characteristic functions $`\varphi _n^{(a)}(\tau )`$ and $`\varphi _n^{(b)}(\tau )`$ for the PDFs of $`a_n(u)\mathrm{cos}(2\pi nx/L)`$ and $`b_n(u)\mathrm{sin}(2\pi nx/L)`$ at any fixed spatial position $`x`$ KS77 :

$`\varphi _n^{(a)}(\tau )`$ $`=`$ $`\mathrm{exp}\left[{\displaystyle \frac{\tau ^2}{2!}}E_n\mathrm{cos}^2\left({\displaystyle \frac{2\pi nx}{L}}\right)\right]`$ (3a)
$`\varphi _n^{(b)}(\tau )`$ $`=`$ $`\mathrm{exp}\left[{\displaystyle \frac{\tau ^2}{2!}}E_n\mathrm{sin}^2\left({\displaystyle \frac{2\pi nx}{L}}\right)\right].`$ (3b)
Here $`E_n\mathrm{cos}^2(2\pi nx/L)`$ and $`E_n\mathrm{sin}^2(2\pi nx/L)`$ are the variances of $`a_n(u)\mathrm{cos}(2\pi nx/L)`$ and $`b_n(u)\mathrm{sin}(2\pi nx/L)`$. From the expansion formula (1) and the independence of the Fourier transforms $`a_n(u)`$ and $`b_n(u)`$, it follows that the sum of the logarithms of the characteristic functions $`\varphi _n^{(a)}(\tau )`$ and $`\varphi _n^{(b)}(\tau )`$ is equal to the logarithm of the characteristic function $`\varphi (\tau )`$ for the PDF of the velocity fluctuations $`u(x)`$:

$$\mathrm{ln}\varphi (\tau )=\underset{n=1}{\overset{\mathrm{}}{}}\mathrm{ln}\varphi _n^{(a)}(\tau )+\mathrm{ln}\varphi _n^{(b)}(\tau )=\frac{\tau ^2}{2!}\underset{n=1}{\overset{\mathrm{}}{}}E_n.$$
(4)
Thus the velocity fluctuations $`u(x)`$ have a Gaussian PDF with the variance $`u(x)^2=_{n=1}^{\mathrm{}}E_n`$. The independence of the Fourier transforms also leads to the statistical independence of the velocity fluctuations $`u(x)`$ at different spatial positions.

Therefore, the necessary and sufficient condition for the velocity fluctuations to have a Gaussian PDF is the independence of the Fourier transforms. This is a good approximation for fully developed turbulence, where large-scale motions of energy-containing eddies are random and independent. Although turbulence contains small-scale coherent structures such as vortex tubes SA97 , their importance to the velocity fluctuations is negligible, i.e., being as small as the energy ratio of the dissipation range to the energy-containing range.

## III EXPERIMENTS

The experiments were done in two wind tunnels of Meteorological Research Institute. Their test sections were of $`0.8\times 0.8\times 3`$ and $`3\times 2\times 18\mathrm{m}`$ in size (hereafter, respectively, the small and large tunnels). The small tunnel was used to study developing and fully developed states of grid turbulence, while the large tunnel was used to study fully developed and decayed states.

Turbulence was produced by placing a grid across the entrance to the test section. The grid consisted of two layers of uniformly spaced rods, the axes of which were perpendicular to each other. We used different grids in the small and large tunnels. The separation of the axes of their adjacent rods were 0.10 and 0.40 m, respectively. The cross sections of the rods were $`0.02\times 0.02`$ and $`0.06\times 0.06\mathrm{m}`$, respectively. The mean wind was set to be $`U10\mathrm{m}\mathrm{s}^1`$ in the small tunnel and $`U20\mathrm{m}\mathrm{s}^1`$ in the large tunnel.

We simultaneously measured the streamwise ($`U+u`$) and transverse ($`v`$) velocities. They are velocity components that are parallel and perpendicular to the mean-wind direction, respectively. The measurements in the small tunnel were done on the tunnel axis from $`d=0.25`$ to 2.00 m downstream of the grid with an interval of 0.25 m. Those in the large tunnel were done from $`d=3.00`$ to 17.00 m with an interval of 1.00 m. The ranges of the measurement positions were to the limit of mechanical constraints of the wind tunnels. Since there was no overlap in the distance $`d`$ between the small- and large-tunnel measurements, the individual data are identified by their $`d`$ values.

We used a hot-wire anemometer, which was composed of a crossed-wire probe and a constant temperature system. The wires were 5 $`\mu `$m in diameter, 1.25 mm in effective length, 1.25 mm in separation, and oriented at $`\pm 45^{}`$ to the mean-wind direction. The wire temperature was 280C, while the air temperature was 29–30C in the small tunnel and 14–19C in the large tunnel. We calibrated the anemometer before and after the measurements.

The signal was low-pass filtered with 24 dB/octave and sampled digitally with 16-bit resolution. In the small-wind measurements, the filtering was at 8 kHz and the sampling was at 16 kHz. In the large-wind measurements, the filtering was at 20 kHz and the sampling was at 40 kHz. The entire length of the signal was as long as $`5\times 10^6`$ points. We obtained longer data of $`2\times 10^7`$ points at the positions $`d=0.25`$, 2.00, 8.00, and 12.00 m.

The turbulence levels, i.e., the ratios of the root-mean-square values of the velocity fluctuations $`u^2^{1/2}`$ and $`v^2^{1/2}`$ to the mean streamwise velocity $`U`$, were always low ($`0.2`$; see below). This is a good characteristic of grid turbulence and allows us to rely on the frozen-eddy hypothesis of Taylor, $`/t=U/x`$, which converts temporal variations into spatial variations in the mean-wind direction.

The resultant spatial resolution is comparable to the probe size, $`1\mathrm{mm}`$. Since the probe is larger than the Kolmogorov length, $`0.1`$$`0.2\mathrm{mm}`$, the smallest-scale motions of the flow were filtered out. The present resolution is nevertheless typical of hot-wire anemometry VC69 ; SD98 ; FK67 ; CG99 ; AGHA84 .

Fig. 1 shows the mean streamwise velocity $`U`$, the root-mean-square fluctuations $`u^2^{1/2}`$ and $`v^2^{1/2}`$, the correlation lengths $`l_u`$ and $`l_v`$, the Taylor microscale $`\lambda =(u^2/(u/x)^2)^{1/2}`$, and the turbulence levels $`u^2^{1/2}/U`$ and $`v^2^{1/2}/U`$. For these and the other similar diagrams, open symbols denote the streamwise velocity while filled symbols denote the transverse velocity. The flow parameters change systematically with the distance from the grid, indicating a systematical change of the turbulence. This is especially the case in the turbulence levels. Our following studies suggest $`u^2^{1/2}/U0.1`$ for the developing state, $`u^2^{1/2}/U0.04`$–0.1 for the fully developed state, and $`u^2^{1/2}/U0.04`$ for the decayed state.

## IV RESULTS AND DISCUSSION

### IV.1 Overview

Fig. 2 demonstrates that Fourier transforms of the velocity fluctuations have Gaussian PDFs at the positions $`d=0.25`$, 2.00, 8.00, and 12.00 m. The individual data were divided into 4864 segments of 212 points, which were regarded as independent realizations of the turbulence. They were windowed by a flat-topped function, which rises from zero to unity in the first small fraction of the data and falls back to zero in the last small fraction. The PDFs shown in Fig. 2 are those for the wave number $`k=n/L=7m^1`$, where the energy spectra at $`d=0.25`$ m have a peak (see below). Since the PDFs of $`a_n`$ and $`b_n`$ should be the same, they were put together in order to minimize statistical uncertainties. We also obtained Gaussian PDFs of the Fourier transforms at the other wave numbers and at the other positions in the wind tunnels.

However, velocity fluctuations do not necessarily have Gaussian PDFs. Fig. 3 shows the PDFs at the positions $`d=0.25`$, 2.00, 8.00, and 12.00 m. The transverse-velocity PDF is sub-Gaussian at $`d=0.25`$ m, Gaussian at $`d=2.00`$ and 8.00 m, and hyper-Gaussian at $`d=12.00`$ m. The streamwise-velocity PDF tends to be skewed owing to a shear flow.

Fig. 4 shows the flatness factors $`F_u=u^4/u^2^2`$ and $`F_v`$, the skewness factor $`S_u=u^3/u^2^{3/2}`$, and the streamwise-transverse correlation $`C_{uv}=(u^2v^2u^2v^2)/[(u^4u^2^2)(v^4v^2^2)]^{1/2}`$. With an increase of the distance from the grid, the transverse-velocity PDF changes from sub-Gaussian ($`F<3`$) to Gaussian ($`F=3`$), and to hyper-Gaussian ($`F>3`$). The streamwise-velocity PDF tends to be skewed ($`S0`$) and changes from hyper-Gaussian to Gaussian, and to hyper-Gaussian. Also at large distances from the grid, the streamwise and transverse velocities have a significant correlation ($`C_{uv}0`$).

Since the streamwise fluctuations suffer from a shear, we are interested mainly in the transverse fluctuations. Their sub-Gaussian, Gaussian, and hyper-Gaussian PDFs are studied separately in the following subsections.

### IV.2 Sub-Gaussian PDF in developing turbulence

The transverse velocity has a sub-Gaussian PDF at the smallest distance from the grid, $`d=0.25\mathrm{m}`$ (Figs. 3 and 4). Fig. 5(a) shows energy spectra of the streamwise and transverse velocities. They have peaks at the wave numbers $`k7`$ and 14 m-1. The flow is in a transition state from quasi-periodic motions due to wavy wakes of the grid rods to weak turbulence. Fig. 5(b) shows the correlation coefficient between the Fourier transforms at adjacent wave numbers, $`C_{nn^{}}=a_na_n^{}/(a_n^2a_n^{}^2)^{1/2}`$ with $`n^{}=n+1`$. At around the energy peaks, we observe significant correlations. The Fourier transforms are not mutually independent.

The quasi-periodic motions should have finite correlation lengths. Namely, the motions should not be exactly periodic. If this were not the case, the central limit theorem would not be applicable to the calculation of the Fourier transforms (Sec. II). Then the transforms would not have a Gaussian PDF. This is inconsistent with our result (Fig. 2).

For velocity differences $`\delta u=u(x+\delta x)u(x)`$ and $`\delta v=v(x+\delta x)v(x)`$, where $`\delta x`$ is the separation in the mean-wind direction, Fig. 6 shows the variances $`\delta u^2`$ and $`\delta v^2`$, the flatness factors $`F_{\delta u}`$ and $`F_{\delta v}`$, the skewness factors $`S_{\delta u}`$ and $`S_{\delta v}`$, and the streamwise-transverse correlation $`C_{\delta u\delta v}`$. At small separations, we observe enhancements of $`F_{\delta u}`$, $`F_{\delta v}`$, $`S_{\delta u}`$ and $`C_{\delta u\delta v}`$. They are due to small-scale coherent structures NWLMF97 ; SA97 ; CG99 ; MOURI  and are not of our interest. At large separations, we observe oscillations of $`\delta u^2`$, $`\delta v^2`$, $`F_{\delta u}`$, $`F_{\delta v}`$, $`S_{\delta u}`$, $`S_{\delta v}`$, and $`C_{\delta u\delta v}`$. Their wavelengths roughly correspond to the wave number $`k7\mathrm{m}^1`$ of the spectral energy peak. The oscillations of $`\delta u^2`$ and $`\delta v^2`$ are in phase, while those of $`S_{\delta u}`$ and $`S_{\delta u}`$ are $`\pm 90^{}`$ out of phase and those of $`F_{\delta u}`$, $`F_{\delta v}`$, and $`C_{\delta u\delta v}`$ are $`180^{}`$ out of phase.

Our results are explained if the velocity field is a superposition of a few quasi-periodic motions and a random background. It is actually possible to reproduce the oscillations and their phases qualitatively with a few sinusoidal functions satisfying the solenoidal condition and a random Gaussian noise.

The quasi-periodic motions are predominant at local maxima of the variances $`\delta u^2`$ and $`\delta v^2`$. Since the flatness factor $`F_{\delta v}`$ is locally minimal and less than 3 there, the transverse-velocity amplitudes of those quasi-periodic motions lie in a limited range (see also Ref. J98 ). This discussion would apply to the streamwise velocity as well. The oscillations of $`S_{\delta u}`$ and $`S_{\delta v}`$ imply the presence of quasi-periodic motions with different wave numbers that are coupled with each other (see also Fig. 5(b)).

The background flow is predominant at local minima of the variances. Since the flatness factor $`F_{\delta u}`$ is greater than 3 there, the streamwise background tends to be intermittent. This tendency is not significant in the transverse background, which exhibits $`F_{\delta v}3`$ even at its local maxima. The slight enhancement of $`C_{\delta u\delta v}`$ implies that the streamwise and transverse components of the background flow tends to have a correlation.

Overall, the observed sub-Gaussianity of the transverse fluctuations is due to quasi-periodic motions that have finite amplitudes. The hyper-Gaussianity of the streamwise fluctuations is due to a background flow that tends to be intermittent.

### IV.3 Gaussian PDF in fully developed turbulence

The transverse velocity has a Gaussian PDF at intermediate distances from the grid, i.e., $`d1\mathrm{m}`$ in the small tunnel and $`d10\mathrm{m}`$ in the large tunnel (Figs. 3 and 4). This is because turbulence is fully developed there. Wavy wakes of the grid rods have evolved to energy-containing eddies that are random and independent. The corresponding Fourier transforms are thus independent. Actually we did not find spectral energy peaks such as those observed at $`d=0.25\mathrm{m}`$. The wave-number–wave-number correlation in the energy-containing range $`k10\mathrm{m}^1`$ is absent within the statistical uncertainty, $`|C_{nn^{}}|0.1`$.

The grid turbulence at the intermediate distances does not exhibit large-separation oscillations of $`\delta u^2`$, $`\delta v^2`$, $`F_{\delta u}`$, $`F_{\delta v}`$, $`S_{\delta u}`$, $`S_{\delta v}`$, and $`C_{\delta u\delta v}`$. The flatness factors $`F_{\delta u}`$ and $`F_{\delta v}`$ as well as the skewness factors $`S_{\delta u}`$ and $`S_{\delta v}`$ at large separations are identical to the Gaussian values. The streamwise-transverse correlation $`C_{\delta u\delta v}`$ at large separations is absent. These are observed usually in fully developed turbulence B53 ; VC69 ; FK67 ; MOURI .

### IV.4 Hyper-Gaussian PDF in decayed turbulence

The streamwise and transverse velocities have hyper-Gaussian PDFs at largest distances from the grid, $`d10\mathrm{m}`$. The streamwise-transverse correlation $`C_{uv}`$ is also enhanced there (Figs. 3 and 4). For the separations $`\delta x=0.20\mathrm{m}l_u`$ and $`\delta x=0.40\mathrm{m}2l_u`$, Fig. 7 shows the flatness factors $`F_{\delta u}`$ and $`F_{\delta v}`$ and the streamwise-transverse correlation $`C_{\delta u\delta v}`$. They increase with the distance. Since the grid turbulence has decayed (Fig. 1), there remain only strong energy-containing eddies B53 , which intermittently pass the probe. The enhancement of the flatness factors $`F_{\delta u}`$ and $`F_{\delta v}`$ is due to enhanced values of $`\delta u`$ and $`\delta v`$, which should be associated with the surviving strong eddies. The streamwise-transverse correlation $`C_{\delta u\delta v}`$ is enhanced if $`\delta u`$ and $`\delta v`$ are enhanced simultaneously at the positions of those eddies. In this case, velocity fluctuations $`u(x)`$ and $`v(x)`$ are also enhanced at the eddy positions, resulting in the hyper-Gaussian PDFs as well as the significant streamwise-transverse correlation $`C_{uv}`$.

The presence of large-scale spatial structures implies a correlation among the corresponding Fourier transforms. Since the spatial structures consist of many Fourier modes, the correlation is not local in the wave number space. Although we failed to detect any significant correlation, we found moderate correlations with $`|C_{nn^{}}|0.2`$ between many wave numbers in the energy-containing range $`k10\mathrm{m}^1`$.

Roughly at the position where the flatness factors $`F_u`$ and $`F_v`$ begin to differ from the Gaussian value, the skewness factor for the streamwise velocity $`S_u`$ changes its sign (Fig. 4). This is probably because turbulence becomes weak. The positive skewness is more significant at a smaller distance from the grid and is attributable to the turbulence itself. On the other hand, the negative skewness is more significant at a larger distance. Using a low-pass filtering technique, we ascertained the presence of long-wavelength motions along the mean-flow direction ($`k1\mathrm{m}^1`$). Their amplitude is larger and their PDF is more negatively skewed at a larger distance from the grid. These long-wavelength motions are attributable to effects of the wind tunnel, e.g., wall effects and deceleration of the mean flow.

### IV.5 Small scale statistics

Small-scale coherent structures such as vortex tubes in fully developed turbulence have attracted much interest NWLMF97 ; SA97 ; CG99 ; MOURI . The statistics that are studied most often are the flatness factors $`F_{u/x}`$ and $`F_{v/x}`$ and the skewness factor $`S_{u/x}`$ of the velocity derivatives $`u/x`$ and $`v/x`$. An increase of the microscale Reynolds number Re$`{}_{\lambda }{}^{}=u^2^{1/2}\lambda /\nu `$ is known to cause the increase of $`F_{u/x}`$ and $`F_{v/x}`$ and the decrease of $`S_{u/x}`$ SA97 . Here $`\nu `$ is the kinematic viscosity. We briefly summarize overall trends of Reλ, $`F_{u/x}`$, $`F_{v/x}`$, and $`S_{u/x}`$ over distances from the grid. They are not of our interest but are expected to be helpful in a future experimental research.

The microscale Reynolds number Reλ decreases with the distance as shown in Fig. 8(a). The flatness factors $`F_{u/x}`$ and $`F_{v/x}`$ increase as shown in Fig. 8(b). The skewness factor $`S_{u/x}`$ decreases as shown in Fig. 8(c). Thus the dependences of $`F_{u/x}`$, $`F_{v/x}`$, and $`S_{u/x}`$ on Reλ are apparently opposite to those mentioned above. This is because the turbulence state changes from developing to fully developed, and to decayed with an increase of the distance. In order to compare with other experimental data, the measurement is required to be done at a position where the grid turbulence is in the fully developed state.

## V CONCLUSION

The PDF of velocity fluctuations was studied systematically for grid turbulence. At small distances from the grid, where the turbulence is still developing, there are quasi-periodic motions having finite amplitudes, and hence the PDF is sub-Gaussian. At intermediate distances from the grid, where the turbulence is fully developed, motions of energy-containing eddies are random, and hence the PDF is Gaussian. At large distances from the grid, where the turbulence has decayed, there remain only strong eddies, and hence the PDF is hyper-Gaussian. The Fourier transforms of the velocity fluctuations always have Gaussian PDFs, in accordance with the central limit theorem. At intermediate distances from the grid, the Fourier transforms are statistically independent of each other. This is the necessary and sufficient condition for Gaussianity of the velocity fluctuations. At small and large distances, the Fourier transforms are dependent.

Our result serves as an example that the velocity fluctuations could have a sub-Gaussian PDF if there exist strong finite-amplitude motions. We suspect that this is the case in experiments of Noullez et al. NWLMF97  and Sreenivasan and Dhruva SD98 . They obtained $`F_v2.85`$ in a free air jet and $`F_u=2.66`$ in an atmospheric boundary layer, respectively (see also Refs. CG99 ; AGHA84 ), for which no explanation has been proposed yet. The signal could suffer from finite-amplitude motions generated by the jet nozzle or the surface. It is of interest to analyze such data in the same manner as in our present work.

###### Acknowledgements.

The authors are grateful to T. Umezawa for his help in our experiments and also to the referee for helpful comments.