text:
# The terahertz phonon laser: a full quantum treatment

## I Introduction

Since the pioneering experiments of Chang, Esaki and Tsu , the tunneling processes in semiconductor heterostructures have attracted great interest. This has been increasing with the development of modern deposition techniques such as metal-organic chemical-vapor deposition (MOCVD) and molecular beam epitaxy (MBE) that can provide artificial nanostructures. Such structures, coupled to electrical leads, can display many unusual behaviors.

Most of works have focalized on electron transport phenomena through heterostructure devices. The double barrier heterostructure (DBH) permits the creation of several kinds of electronic devices such as semiconductor laser diodes, ultra-high-frequency oscillators, and many others .

One of the most important problems in the study of resonant tunneling in DBH is the e-ph interaction. After the works of Goldman, Tsui and Cunningham , it was realized the importance of electron-phonon interaction on the electronic properties of these structures. However, little importance has been given to the study of the phonons generated in this process, the way they propagate, their decay processes, etc.

Our purpose is to study the phonon generation in a DBH under the effect of an external applied bias. Our device consists in a DBH made of GaAs-AlGaAs. It was designed in such a way that the difference $`\mathrm{\Delta }\epsilon `$ between the first excited level $`\epsilon _1`$ and the ground state $`\epsilon _0`$ in the well is of the order of the $`LO`$ phonon energy $`\mathrm{}\omega _1`$, for a small applied bias $`V`$. Here, $`\omega _1`$ is the $`\mathrm{\Gamma }`$-point $`LO`$ phonon frequency. The Fermi level $`\epsilon _F^L`$ in the emitter is such that the excited level is above it. With a further increase in $`V`$, the levels in the well are lowered relative to $`\epsilon _F^L`$. When the ground state comes under the bottom of the conduction band, the current is almost suppressed until the level $`\epsilon _1`$ reaches $`\epsilon _F^L`$. To continue increasing $`V`$, the current begins to flow through the excited level, but since $`\mathrm{\Delta }\epsilon `$ remains less than $`\mathrm{}\omega _1`$, the phonon emission is inhibited. For a given bias $`V`$ the resonant condition $`\mathrm{\Delta }\epsilon \mathrm{}\omega _1`$ is achieved and the electrons begin to decay to the ground state by emitting primary$`LO`$ phonons. The potential profile and the level positions at this resonant condition are shown in figure 1.

For an $`Al`$ concentration greater than $`0.25`$  or $`0.3`$ , these phonons are confined inside the well (they can be also absorbed by exciting electrons from $`\epsilon _0`$ to $`\epsilon _1`$).

The process described above, acts in parallel with the decay of primary $`LO`$ phonons due to anharmonicity. One product of this decay is a secondary longitudinal optical phonon ($`\stackrel{~}{LO}`$), the other is a transversal acoustic ($`TA`$) phonon . The $`\stackrel{~}{LO}TA`$ pair is produced by stimulated emission. Therefore these $`TA`$ phonons, for a bias greater than the operation threshold, are coherent and form a beam that we call saser by analogy with a laser. The coherence of the device presented here was studied previously  by developing a formalism similar to that employed usually to study the coherence in lasers  which take into account the competition among the several emitted modes.

In the early 60’s , generation of coherent phonons was discussed theoretically and observed experimentally . Several experiments and applications were suggested . These phonons can be generated by intense pulses of laser . More details were described in the review of Kurz . In a recent review, Merlin  discusses the generation of pulses of coherent $`THz`$ optical phonons by intense femtosecond-pulsed lasers.

Coherent acoustic phonons can be also obtained with a wavelength given by the period of the superlattice where phonons are optically generated (i.e., $`500`$ Å) . In our case, the system generates coherent sound by electronic means with much shorter wavelength. The generation of high-frequency monochromatic acoustic waves by laser-induced thermomodulation was reported by Damen . The phonon frequency ranged from $`2`$ to $`4GHz`$ comparing with $`2THz`$ for our device.

Generation of phonons by stimulated emission was theoretically considered  and observed experimentally. Bron and Grill  reported the direct observation of the stimulated emission of $`24.7cm^1`$ phonons in a three level system of $`V^{4+}`$ ions in $`Al_2O_3`$. Prieur et al.  studied the phonon emission in a two level system of a glass. These phonons had a wide frequency distribution around $`\omega =0.34GHz`$. Zavtrak and collaborators  proposed a device consisting in a dielectric liquid with small particles or gas bubbles working in such a way similar to a free-electron laser. That device could produce coherent acoustic phonons by stimulated emission at a frequency $`\omega 2KHz`$. Stimulated emission of phonons in an acoustic cavity was obtained experimentally by Fokker et al. . They use the metastable Zeeman-split doublet in ruby. These phonons have frequencies of about $`60MHz`$.

At present, several kinds of amplified coherent beams of bosons are produced by stimulated emission. Besides the traditional visible or near infrared photons (lasers) and microwave photons (masers), other coherent photon beams are produced. High energy lasing processes like X ray lasers (called also xasers)  and gamma ray lasers (GRASERS)  were studied.

A very exciting field has been recently opened with the development of atoms lasers . Several applications of these novel devices, like high precision gravimetry, are proposed. Strongly related, it is the study of vibrational amplification by stimulated emission of radiation (vaser) .

Kleppner  claims against that complicated jargon of acronyms and proposes to simplify to atom laser, X ray laser, etc.

In summary, a great effort has been made to produce terahertz phonons using other mechanisms different of stimulated emission. On the other side, acoustic phonons have been produced by stimulated emission but at frequencies lower than $`2THz`$. The DBH resonator proposed here could produce a continuous beam of terahertz coherent phonons. The ultra short wavelength of these phonons permits potential applications.

Acoustic imaging is a very well known procedure, mainly in medicine, where it is exploited the fact that different water contents on several tissues produce different attenuations for a signal of about $`4MHz`$. Today are available commercial equipments that can produce three-dimensional images of very high resolution.

Equipments for acoustic microscopy are less known but they are also commercially available. The most used technique is Scanning Acoustic Microscopy (SAM)  that can attain resolutions higher than $`0.4\mu m`$, taken with phonons of $`2GHz`$.

Besides the scanning acoustic microscopy, several techniques have been presented in a recent review : photoacoustic microscopy (PAM), scanning electron acoustic microscopy (SEAM) and scanning laser acoustic microscopy (SLAM). The theory of acoustic imaging was also reviewed, schematic diagrams of the different microscopes are presented and several industrial applications are discussed in reference .

The equipments mentioned above work at up to $`10GHz`$, three orders of magnitude lower than the frequency of the phonon beam considered in the present paper. Therefore, they have a limited theoretical resolution given by its wavelength that is $`1000`$ times greater than the generated by a DBH device.

To make the acoustic nanoscopy possible, it is necessary to have a coherent beam of ultra short wavelength phonons. To obtain images of a small structure, one would need waves shorter than the size of the details of the system.

H. Maris  describes the process of acoustic nanoscopy by mean of echoes of ultra-short phonon pulses optically generated. They are detected by measuring the variation of optical reflectivity due to the scattered phonon pulses. They work with coherent phonons of about $`2THz`$ like those produced by our device, but they were not produced by a phonon laser.

Evidently, a laser is not necessary to produce optical images. Nevertheless, lasers are a very useful tool for imaging and up today they are indispensable to make holograms. On the other hand, many optical experiments that can be made with incoherent light are done more easily and more precisely with the aid of a laser. Therefore many experiments that are made today with incoherent phonons could be done better with the aid of a phonon laser.

The phonon laser described in this paper has the wavelength smaller than $`25`$ Å (what is about 20 times less than for the device described in ). Therefore, it is possible to use it to see structures of about $`50nm`$. A coherent source of X rays of about $`20nm`$ begins to be available now . Even if it will be possible to make images of structures at the submicron scale by using X rays of a wavelength of about $`25`$ Å, they would have energies of about $`0.5KeV`$, what could affect the system to be studied.

The device presented here could be used to get a hologram at nano-scale. To do this, it is possible to build in a GaAs an acoustic analogous of a Michelson interferometer. It consists of a thin barrier of AlGaAs designed as a semimirror, placed at $`45^o`$ to split the beam and to force it to pass one half through the part of the sample we are interested in and the other across a very pure reference GaAs. After that, we can use mirrors (made of thick AlGaAs barriers) to superpose them and record the interference pattern with a detector. This pattern could be processed in a computer to get a hologram of the sample showing the three-dimensional shape of the nanostructure. In a previous paper  we described, as an alternative, the acoustic version of a Sagnac interferometer.

Because of the high sensitivity to all types of lattice defects, the phonon beam could be used to give information about the inner lattice structure ,  and form part of a phonon imaging setup as a source .

We will show that the amplitude of the phonon beam is with a great precision directly proportional to the incoming current. Therefore, the saser amplitude could be modulated by modulating the amplitude of the current. This will permit the use of the phonon laser as an information carrier between circuit components at very short distances. These “phonoelectronic” devices could work at smaller distances and lower energies than optoelectronic devices.

If we achieve to make a “phonoresist”, i.e., a material with a soft mode close to the frequency of our phonon laser, the saser could be used for nanolithography. Due to the short wavelength and low energy of the phonon beam, “phonolithography” would be a better option than the methods currently in used.

The rest of this paper is organized as follow. The Hamiltonian which describes our system is presented in Sec. II. The solutions of the Hamiltonian are obtained is Sec. III. In the Sec. IV, we show the procedure used to obtain the kinetic equations. We devote Sec. V to present our results and conclusions.

## II The Hamiltonian

To describe our system, we considered a Hamiltonian composed by a single particle part and one representing the interactions. The single particle part can be written as

$$_o=_e+_1+_2+_3,$$
(1)
where $`_e`$ describes the electrons, and $`_1`$, $`_2`$, $`_3`$ the $`LO`$, $`\stackrel{~}{LO}`$ and $`TA`$ phonons respectively.

Assuming the grown direction as the $`z`$ axis, the Hamiltonians $`_1`$, $`_2`$, $`_3`$ for the phonons can be written

$$_1=\underset{𝐪_1}{}\mathrm{}\omega _1b_{1_{𝐪_1}}^{}b_{1_{𝐪_1}},$$
(2)
$$_2=\underset{𝐪_2}{}\left(\mathrm{}\omega _2i\mathrm{}\kappa _2\right)b_{2_{𝐪_2}}^{}b_{2_{𝐪_2}},$$
(3)
$$_3=\underset{𝐪_3}{}\left(\mathrm{}\omega _3i\mathrm{}\kappa _3\right)b_{3_{𝐪_3}}^{}b_{3_{𝐪_3}}$$
(4)
where $`\mathrm{}\omega _1`$, $`\mathrm{}\omega _2`$, $`\mathrm{}\omega _3`$ are the energies for $`LO`$, $`\stackrel{~}{LO}`$ and $`TA`$ phonons respectively and $`b_{i_{𝐪_i}}^{}`$($`b_{i_{𝐪_i}}`$) are the creation (annihilation) operators for phonons with momentum $`𝐪_i`$ perpendicular to the $`z`$ direction. Two imaginary terms $`\left(i\mathrm{}\kappa _2\text{}i\mathrm{}\kappa _3\right)`$ were introduced to take into account the decay by anharmonicity of the $`\stackrel{~}{LO}`$ phonons and the escape of the $`TA`$ phonons, respectively.

To describe the electrons, we use a tight binding approach with hoppings $`v`$ between nearest neighbors. Writing the Hamiltonian in a Wannier basis, in terms of the electronic creation (annihilation) operators $`c_𝐫^{}`$($`c_𝐫`$) at site $`𝐫`$, we get

$$_e=\underset{𝐫}{}\epsilon _𝐫c_𝐫^{}c_𝐫+v\underset{<𝐫,𝐫^{}>}{}\left(c_𝐫^{}c_𝐫^{}+c_𝐫^{}^{}c_𝐫\right)$$
(5)
where $`𝐫(l,m,j)`$ with $`l,`$ $`m,`$ $`j=\mathrm{},`$ $`\mathrm{},`$ $`\mathrm{}`$. As the system has not magnetic behavior the subscript $`\sigma `$ (for spin) is left implicit.

Our system has translational symmetry in the perpendicular directions ($`xy`$). Therefore the Hamiltonian can be uncoupled by expanding the operators $`c_𝐫`$ in plane waves in these directions

$$c_𝐫=\underset{𝐤}{}c_{j𝐤}e^{i𝐤𝐱_{lm}}.$$
(6)
We can treat the system as a sum over 1D Hamiltonians for each wave vector $`𝐤`$ perpendicular to the current $`\left(z\right)`$ direction

$$_e=\underset{𝐤}{}\left\{\underset{j}{}\stackrel{~}{\epsilon }_{j𝐤}c_{j𝐤}^{}c_{j𝐤}+v\underset{<jj^{}>}{}\left(c_{j𝐤}^{}c_{j^{}𝐤}+c_{j^{}𝐤}^{}c_{j𝐤}\right)\right\}$$
(7)
where $`\stackrel{~}{\epsilon }_{j𝐤}=\epsilon _{j𝐤}4v`$. The energies $`\stackrel{~}{\epsilon }_{j𝐤}`$ are measured from the bottom of the conduction band of the emitter and $`\epsilon _{j𝐤}=\epsilon _j+\epsilon _𝐤`$. The energies $`\epsilon _j`$ are chosen to describe the energy profile of the DBH. For the sake of simplicity, we will leave implicit the dependence on $`𝐤`$.

For the $`z`$ direction, we separate the space in three regions: the dispersion region and two semi-infinite one-dimensional chains. On the left we will have planes with energy $`\epsilon _j=0`$ for $`j0`$ and to the right $`\epsilon _j=V`$ for $`jL+1`$, $`L`$ is the length of the DBH. The corresponding eigenstates of these two regions are planes waves. If we disconnect the DBH from the left and right chains we get for the dispersion region the profile of an infinite (even not rectangular) well as it is shown in figure 2.

For the dispersion region, we diagonalize a three-diagonal matrix of order $`L`$ corresponding to the profile showed above getting the eigenvalues $`\epsilon _m`$ and the eigenvectors $`|m`$ of equation

$$_{ez}^{}|m=\epsilon _m|m.$$
(8)
Here, $`_{ez}^{}`$ is the part of $`_{ez}`$ that goes from the beginning of the left barrier to the end of the right one. Written in the basis of planes this is

$$_{ez}^{}=\left(\begin{array}{cccccc}\stackrel{~}{\epsilon }_1& v& 0& & & \\ v& \stackrel{~}{\epsilon }_2& v& & & \\ 0& v& \stackrel{~}{\epsilon }_3& \mathrm{}& & \\ & & \mathrm{}& \mathrm{}& & \\ & & & & & v\\ & & & & v& \stackrel{~}{\epsilon }_L\end{array}\right),$$
(9)
where $`\stackrel{~}{\epsilon }_j=\epsilon _j2v`$.

After the diagonalization we get $`L`$ discrete levels. We label these levels through the index $`m`$ as $`m=0`$, $`1^{}`$, $`2^{}\mathrm{}`$. For the right chain we rename the planes, from $`L+1`$, $`L+2`$, $`\mathrm{}`$ to $`1`$, $`2`$, $`\mathrm{}`$; and for the left one, from $`0`$, $`1`$, $`2`$,$`\mathrm{}`$, $`\mathrm{}`$ to $`\overline{1}`$, $`\overline{2}`$, $`\mathrm{}`$, $`\overline{\mathrm{}}`$. Therefore we get a new picture where we have two semi-infinite chains, the first one corresponding to the planes $`j=\overline{\mathrm{}}`$, $`\mathrm{}`$, $`\overline{1}`$ and the other one corresponding to $`j=1`$, $`\mathrm{}`$, $`\mathrm{}`$. Among the $`L`$ levels obtained from the diagonalization of the matrix (9), only the two states with lowest energies will participate significantly in the electronic transport. The others are far over the Fermi level and their effects will be considered later.

To connect the DBH with the left and right chains we compute the matrix elements

$$v_{jm}=j|_{ez}|m.$$
(10)
The electronic Hamiltonian is

$`_e`$ $`=`$ $`{\displaystyle \underset{j}{}}\epsilon _jc_j^{}c_j+{\displaystyle \underset{m}{}}\epsilon _mc_m^{}c_m+{\displaystyle \underset{j\overline{1},m}{}}v\left(c_j^{}c_{j+1}+c_{j+1}^{}c_j\right)+`$ (12)
$`{\displaystyle \underset{m}{}}\left\{v_{\overline{1}m}\left(c_{\overline{1}}^{}c_m+c_m^{}c_{\overline{1}}\right)+v_{m1}\left(c_m^{}c_1+c_1^{}c_m\right)\right\}.`$
A diagram representing the electronic Hamiltonian (12) is shown in figure 3.

As it is well known , the dominant dispersion process in polar semiconductors is due to the coupling between electrons and $`LO`$ phonons.

The interaction Hamiltonian $`_{int}`$ has three terms

$$_{int}=_{ee}+_{eph}+_{phph}$$
(13)
where $`_{ee}`$, $`_{eph}`$ and $`_{phph}`$ describe the electron-electron, the electron-phonon, and the phonon-phonon interactions respectively.

$`_{ee}`$ is taken in the Hartree approximation. That means that the energies $`\epsilon _m`$ of the levels in the well depend on the accumulated charge.

The electron-phonon interaction Hamiltonian has the form,

$$_{eph}=\underset{𝐤,𝐪_1}{}g_{_{\mathrm{𝐤𝐪}_1}}(c_{0_{𝐤𝐪_1}}^{}c_{1^{}𝐤}b_{1_{𝐪_1}}^{}+c_{1^{}𝐤}^{}c_{0_{𝐤𝐪_1}}b_{1_{𝐪_1}}).$$
(14)
In order to recover the translational symmetry in the $`xy`$ plane we take the average on $`𝐪_1`$

$$_{eph}=\underset{𝐤}{}g__𝐤\left(c_{0_𝐤}^{}c_{1^{}𝐤}b_{1_𝐤}^{}+c_{1^{}𝐤}^{}c_{0_𝐤}b_{1_𝐤}\right)$$
(15)
where $`g_𝐤`$ measure the effective strength of the interaction and was calculated as in .

The phonon-phonon interaction is given by

$$_{phph}=\underset{𝐪_1}{}\gamma _{𝐪_1}(b_{1_{𝐪_1}}b_{2_{𝐪_1}}^{}b_3^{}+b_{1_{𝐪_1}}^{}b_{2_{𝐪_1}}b_3).$$
(16)
This Hamiltonian describes the decay $`LO\stackrel{~}{LO}+TA`$ and its inverse process (recombination).

As in the GaAs the $`TA`$ phonons are emitted fundamentally in the \[$`111`$\] direction , we considered here a device that was grown in this direction. The inside barrier walls act as mirrors of a Fabry-Perot interferometer selecting discrete cavity modes. This mirrors also serve for a selection of the saser modes: those running in the grown direction are reflected several times between the mirrors and stay longer in the well, while all others are lost.

With the increasing of the applied bias $`V`$, from a certain value, called saser threshold, the $`TA`$ phonons associated with the cavity resonant mode, i.e., $`q_3=0`$, stimulate strongly the $`LO`$ phonon decay and their number grow sharply. Once this mode is sufficiently populated, the probability of stimulated emission for that mode exceeds the sum of the emission probabilities for all other modes. Thus, the $`TA`$ phonons emitted in this mode will slave the other ones, giving rise to a coherent acoustic phonon beam. This process is similar to what happens in a laser. Due to the momentum conservation $`𝐪_1=𝐪_2+𝐪_3`$ and because of $`𝐪_3=0`$  we get $`𝐪_2=𝐪_1`$.

The parameter $`\gamma _{𝐪_1}`$ appearing in (16) represents the intensity of the interaction and was estimated from ,  and .

Finally, the total Hamiltonian is the sum of both terms: the single particle and the interaction parts

$$=\underset{𝐤}{}_o\left(𝐤\right)+_{int}\left(𝐤\right)$$
(17)
where we had made the approximation that the momentum $`𝐪_1`$ (appearing in expressions (2), (3) and (16)) varies in the same range as $`𝐤`$. We can see at expression (17) that the Hamiltonian became, in this approximation, fully uncoupled for each $`𝐤`$. This simplifies strongly the obtention of its solutions.

## III The resolution of the Hamiltonian

In order to solve de Hamiltonian described in the former section we introduce the operators

$$𝒪_{j𝐤}^{n_1,n_2,n_3}=c_{j𝐤}\frac{b_{1𝐤}^{n_1}b_{2𝐤}^{n_2}b_3^{n_3}}{\sqrt{n_1!}\sqrt{n_2!}\sqrt{n_3!}}.$$
(18)
The operators $`𝒪_{j𝐤}^{n_1,n_2,n_{3}^{}{}_{}{}^{}}`$ creates the eigenstates of $`_o`$ (without the hopping terms)

$$|j𝐤n_1n_2n_3𝒪_{j𝐤}^{n_1,n_2,n_{3}^{}{}_{}{}^{}}|0$$
(19)
formed by an electron at plane $`j`$ and $`n_1`$ $`LO`$ phonons, $`n_2`$ $`\stackrel{~}{LO}`$ phonons with parallel momentum $`𝐤`$; and $`n_3`$ $`TA`$ phonons.

As it was done before, we will leave implicit the dependence on $`𝐤`$.

Using the operators (18) and the Hamiltonian (17) we can calculate their equations of motion. In order to simplify the notation we call $`𝐧`$ for the set of $`\{n_1,n_2,n_3\}`$ phonons

$$i\mathrm{}\frac{d𝒪_j^𝐧}{dt}=[𝒪_j^𝐧,].$$
(20)
By expanding the eigenstates $`|\mathrm{\Phi }`$ of $``$ in the eigenstates of $`_o`$

$$|\mathrm{\Phi }=\underset{j,𝐧}{}p_j^𝐧|j𝐧,$$
(21)
and due to the orthogonality of $`|j𝐧,`$ the amplitudes can be calculated as

$$p_j^𝐧=jn|\mathrm{\Phi }=0|𝒪_j^𝐧|\mathrm{\Phi }.$$
(22)
Making use of the equations turning out from (20) and the previous definition for the amplitudes, we get for $`j\overline{1},1,m\left(m=0,1^{},2^{},3^{}\mathrm{}L\right)`$

$$i\mathrm{}\frac{dp_j^𝐧}{dt}=\epsilon _j^𝐧p_j^𝐧+\mathrm{{\rm Y}}_j+vp_{j+1}^𝐧+vp_{j1}^𝐧$$
(23)
and for $`j=\overline{1},1,m`$

$`i\mathrm{}{\displaystyle \frac{dp_{\overline{1}}^𝐧}{dt}}`$ $`=`$ $`\epsilon _{\overline{1}}^𝐧p_{\overline{1}}^𝐧+\mathrm{{\rm Y}}_{\overline{1}}+v_{\overline{1}1^{}}p_1^{}^𝐧+v_{\overline{1}0}p_0^𝐧+vp_{\overline{2}}^𝐧+{\displaystyle \underset{m2^{}}{}}v_{\overline{1}m}p_m^𝐧,`$ (25)
$`i\mathrm{}{\displaystyle \frac{dp_0^𝐧}{dt}}`$ $`=`$ $`\epsilon _0^𝐧p_0^𝐧+g\sqrt{n_1}p_1^{}^𝐧^{}+\mathrm{{\rm Y}}_0+v_{\overline{1}0}p_{\overline{1}}^𝐧+v_{01}p_1^𝐧,`$ (26)
$`i\mathrm{}{\displaystyle \frac{dp_1^{}^𝐧}{dt}}`$ $`=`$ $`\epsilon _1^{}^𝐧p_1^{}^𝐧+g\sqrt{n_1+1}p_0^{𝐧^+}+\mathrm{{\rm Y}}_1^{}+v_{\overline{1}1^{}}p_{\overline{1}}^𝐧+v_{1^{}1}p_1^𝐧,`$ (27)
$`i\mathrm{}{\displaystyle \frac{dp_1^𝐧}{dt}}`$ $`=`$ $`\epsilon _1^𝐧p_1^𝐧+\mathrm{{\rm Y}}_1+v_{1^{}1}p_1^{}^𝐧+v_{01}p_0^𝐧+vp_2^𝐧+{\displaystyle \underset{m2^{}}{}}v_{1m}p_m^𝐧,`$ (28)
$`i\mathrm{}{\displaystyle \frac{dp_m^𝐧}{dt}}`$ $`=`$ $`\epsilon _m^𝐧p_m^𝐧+v_{1m}p_1^𝐧+v_{\overline{1}m}p_{\overline{1}}^𝐧+\mathrm{{\rm Y}}_m,m2^{}`$ (29)
where

$`\epsilon _j^𝐧`$ $`=`$ $`\epsilon _j+n_1\mathrm{}\omega _1+n_2\left(\mathrm{}\omega _2i\mathrm{}\kappa _2\right)+n_3\left(\mathrm{}\omega _3i\mathrm{}\kappa _3\right),`$
$`𝐧^{}`$ $``$ $`\{n_11,n_2,n_3\},`$
$`𝐧^+`$ $``$ $`\{n_1+1,n_2,n_3\},`$
$`\mathrm{{\rm Y}}_j`$ $`=`$ $`\gamma \left\{\sqrt{N^{}}p_j^𝐧^{}+\sqrt{N^+}p_j^{𝐧^{\prime \prime }}\right\},`$
$`𝐧^{}`$ $``$ $`\{n_11,n_2+1,n_3+1\},`$
$`𝐧^{\prime \prime }`$ $``$ $`\{n_1+1,n_21,n_31\},`$
$`N^{}`$ $`=`$ $`n_1\left(n_2+1\right)\left(n_3+1\right),`$
$`N^+`$ $`=`$ $`\left(n_1+1\right)n_2n_3.`$
The amplitudes $`p_j^𝐧`$ are time dependent. We transform them

$$p_j^𝐧\left(t\right)=\frac{1}{\sqrt{2\pi }}_{\mathrm{}}^{\mathrm{}}𝑑\omega a_j^𝐧\left(\omega \right)e^{i\omega t}.$$
(30)
By doing that, we get the following set of equations

$`\mathrm{}\omega a_j^𝐧=\epsilon _j^𝐧a_j^𝐧+\mathrm{\Gamma }_j+va_{j+1}^𝐧+va_{j1}^𝐧,\text{for }j\overline{1},1,m`$ (32)
$`\mathrm{}\omega a_{\overline{1}}^𝐧=\epsilon _{\overline{1}}^𝐧a_{\overline{1}}^𝐧+\mathrm{\Gamma }_{\overline{1}}+va_{\overline{2}}^𝐧+v_{\overline{1}1^{}}a_1^{}^𝐧+v_{\overline{1}0}a_0^𝐧+{\displaystyle \underset{m2^{}}{}}v_{\overline{1}m}a_m^𝐧,`$ (33)
$`\mathrm{}\omega a_0^𝐧=\epsilon _0^𝐧a_0^𝐧+g\sqrt{n_1}a_1^{}^𝐧^{}+\mathrm{\Gamma }_0+v_{\overline{1}0}a_{\overline{1}}^𝐧+v_{01}a_1^𝐧,`$ (34)
$`\mathrm{}\omega a_1^{}^𝐧=\epsilon _1^{}^𝐧a_1^{}^𝐧+g\sqrt{n_1+1}a_0^{𝐧^+}+\mathrm{\Gamma }_1^{}+v_{\overline{1}1^{}}a_{\overline{1}}^𝐧+v_{1^{}1}a_1^𝐧,`$ (35)
$`\mathrm{}\omega a_1^𝐧=\epsilon _1^𝐧a_1^𝐧+\mathrm{\Gamma }_1+va_2^𝐧+v_{1^{}1}a_1^{}^𝐧+v_{01}a_0^𝐧+{\displaystyle \underset{m2^{}}{}}v_{1m}a_m^𝐧,`$ (36)
$`\mathrm{}\omega a_m^𝐧=\epsilon _m^𝐧a_m^𝐧+\mathrm{\Gamma }_m+v_{\overline{1}m}a_{\overline{1}}^𝐧+v_{1m}a_1^𝐧,m2^{}.`$ (37)
Here $`\mathrm{\Gamma }_j`$ is the Fourier transform of $`\mathrm{{\rm Y}}_j`$.

In order to obtain the solutions of (III), we considered only three channels through which the electrons can tunnel. Such channels are the only that are between the bottom of the conduction band ($`\epsilon =0`$) and the Fermi energy $`\epsilon _F^L`$ at the emitter, near the resonant condition. Therefore, they are the only ones that contribute significantly to the current. These channels are: one electron through the excited level and no phonons (channel $`0`$, i.e., $`𝐧=\{0,0,0\}`$); one electron through the ground level and one $`LO`$ phonon (channel $`100`$) and one electron through the ground level and a pair $`\stackrel{~}{LO}TA`$ phonons respectively (channel $`011`$).

To solve the system above, we drop the imaginary part of $`\epsilon _j^𝐧`$. This cope with the problems caused by the non-Hermitian Hamiltonian.

We have to note that the amplitudes $`a_j^0`$ (for $`j1`$ or $`j\overline{1}`$) remain uncoupled (this tunneling channel does not mix with any other). On the other hand, the other two channels are coupled due to the ph-ph interaction. To uncouple them we propose the following transformation for the amplitudes

$`a_j^{100}`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{2}}}\left(a_j^++a_j^{}\right),`$ (38)
$`a_j^{011}`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{2}}}\left(a_j^+a_j^{}\right).`$ (39)
Making the above transformations for the amplitudes $`a_j^{100}`$ and $`a_j^{011}`$ we get two new uncoupled channels (called channels $`+`$ and $``$).

Considering only the channels $`0`$ and $`\pm `$, the system (III) can be rewritten as

$`va_{j1}^0+\left(\epsilon _j^0\mathrm{}\omega \right)a_j^0+va_{j+1}^0=0`$ (41)
$`va_{j1}^\pm +\left(\epsilon _j^{100}\pm \gamma \mathrm{}\omega \right)a_j^\pm +va_{j+1}^\pm =0\text{for }j\overline{1},m,1;`$ (42)
$`va_{\overline{2}}^0+\left(\epsilon _{\overline{1}}^0\mathrm{}\omega \right)a_{\overline{1}}^0+v_{\overline{1}1^{}}a_1^{}^0+{\displaystyle \underset{m1^{}}{}}v_{\overline{1}m}a_m^0=0`$ (43)
$`va_{\overline{2}}^\pm +\left(\epsilon _{\overline{1}}^{100}\pm \gamma \mathrm{}\omega \right)a_{\overline{1}}^\pm +v_{\overline{1}0}a_0^\pm +{\displaystyle \underset{m0}{}}v_{\overline{1}m}a_m^\pm =0`$ (44)
$`v_{1^{}1}a_1^{}^0+{\displaystyle \underset{m1^{}}{}}v_{m1}a_m^0+\left(\epsilon _1^0\mathrm{}\omega \right)a_1^0+va_2^0=0`$ (45)
$`v_{01}a_0^\pm +{\displaystyle \underset{m1^{}}{}}v_{m1}a_m^\pm +\left(\epsilon _1^{100}\pm \gamma \mathrm{}\omega \right)a_1^\pm +va_2^\pm =0`$ (46)
$`v_{\overline{1}m}a_{\overline{1}}^\pm +\left(\epsilon _m^{100}\pm \gamma \mathrm{}\omega \right)a_m^\pm +v_{m1}a_1^\pm =0\text{for }m0;`$ (47)
$`v_{\overline{1}m}a_{\overline{1}}^0+\left(\epsilon _m^0\mathrm{}\omega \right)a_m^0+v_{m1}a_1^0=0\text{for }m1^{};`$ (48)
$`v_{\overline{1}0}a_{\overline{1}}^\pm +\left(\epsilon _0^{100}\pm \gamma \mathrm{}\omega \right)a_0^\pm +{\displaystyle \frac{g}{\sqrt{2}}}a_1^{}^0+v_{01}a_1^\pm =0`$ (49)
$`v_{\overline{1}1^{}}a_{\overline{1}}^0+\left(\epsilon _1^{}^0\mathrm{}\omega \right)a_1^{}^0+{\displaystyle \frac{g}{\sqrt{2}}}a_0^++{\displaystyle \frac{g}{\sqrt{2}}}a_0^{}+v_{1^{}1}a_1^0=0`$ (50)
The amplitudes $`a_m^\pm `$ and $`a_m^0`$ appearing in equations (III) take into account the electronic energy levels different from the fundamental and excited one. To eliminate these amplitudes (the fundamental and excited levels are the only ones that participate directly in the $`LO`$ phonon generation) we use the equations (43)-(48). The resultant system is

$`va_{j1}^0+\left(\epsilon _j^0\mathrm{}\omega \right)a_j^0+va_{j+1}^0=0`$ (52)
$`va_{j1}^\pm +\left(\epsilon _j^{100}\pm \gamma \mathrm{}\omega \right)a_j^\pm +va_{j+1}^\pm =0\text{for }j\overline{1},m,1;`$ (53)
$`va_{\overline{2}}^0+\left(\stackrel{~}{\epsilon }_{\overline{1}}^0\mathrm{}\omega \right)a_{\overline{1}}^0+v_{\overline{1}1^{}}a_1^{}^0+\stackrel{~}{v}_{\overline{1}1}^0a_1^0=0`$ (54)
$`va_{\overline{2}}^\pm +\left(\stackrel{~}{\epsilon }_{\overline{1}}^\pm \pm \gamma \mathrm{}\omega \right)a_{\overline{1}}^\pm +v_{\overline{1}0}a_0^\pm +\stackrel{~}{v}_{\overline{1}1}^\pm a_1^\pm =0`$ (55)
$`\stackrel{~}{v}_{\overline{1}1}^0a_{\overline{1}}^0+v_{1^{}1}a_1^{}^0+\left(\stackrel{~}{\epsilon }_1^0\mathrm{}\omega \right)a_1^0+va_2^0=0`$ (56)
$`\stackrel{~}{v}_{\overline{1}1}^\pm a_{\overline{1}}^\pm +v_{01}a_0^\pm +\left(\stackrel{~}{\epsilon }_1^\pm \pm \gamma \mathrm{}\omega \right)a_1^\pm +va_2^\pm =0`$ (57)
$`v_{\overline{1}1^{}}a_{\overline{1}}^0+\left(\epsilon _1^{}^0\mathrm{}\omega \right)a_1^{}^0+{\displaystyle \frac{g}{\sqrt{2}}}a_0^++{\displaystyle \frac{g}{\sqrt{2}}}a_0^{}+v_{1^{}1}a_1^0=0`$ (58)
$`v_{\overline{1}0}a_{\overline{1}}^\pm +\left(\epsilon _0^{100}\pm \gamma \mathrm{}\omega \right)a_0^\pm +{\displaystyle \frac{g}{\sqrt{2}}}a_1^{}^0+v_{01}a_1^\pm =0`$ (59)
where

$$\stackrel{~}{v}_{\overline{1}1}^0=\underset{m1^{}}{}\frac{v_{\overline{1}m}v_{m1}}{\left(\mathrm{}\omega \epsilon _m^0\right)}$$
(60)
and

$$\stackrel{~}{v}_{\overline{1}1}^\pm =\underset{m0}{}\frac{v_{\overline{1}m}v_{m1}}{\left(\mathrm{}\omega \epsilon _m^{100}\gamma \right)}$$
(61)
are renormalized hoppings joining sites $`j=\overline{1},1`$, and

$`\stackrel{~}{\epsilon }_{\overline{1}}^0`$ $`=`$ $`\epsilon _{\overline{1}}^0+\mathrm{\Sigma }_{\overline{1}}^0`$ (63)
$`\stackrel{~}{\epsilon }_{\overline{1}}^\pm `$ $`=`$ $`\epsilon _{\overline{1}}^{100}+\mathrm{\Sigma }_{\overline{1}}^\pm `$ (64)
$`\stackrel{~}{\epsilon }_1^0`$ $`=`$ $`\epsilon _1^0+\mathrm{\Sigma }_1^0`$ (65)
$`\stackrel{~}{\epsilon }_1^\pm `$ $`=`$ $`\epsilon _1^{100}+\mathrm{\Sigma }_1^\pm `$ (66)
are renormalized energies at the sites $`j=\overline{1},1`$ for each channel, where

$`\mathrm{\Sigma }_{\overline{1}}^0`$ $`=`$ $`{\displaystyle \underset{m1^{}}{}}{\displaystyle \frac{v_{\overline{1}m}^2}{\left(\mathrm{}\omega \epsilon _m^0\right)}}`$ (68)
$`\mathrm{\Sigma }_{\overline{1}}^\pm `$ $`=`$ $`{\displaystyle \underset{m0}{}}{\displaystyle \frac{v_{\overline{1}m}^2}{\left(\mathrm{}\omega \epsilon _m^{100}\gamma \right)}}`$ (69)
$`\mathrm{\Sigma }_1^0`$ $`=`$ $`{\displaystyle \underset{m1^{}}{}}{\displaystyle \frac{v_{m1}^2}{\left(\mathrm{}\omega \epsilon _m^0\right)}}`$ (70)
$`\mathrm{\Sigma }_1^\pm `$ $`=`$ $`{\displaystyle \underset{m0}{}}{\displaystyle \frac{v_{m1}^2}{\left(\mathrm{}\omega \epsilon _m^{100}\gamma \right)}}.`$ (71)
A diagram representing equations (III) is shown in figure 4.

To solve the system (III) we proceed as follows.

For the channel $`0`$ and $`\pm `$, in the region $`j1`$ the solutions of (52) and (53) are respectively

$`a_j^0`$ $`=`$ $`A_T^0e^{ik_0^{^{}}z_j}`$ (72)
$`a_j^\pm `$ $`=`$ $`A_T^\pm e^{ik^^\pm z_j}`$ (73)
where $`z_j`$ is the coordinate of the plane $`j`$ ($`z_j=ja`$; $`a`$ is the lattice parameter), $`A_T^0`$, $`A_T^\pm `$ are the amplitudes of the transmitted waves, and $`k_0^{^{}}`$ and $`k^^\pm `$ are the $`z`$ component of the waves vectors that fulfills the dispersion relation

$`\mathrm{}\omega `$ $`=`$ $`\epsilon _1^0+2v\mathrm{cos}\left(k_0^{^{}}a\right),`$ (74)
$`\mathrm{}\omega `$ $`=`$ $`\epsilon _1^{100}+2v\mathrm{cos}\left(k^\pm a\right)\pm \gamma `$ (75)
We did not consider the terms $`e^{ik_𝐧^{^{}}z_j}`$ because we assume that there are no electrons coming in from the right of the DBH.

In the region $`j1`$ we have two kinds of solutions. The first one is obtained for the excited level that is above the bottom of the conduction band, in this case the solution is

$$a_j^0=A_Re^{ik_0z_j}+A_Ie^{ik_0z_j}$$
(76)
where $`A_R`$ and $`A_I`$ are the amplitudes of the reflected and incident waves and $`k_0`$ satisfies

$$\mathrm{}\omega =\epsilon _{\overline{1}}^0+2v\mathrm{cos}\left(k_0a\right).$$
(77)
The second kind of solution is obtained for the ground state that is below the bottom of the conduction band then, in this case, the solutions are evanescent modes with amplitudes $`A_E^+`$ and $`A_E^{}`$

$$a_j^\pm =A_E^\pm e^{\kappa ^^\pm z_j}$$
(78)
where the $`\kappa ^\pm `$ are positive and fulfill the relations

$$\mathrm{}\omega =\epsilon _{\overline{1}}^{100}+2v\mathrm{cosh}\left(\kappa ^\pm a\right)\pm \gamma .$$
(79)
After replacing (72), (73), (76), and (78) in (III) we get a system of nine equations with the following nine unknown quantities

$`a_1^{}^0,a_0^+,a_0^{},A_R,A_T^0,A_E^{},A_E^+,A_T^{}\text{and }A_T^+.`$

The solution of this system permit to calculate the electronic current taking into account the e-ph and ph-ph interactions and to obtain a set of five kinetic equations as it is shown in the next section.

## IV The kinetic equations

To get the kinetic equations that describe the population of electrons and phonons we need to obtain the equations of motion for the electron number operators $`c_j^{}c_j`$ and the phonon number operators $`b_i^{}b_i`$

$$i\mathrm{}\frac{d}{dt}\left(c_j^{}c_j\right)=[c_j^{}c_j,],i\mathrm{}\frac{d}{dt}\left(b_i^{}b_i\right)=[b_i^{}b_i,].$$
(80)
We define the average populations for electrons in the excited and ground states as

$`\overline{n_j}={\displaystyle \underset{occ}{}}\mathrm{\Phi }|c_j^{}c_j|\mathrm{\Phi }\text{ for }j=0,1^{}.`$

From the definition (21) of $`|\mathrm{\Phi }`$, the electrons average populations can be calculated as

$$\overline{n_j}=\underset{𝐧,occ}{}\left|p_j^𝐧\right|^2$$
(81)
For the $`LO`$, $`\stackrel{~}{LO}`$, and $`TA`$ phonons, the average populations are

$$\overline{n_i}=\underset{occ}{}\mathrm{\Phi }|b_i^{}b_i|\mathrm{\Phi }=\underset{𝐧,occ,j}{}\left|p_j^𝐧\right|^2n_i\text{ for }i=1,2,3$$
(82)
Using the equations (80) we get the set of quantum kinetic equations

$`{\displaystyle \frac{d\overline{n_1^{}}}{dt}}`$ $`=`$ $`G_1^{}^IG_1^{}^OG_1^{}^E`$ (84)
$`{\displaystyle \frac{d\overline{n_0}}{dt}}`$ $`=`$ $`G_0^IG_0^O+G_1^{}^E`$ (85)
$`{\displaystyle \frac{d\overline{n_1}}{dt}}`$ $`=`$ $`G_1^{}^EG_{LO}^E`$ (86)
$`{\displaystyle \frac{d\overline{n_2}}{dt}}`$ $`=`$ $`G_{LO}^ED_{\stackrel{~}{LO}}`$ (87)
$`{\displaystyle \frac{d\overline{n_3}}{dt}}`$ $`=`$ $`G_{LO}^EE_{TA},`$ (88)
where

$$G_1^{}^I=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}v_{\overline{1}1^{}}Im\left\{p_1^{}^𝐧^{}p_{\overline{1}}^𝐧\right\},$$
(89)
is the electronic input rate from the emitter to the excited level,

$$G_1^{}^O=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}v_{1^{}1}Im\left\{p_1^𝐧^{}p_1^{}^𝐧\right\},$$
(90)
is the output rate from there to the collector.

$$G_1^{}^E=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}gIm\left\{p_0^{𝐧^+}p_1^{}^𝐧\right\}\sqrt{n_1+1},$$
(91)
represents the net balance between emission and absorption of $`LO`$ phonons via electron transitions,

$$G_0^I=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}v_{\overline{1}0}Im\left\{p_0^𝐧^{}p_{\overline{1}}^𝐧\right\},$$
(92)
is the input rate from the emitter directly to the ground state and

$$G_0^O=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}v_{01}Im\left\{p_1^𝐧^{}p_0^𝐧\right\}.$$
(93)
stands for the output rate from the ground state to the collector.

The balance between the decay of one phonon $`LO`$ by emission of a pair of $`\stackrel{~}{LO}`$ and $`TA`$ phonons, and its inverse process (i.e., recombination) is given by

$$G_{LO}^E=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}\gamma Im\left\{p_j^{𝐧^+}p_j^𝐧^{}\right\}\sqrt{N},$$
(94)
where $`N\left(n_1+1\right)\left(n_2+1\right)\left(n_3+1\right)`$.

The decay rate of the $`\stackrel{~}{LO}`$ population turns out to be

$$D_{\stackrel{~}{LO}}=2\kappa _2\underset{𝐧,j}{}\left\{p_j^𝐧^{}p_j^𝐧\right\}n_2=2\kappa _2\overline{n_2}.$$
(95)
Finally the escape rate of the $`TA`$ phonons is given by

$$E_{TA}=2\kappa _3\underset{𝐧,j}{}\left\{p_j^𝐧^{}p_j^𝐧\right\}n_3=2\kappa _3\overline{n_3}.$$
(96)
A diagram for the input and output electron rates is shown in figure 5.

The kinetic equations (IV) can be approximated by a system involving only populations. To do that we have to replace the amplitudes $`p_j^𝐧`$ from (23) and (III) (See the appendix for a detailed calculation). Then, we get a new set of kinetic equations

$`{\displaystyle \frac{d\overline{n_1^{}}}{dt}}=GR_1^{}\overline{n_1^{}}R_w`$ (98)
$`{\displaystyle \frac{d\overline{n_0}}{dt}}=R_0\overline{n_0}+R_w`$ (99)
$`{\displaystyle \frac{d\overline{n_1}}{dt}}=R_wR_\gamma `$ (100)
$`{\displaystyle \frac{d\overline{n_2}}{dt}}=\kappa _2\overline{n_2}+R_\gamma `$ (101)
$`{\displaystyle \frac{d\overline{n_3}}{dt}}=\kappa _3\overline{n_3}+R_\gamma `$ (102)
where

$`R_w`$ $`=`$ $`w\left[\overline{n_1^{}}\left(\overline{n_1}+1\right)\overline{n_0}\overline{n_1}\right]`$ (103)
$`R_\gamma `$ $`=`$ $`\gamma _o[\overline{n_1}\left(\overline{n_2}+1\right)\left(\overline{n_3}+1\right)\left(\overline{n_1}+1\right)\overline{n_2}\overline{n_3}].`$ (104)
$`G`$ is the input current, $`R_1^{}`$ and $`R_0`$ are the escape rates for one electron through the right barrier and $`w`$ is the transition rate due to $`LO`$ phonon emission estimated in the appendix (and more roughly in ). The other parameters are: the decay rate of a $`LO`$ phonon $`\gamma _o`$, the decay rate of a $`\stackrel{~}{LO}`$ phonon $`\kappa _2`$ and the $`TA`$ phonon escape rate $`\kappa _3`$.

It can be seen in the expressions (103) and (104) that the Bose factors were obtained. That means that the phonons are produced in a stimulated regime. For the case of $`TA`$ phonons, the Bose factor $`\overline{n_3}+1`$ is thousands times greater than the others.

## V Results and conclusions

The results presented here correspond only to the stationary case, i.e., $`p_j^𝐧\left(t\right)=a_j^𝐧\left(\omega \right)e^{i\omega t}`$. In this case, all the populations (81) and (82), and the rates (89)-(96) are time independent and they can be written as the same expressions of $`a_j^𝐧\left(\omega \right)`$ instead of $`p_j^𝐧\left(t\right)`$. The amplitudes $`a_j^𝐧\left(\omega \right)`$ are the solutions of the system (III).

The following parameters were used in our calculations. The barriers are $`42.4`$Å wide and $`300meV`$ high ($`38\%`$ of Al concentration). The well has a width of $`203.4`$Å. The Fermi energy at the emitter was fixed at $`15meV`$. The values estimated for the decay rate $`\kappa _2`$ and the escape rate $`\kappa _3`$ were $`30ps^1`$ and $`0.05ps^1`$ respectively. The other parameters used are: $`w`$, that depends on the applied bias but it is of the order of $`2.5ps^1`$, $`\gamma _o=0.11ps^1`$ taken from , and $`g`$, that also depends on the applied bias and it is of the order of $`3.2ps^1`$. The value of $`\gamma `$ was taken equal to $`0.65meV`$ and was calculated from the experimental value of $`\gamma _o`$ using  and .

The transmittance through the DBH is calculated from

$$T=\frac{\left|A_T^0\right|^2+\left|A_T^+\right|^2+\left|A_T^{}\right|^2}{\left|A_I\right|^2}.$$
(105)
A plot of the transmittance as a function of energy is shown in figure 6. It can be seen that when the ph-ph interaction is considered (taking $`\gamma 0`$), a central peak appears. The other two peaks correspond to the polaronic branches already shown in previous works , .

Once the transmission probabilities through each channel are obtained, the output electronic current (equal to the input current at the stationary regime) through the structure can be calculated by summing over the occupied states at the emitter. This sum can be easily transformed into the integral

$$G\left(V\right)=\frac{S}{2\pi ^2}\frac{m^{}}{\mathrm{}^3}\underset{0}{\overset{\epsilon _F^L}{}}\frac{\left(\epsilon _F^L\epsilon \right)j\left(\epsilon \right)\rho \left(\epsilon \right)}{1+R+T}𝑑\epsilon $$
(106)
where $`S=0.510^3mm^2`$ is the device area, $`m^{}=0.067m_e`$ is the effective electron mass, $`R`$ is the reflectance and

$`j\left(\epsilon \right)`$ $`=`$ $`{\displaystyle \frac{4v}{\mathrm{}}}{\displaystyle \frac{1}{\left|A_I\right|^2}}\left\{\left|A_T^0\right|^2\mathrm{sin}\left(k_0^{^{}}a\right)+\left|A_T^+\right|^2\mathrm{sin}\left(k^+a\right)+\left|A_T^{}\right|^2\mathrm{sin}\left(k^{}a\right)\right\},`$
$`\rho \left(E\right)`$ $`=`$ $`{\displaystyle \frac{1}{2v\mathrm{sin}\left(k_0a\right)}}.`$
Instead of using the expression (106) to calculate the electronic current $`G`$, it is calculated from the sum of the terms (89) and (92).

The difference between the currents calculated for $`\gamma 0`$ and $`\gamma =0`$ is appreciable only in the dashed region of figure 7, that it is amplified in the inset.

Finally, the saser intensity and the $`TA`$ phonon population are shown in figure 8. This curve is obtained by solving the kinetic equations (IV).

In conclusion it was obtained here a set of kinetic equations using a full quantum treatment. This set is different from those presented in a previous work  where the kinetic equations were derived phenomenologically. Besides, we consider for the first time all the relevant phonons and its interactions in the calculation of the current.

It was shown that when the phonon-phonon interaction is considered, only small modifications are observed for the current.

The accurate full quantum calculations for a DBH phonon laser confirm the results obtained using rougher phenomenological methods .

## acknowledgments

One of us (IC) wants to acknowledge Latin American agency CLAF for financial support.

##

In this appendix, we will show the way to obtain each term of equations (IV) from (IV). To illustrate the procedure, we will limit us to show only two examples.

We begin with the term

$$G_1^{}^O=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}v_{1^{}1}Im\left\{p_1^𝐧^{}p_1^{}^𝐧\right\}\text{.}$$
(107)
The channel $`000`$ is the only one that contribute significantly to the sum appearing in (107) (see Fig. 3 and Fig. 4). Taking this into account, the expression (107) can be rewriten as

$$G_1^{}^O=\frac{2}{\mathrm{}}\underset{occ}{}v_{1^{}1}Im\left\{p_1^0^{}p_1^{}^0\right\}\text{.}$$
(108)
As the amplitudes appearing in (108) are time dependent, we substitute one of them by their Fourier transform using (30)

$$Im\left\{p_1^0^{}p_1^{}^0\right\}=Im\left\{\frac{1}{\sqrt{2\pi }}_{\mathrm{}}^{\mathrm{}}𝑑\omega a_1^0^{}\left(\omega \right)e^{i\omega t}p_1^{}^0\right\}.$$
(109)
To obtain an expression for the amplitude $`a_1^0^{}\left(\omega \right)`$, we use (56) and neglect the term containing the renormalized hopping $`\stackrel{~}{v}_{\overline{1}1}^0`$ in it. Taken into account the dispersion relation (74), (65) and using (72) to calculate the amplitude in the plane $`j=2`$, we get

$$a_1^0^{}\left(\omega \right)=\mathrm{\Xi }_{1^{}1}a_1^{}^0^{}\left(\omega \right)$$
(110)
where

$`\mathrm{\Xi }_{1^{}1}=v_{1^{}1}{\displaystyle \frac{\left(\mathrm{\Sigma }_1^0ve^{ik_0^{^{}}a}\right)}{\left[\mathrm{\Sigma }_1^{0^2}2v\mathrm{\Sigma }_1^0\mathrm{cos}\left(k_0^{^{}}a\right)+v^2\right]}}.`$

Replacing (110) in (109)

$$Im\left\{p_1^0^{}p_1^{}^0\right\}=\frac{1}{\sqrt{2\pi }}Im\left\{_{\mathrm{}}^{\mathrm{}}\mathrm{\Xi }_{1^{}1}a_1^{}^0^{}\left(\omega \right)𝑑\omega e^{i\omega t}p_1^{}^0\right\}.$$
(111)
The term $`\mathrm{\Xi }_{1^{}1}`$ inside the integral varies very slowly in the range where the amplitude $`a_1^{}^0^{}\left(\omega \right)`$ is not null. Therefore, we can take it out of the integral. Now, performing the integral over $`\omega `$ and taking the imaginary part we get

$$G_1^{}^O=\frac{2}{\mathrm{}}\underset{occ}{}v_{1^{}1}^2\frac{v}{\left[\mathrm{\Sigma }_1^{0^2}2v\mathrm{\Sigma }_1^0\mathrm{cos}\left(k_0^{^{}}a\right)+v^2\right]}\mathrm{sin}\left(k_0^{^{}}a\right)\left|p_1^{}^0\right|^2.$$
(112)
Defining $`\overline{n_1^{}}_{occ}\left|p_1^{}^0\right|^2`$ we get

$`G_1^{}^O=R_1^{}\overline{n_1^{}}`$

where

$`R_1^{}={\displaystyle \frac{2}{\mathrm{}}}v_{1^{}1}^2{\displaystyle \frac{v}{\left[\mathrm{\Sigma }_1^{0^2}2v\mathrm{\Sigma }_1^0\mathrm{cos}\left(k_0^{^{}}a\right)+v^2\right]}}\mathrm{sin}\left(k_0^{^{}}a\right).`$

When the renormalization efects are not taken into account, i.e., $`\mathrm{\Sigma }_1^0=0`$, the expression obtained here for $`R_1^{}`$ is the same as that appearing in .

The other term to be analyzed here is the rate of $`LO`$ phonon emission (91)

$$G_1^{}^E=\frac{2}{\mathrm{}}\underset{𝐧,occ}{}gIm\left\{p_0^{𝐧^+}p_1^{}^𝐧\right\}\sqrt{n_1+1}.$$
(113)
First, we expand the term inside the sum and change the range over which the $`n_1`$ runs to obtain (103). That is

$$\underset{𝐧}{}Im\left\{p_0^{𝐧^+}p_1^{}^𝐧\right\}\sqrt{n_1+1}=\underset{𝐧}{}i\left(p_0^{\left(n_1+1\right)n_2n_3^{}}p_1^{}^𝐧\sqrt{n_1+1}p_1^{}^{\left(n_11\right)n_2n_3^{}}p_0^𝐧\sqrt{n_1}\right).$$
(114)
Writing $`p_0^{\left(n_1+1\right)n_2n_3^{}}`$ in the first term of (114) in terms of its Fourier transform $`a_0^{\left(n_1+1\right)n_2n_3^{}}`$, doing the same with $`p_1^{}^{\left(n_11\right)n_2n_3^{}}`$ in the second term and using equations (34), (35) (only the terms that contain the factor $`g`$) we get

$$a_0^{\left(n_1+1\right)n_2n_3}=\frac{g}{\mathrm{\Delta }_0^{\left(n_1+1\right)n_2n_3}}\sqrt{n_1+1}a_1^{}^𝐧$$
(115)
and

$$a_1^{}^{\left(n_11\right)n_2n_3}=\frac{g}{\mathrm{\Delta }_1^{}^{\left(n_11\right)n_2n_3}}\sqrt{n_1}a_0^𝐧,$$
(116)
where

$`\mathrm{\Delta }_0^{\left(n_1+1\right)n_2n_3}=\mathrm{}\omega \left(n_1+1\right)\mathrm{}\omega _1n_2\mathrm{}\omega _2n_3\mathrm{}\omega _3\epsilon _0`$

and

$`\mathrm{\Delta }_1^{}^{\left(n_11\right)n_2n_3}=\mathrm{}\omega \left(n_11\right)\mathrm{}\omega _1n_2\mathrm{}\omega _2n_3\mathrm{}\omega _3\epsilon _1^{}.`$

Replacing (115) and (116) in (113) and doing the inverse Fourier transform for the amplitudes with the same approximation as those made to get equation (112)

$$G_1^{}^E=\underset{𝐧,occ}{}\left\{\frac{g^2}{\mathrm{}\mathrm{\Delta }_0^{\left(n_1+1\right)n_2n_3}}\left(n_1+1\right)\left|p_1^{}^𝐧\right|^2\frac{g^2}{\mathrm{}\mathrm{\Delta }_1^{}^{\left(n_11\right)n_2n_3}}n_1\left|p_0^𝐧\right|^2\right\}.$$
(117)
If we do now a mean field approximation to take the phonon terms out of the sums, we get the following expression

$$G_1^{}^E=\left(\overline{n_1}+1\right)\underset{𝐧,occ}{}\frac{g^2}{\mathrm{}\mathrm{\Delta }_0^{\left(n_1+1\right)n_2n_3}}\left|p_1^{}^𝐧\right|^2\overline{n_1}\underset{𝐧,occ}{}\frac{g^2}{\mathrm{}\mathrm{\Delta }_1^{}^{\left(n_11\right)n_2n_3}}\left|p_0^𝐧\right|^2.$$
(118)
The amplitudes $`p_1^{}^{000}`$ and $`p_0^{100}`$ are the only ones that contribute significantly to the first and second sum, respectively. Taking this into account and the fact that close to the resonance ($`\epsilon _1^{}\epsilon _0`$ $`\mathrm{}\omega _1`$) $`\mathrm{\Delta }_1^{}^0\mathrm{\Delta }_0^{100}=\mathrm{\Delta }`$, equation (118) can be written as

$$G_1^{}^E=\left(\overline{n_1}+1\right)\underset{occ}{}\frac{g^2}{\mathrm{}\mathrm{\Delta }}\left|p_1^{}^0\right|^2\overline{n_1}\underset{occ}{}\frac{g^2}{\mathrm{}\mathrm{\Delta }}\left|p_0^{100}\right|^2$$
(119)
Finally the expression (119) can be written as

$`G_1^{}^Ew\left[\overline{n_1^{}}\left(\overline{n_1}+1\right)\overline{n_0}\overline{n_1}\right]R_w`$ (120)
$`\text{where }w={\displaystyle \underset{occ}{}}{\displaystyle \frac{g^2}{\mathrm{}\mathrm{\Delta }}},`$ (121)
$`\overline{n_0}{\displaystyle \underset{occ}{}}\left|p_0^{100}\right|^2\text{ and }\overline{n_1^{}}{\displaystyle \underset{occ}{}}\left|p_1^{}^0\right|^2.`$ (122)
The other terms appearing in (IV) can be transformed by the same procedure described above.

In summary, to obtain the kinetic equations similar to that used in , i.e., involving only populations, we did not take into account the crossed terms. For example, if we consider the term with $`v_{\overline{1}0}`$ in (34) to obtain (115), then the expression for the amplitude $`a_0^{\left(n_1+1\right)n_2n_3}`$ would have a term depending on $`v_{\overline{1}0}`$. That term would describe a virtual process. In this process, an electron tunnel through the left barrier ($`j=\overline{1}`$) directly to the ground state ($`j=0`$) emitting one $`LO`$ phonon. Even when this process could occur, its probability is much smaller than that of the indirect process: first tunnel to the excited level, and then decay to the ground state emitting a $`LO`$ phonon. Therefore, we make the sum of each process as being independent. The results obtained from equations (IV) are quite similar to those obtained from (IV).


format: text

text:
# Neutrinos, Large Extra Dimensions and Solar Neutrino Puzzle*footnote **footnote *Invited talk presented by R. N. Mohapatra at the workshop on Lepton Flavor Violation, Honolulu, October (2000)

## I Introduction

One of the important predictions of string theories is the existence more than three space dimensions. For a long time, it was believed that these extra dimensions are small and are therefore practically inconsequential as far as low energy physics is concerned. However, recent progress in the understanding of the nonperturbative aspects of string theories have opened up the possibility that some of these extra dimensions could be large without contradicting observations. In particular, models, where some of the extra dimensions have sizes as large as a milli-meter and where the string scale is in the few TeV range have attracted a great deal of phenomenological attention in the past two years. The basic assumption of these models, largely inspired by generic observations in string theories is that the space time has a brane-bulk structure, where the brane is the familiar (3+1) dimensional space-time, with the standard model particles and forces residing in it and bulk consists of all space dimensions where gravity and other possible gauge singlet particles live. The main interest in these models is of course due to the fact that the low string scale provides an opportunity to test them using existing collider facilities.

A major challenge to these theories comes from the neutrino sector, the first problem being how one understands the small neutrino masses in a natural manner. The conventional seesaw explanation which is believed to provide the most satisfactory way to understand this, requires that the new physics scale (or the scale of $`SU(2)_R\times U(1)_{BL}`$ symmetry) be around $`10^{12}`$ GeV or higher and clearly does not work for these models. So one must look for alternative ways. The second problem is that if one considers only the standard model group in the brane, operators such as $`LHLH/M^{}`$ could be induced by string theory in the low energy effective Lagrangian. For TeV scale strings this would obviously lead to unacceptable neutrino masses.

One mechanism suggested in Ref. is to postulate the existence of one or more gauge singlet neutrinos, $`\nu _B`$ in the bulk which couple to the lepton doublets in the brane. After electroweak symmetry breaking this coupling can lead to neutrino Dirac masses, which are suppressed by the ratio $`M_{}/M_P\mathrm{}`$. This is sufficient to explain small neutrino masses and owes its origin to the large bulk volume that suppresses the effective Yukawa couplings of the Kaluza-Klein (KK) models of the bulk neutrino to the brane fields. In this class of models, naturalness of small neutrino mass requires that one must assume the existence of a global B-L symmetry in the theory, since that will exclude the undesirable higher dimensional operators from the theory.

An alternative possibility is to consider the brane theory to have an extended gauge symmetry which contains B-L symmetry as a subgroup, which will eliminate the higher dimensional operators. This is perhaps more in the spirit of string theories, which generally do not allow global symmetries. Phenomenological considerations, however, require that the local $`BL`$ scale and hence the string scale be of order of $`10^{10}`$ GeV or so. The extra dimensions in these models could of course be large.

Regardless of which path one chooses for understanding small neutrino masses, a very desirable feature of these models is that if the size of extra dimensions is of order of a milli-meter, the KK excitations of the bulk neutrino have masses in the range of $`10^3`$ eV, which is in the range needed for a unified understanding of oscillation data.

In this paper, we will focus on TeV scale models and discuss their implications. In this class of models, one gets for the Dirac masses of the familiar neutrinos $`\nu _{e,\mu ,\tau }`$

$`m_0={\displaystyle \frac{hv_{wk}M_{}}{M_P\mathrm{}}}.`$ (1)
The righthanded neutrino in the Dirac mass is a mode of the bulk neutrino $`\nu _B`$. For $`M_{}10`$ TeV, this leads to $`m_\nu 10^4h`$ eV. It is encouraging that this number is in the right range to be of interest in the discussion of solar neutrino oscillation if the Yukawa coupling $`h`$ is appropriately chosen.

If one takes these models for the neutrinos seriously, some immediate questions arise: first, while the KK modes of the bulk neutrinos are candidates for sterile neutrinos if the bulk radius is in the range of millimeters, are their small masses natural? Secondly, how well do the models of this type describe the observed oscillation data from solar, atmospheric and LSND experiments ? It is clear that they provide a tempting possibility to simultaneously explain all data since the needed sterile neutrino could be one of the KK modes of the bulk neutrino, as already emphasized.

We discuss the answer to the first question in section 2. As far as the second question is concerned, several recent papers have addressed this issue. In particular, in Ref. it has been shown that while the overall features of the solar and atmospheric data can be accomodated in minimal versions of these models with three bulk neutrinos, it is not possible to simultaneously explain the LSND observation for the $`\nu _\mu \nu _e`$ oscillation probability, and one must incorporate new physics in the brane.

On a phenomenological level, a first glance at the values of the parameters of the model such as $`m_0`$ from Eq.(1) and $`R^110^3`$ eV suggests that perhaps one should seek a solution of the solar neutrino data in these models using the small angle MSW mechanism. However, present Super-Kamiokande recoil energy distribution seems to disfavor such an interpretation, although any definitive conclusion should perhaps wait till more data accumulates. In any case if the present trend of the data near the higher energy region of the solar neutrino spectrum from Super-Kamiokande persists, it is likely to disfavor the small angle MSW solution and favor a vacuum oscillation. We therefore pursue the possibility of vacuum oscillation between the $`\nu _e`$ and $`\nu _s`$ to solve the solar neutrino puzzle motivated by a model that leads to desired parameters in a natural manner.

As is known from many discussions in literature, chlorine results play a pivotal role in deciding on the nature of the oscillation solution to the solar neutrino puzzle. To fit water and the Gallium data in conjunction with the Chlorine data, one must “kill” the 7Be neutrinos almost completely. In the vacuum solution, it is in general hard to achieve this without simultaneously suppressing the pp neutrinos. Therefore, the general tactic adopted is to suppress the 7Be neutrinos as much as possible (but not quite to zero) and suppressing the 8B neutrinos to less than 30% or so. This, however, suppresses the higher energy 8B contributions below the level observed by Super-Kamiokande collaboration. In the case of oscillation to active neutrinos, this can be fixed by including the neutral current scattering which contributes about 16% of the charged current data. On the other hand, this disfavors the vacuum oscillation to sterile neutrinos which do not have a neutral current component. As we will show later, if the sterile neutrino is a bulk mode, its KK modes play a significant role in resolving this problem due to the existence of MSW transitions of $`\nu _e`$ to these modes for values of extra dimension sizes of interest.

## II Naturalness of ultra light sterile neutrinos in brane-bulk models

Before getting to the detailed discussion of the model, we make some general comments about the naturalness of the ultralightness of the bulk neutrino. The bulk neutrino “self mass” terms are constrained by the geometry of the bulk and could therefore under certain circumstances be zero. If that happens, the only mass of the KK states of the $`\nu _B`$ will arise from the kinetic energy terms such as $`\overline{\nu }_B\mathrm{\Gamma }^I_I\nu _B`$, where $`I=5,6\mathrm{}`$ and will be given by $`n/R`$ where $`R`$ is the radius of the extra dimensions. In such a situation, an ultralight $`\nu _{B,KK}`$ arises naturally.

The key to naturalness of the ultralight bulk neutrino is the geometry that forbids both Dirac and Majorana mass terms. Let us give a few examples. In five dimensions, if we impose the $`Z_2`$ orbifold symmetry ($`yy`$), then it follows that the Dirac mass vanishes. Now if we impose lepton number symmetry in the brane, the Majorana mass vanishes, leaving us with no mass term for the bulk neutrino in 5-dimensions.

Another interesting example is the 10-dimensional bulk, where the bulk neutrino is a 16-component spinor, which when reduced to 4-dimensions leads to eight 2-component spinors. The interesting point is that for a 16-dimensional spinor, one cannot write a Dirac or Majorana mass term consistent with 10-dimensional Lorentz invariance. In this case, there is no need for assuming lepton number to get an ultralight sterile neutrino. A similar situation is also expected in six dimensions if we choose the bulk neutrino to be a 4-component complex chiral spinor.

## III Basic ideas of the model and the mixing of bulk modes

The model we will present will consist of the standard model in the brane (suitably extended to include some additional Higgs fields) and one bulk neutrino. The bulk may be five, six or higher dimensional; we will assume that only one of those extra dimensions is large. The bulk neutrino will be assumed to couple to $`\nu _e`$In general, the bulk neutrino will couple to all three species of brane neutrinos; however, when the $`\nu _{\mu ,\tau }`$ masses are much bigger than the mass parameters in the $`\nu _e`$ and $`\nu _B`$ sector, one can decouple the $`\nu _{\mu ,\tau }`$ and restrict oneself to only the $`\nu _e,\nu _B`$ sector, to analyse the effect of the bulk neutrino,as we do in the paper.. Let us begin our discussion by focusing on the $`\nu _e`$ in the brane and one bulk neutrino. Obviously, the fields that could propagate in the extra dimensions are chosen to be gauge singlets. Let us denote bulk neutrino by $`\nu _B(x^\mu ,y)`$. It has a five dimensional kinetic energy term and a coupling to the brane field $`L(x^\mu )`$ given by

$$=\kappa \overline{L}H\nu _{BR}(x,y=0)+𝑑y\overline{\nu }_{BL}(x,y)_5\nu _{BR}(x,y)+h.c.,$$
(2)
where from the five dimensional kinetic energy, we have only kept the 5th component that contributes to the mass terms of the KK modes in the brane; $`H`$ denotes the Higgs doublet, and $`\kappa =h\frac{M^{}}{M_P\mathrm{}}`$ the suppressed Yukawa coupling. It is worth pointing out that this suppression is independent of the number and radius hierarchy of the extra dimensions, provided that our bulk neutrino propagates in the whole bulk. For simplicity, we will assume that there is only one extra dimension with radius of compactification as large as a millimiter, and the rest with much smaller compactification radii. The smaller dimensions will only contribute to the relationship between the Planck and the string scale, but their KK excitations will be very heavy and decouple from neutrino spectrum. Thus, all the analysis could be done as in five dimensions.

A second point we wish to note is that we will include new physics in the brane that will generate a Majorana mass matrix for the three standard model neutrinos as follows:

$`=(\begin{array}{ccc}\delta _{ee}& \delta _{e\mu }& \delta _{e\tau }\\ \delta _{e\mu }& \delta _{\mu \mu }& m_0\\ \delta _{e\tau }& m_0& \delta _{\tau \tau }\end{array}.)`$ (6)
The origin of this pattern of brane neutrino masses will be discussed in a latter section. In this section let’s discuss the mixing pattern of the bulk neutrinos with the brane ones. For this we will assume that $`m_0\delta _{ij}`$; as a result the $`\nu _{\mu ,\tau }`$ in effect do not affect the mixing between the bulk neutrino modes and the $`\nu _e`$.

To discuss the effect of the bulk modes, note that the first term in Eq. (2) will be responsible for the neutrino mass once the Higgs field develops its vacuum expectation. The induced Dirac mass parameter will be given by $`m=\kappa v`$, which for $`M^{}=1`$ TeV is about $`h10^5`$ eV. Obviously this value depends only linearly on the fundamental scale. Larger values for $`M^{}`$ will increase $`m`$ proportionally. After introducing the expansion of the bulk field in terms of the KK modes, the Dirac mass terms in (2) could be written as

$$(\nu _e\nu _{0B}\nu _{B,}^{}\nu _{B,+}^{})\left(\begin{array}{cccc}\delta _{ee}& m& \sqrt{2}m& 0\\ m& 0& 0& 0\\ \sqrt{2}m& 0& 0& _5\\ 0& 0& _5& 0\end{array}\right)\left(\begin{array}{c}\nu _e\\ \nu _{0B}\\ \nu _{B,}^{}\\ \nu _{B,+}^{}\end{array}\right),$$
(7)
where our notation is as follows: $`\nu _B^{}`$ represents the KK excitations; the off diagonal term $`\sqrt{2}m`$ is actually an infinite row vector of the form $`\sqrt{2}m(1,1,\mathrm{})`$. The operator $`_5`$ stands for the diagonal and infinite KK mass matrix whose $`n`$-th entry is given by $`n/R`$. This notation was introduced in  to represent the infinite mass matrix in a compact manner.

Using this short hand notation makes it easier to calculate the exact eigenvalues and the eigenstates of this mass matrix . Simple algebra yields the characteristic equation

$$m_n=\delta _{ee}+\frac{\pi m^2}{\mu _0}\mathrm{cot}(\frac{\pi m_n}{\mu _0}),$$
(8)
where $`m_n`$ is the mass eigenvalue . The equation for eigenstates is

$`\stackrel{~}{\nu }_n={\displaystyle \frac{1}{N_n}}\left[\nu _e+{\displaystyle \frac{m}{m_n}}\nu _{0B}+{\displaystyle \frac{\sqrt{2}m\left(m_n\nu _{B,}^{}+_5\nu _{B,+}^{}\right)}{m_n^2_5^2}}\right],`$ (9)
where the sum over the KK modes in the last term is implicit, and $`N_n`$ is the normalization factor given by

$$N_n^2=1+m^2\pi ^2R^2+\frac{(m_n\delta _{ee})^2}{m^2}.$$
(10)
Using the expression (9), we can write down the weak eigenstate $`\nu _L`$ in terms of the massive modes as

$$\nu _L=\underset{n=0}{\overset{\mathrm{}}{}}\frac{1}{N_n}\stackrel{~}{\nu }_{nL}.$$
(11)
Thus, the weak eigenstate is actually a coherent superposition of an infinite number of massive modes. Therefore, even for this single flavour case, the time evolution of the mass-eigenstates involves in principle all mass eigenstates and is very different from the simple oscillatory behaviour familiar from the conventional two or three neutrino case. The survival probability depends strongly on the parameter $`\xi `$, reflecting the universal coupling of all the KK components of $`\nu _B`$ with $`\nu _L`$ in (2). We will be interested in the small $`\xi `$ region of the parameter space. In this case it is easy to see that the coupling of the nth mode is given by $`\xi /n`$.

Note that in the limit of $`\delta _{ee}=0`$, the $`\nu _e`$ and $`\nu _{B,0}`$ are two, two-component spinors that form a Dirac fermion with mass m. Once we include the effect of $`\delta _{ee}0`$, they become Majorana fermions with masses given by: $`m_1\delta _{ee}/2+m`$ and $`m_2\delta _{ee}/2m`$, and they are maximally mixed; i.e., the two mass eigenstates are $`\nu _{1,2}\frac{\nu _e\pm \nu _{B,0}}{\sqrt{2}}`$. Thus as the $`\nu _e`$ produced in a weak interaction process evolves, it oscillates to the $`\nu _{B,0}`$ state with an oscillation length $`LE/(2m\delta _{ee})`$. Later on in this paper we will find that for natural values of $`m,\delta _{ee}`$ in this model, $`L`$ is of order of the Sun-Earth distance so that our model leads to vacuum oscillation of the solar neutrinos.

Furthermore since the $`\nu _e`$ also mixes with the KK modes of the bulk neutrinos with a mass difference square of order $`10^6`$ eV2, there is MSW resonance transition of $`\nu _e`$ to $`\nu _{B,KK}`$ modes. In our discussion of the solar neutrino puzzle, we will take into account both these effects.

In our model, we will also have $`m_0`$ in the brane neutrino mass matrix which is of order eV so that we can account for LSND results. Atmospheric neutrino oscillation is now purely $`\nu _\mu `$-$`\nu _\tau `$ oscillation.

## IV Solar neutrino data by a combination of vacuum and MSW oscillations

Let us now discuss how the solar neutrino data is explained within this model. To discuss the MSW effect for the case of bulk neutrinos, we need to include the matter effect in the diagonalization of the infinite dimensional neutrino mass matrix. The eigenvectors and eigenvalues of this matrix can be found in the presence of matter, when the squared mass matrix, $`M^2`$, is replaced with $`M^2+2EH_1`$, where $`H_1`$ is $`\rho _e=\sqrt{2}G_F(n_e0.5n_n)`$ when acting on $`\nu _e`$, and is zero on sterile neutrinos. Define

$$w_k=\frac{E\rho _e}{m_k\delta _{ee}}+\sqrt{1+\left(\frac{E\rho _e}{m_k\delta _{ee}}\right)^2};$$
(12)
$`w_k=1`$ in vacuum. The characteristic equation becomes

$$m_k=w_k\delta _{ee}+\frac{\pi m^2}{\mu _0}cot\frac{\pi m_k}{\mu _0},$$
(13)
The eigenvectors are as in Eq. 9, except the coefficients of $`\nu _{0B}`$ and $`\nu _{B,}^{}`$ acquire an additional factor $`1/w_n`$, and the normalization becomes $`N_n^2=1+\frac{(1+\frac{1}{w_n^2})}{2}\left(\pi ^2m^2R^2+\frac{(m_nw_n\delta _{ee})^2}{m^2}\right)\frac{(1\frac{1}{w_n^2})}{2}\frac{(m_nw_n\delta _{ee})}{m}`$.

In fitting the solar data the simplest way to reconcile the rates for the three classes of experiments is to “kill” the 7Be neutrinos, reduce the 8B neutrinos by half and leave the pp neutrinos alone. To achieve this in the VO case, one may put a node of the survival probability function $`P_{ee}`$ around $`0.86`$ MeV. However, for an arbitrary node number, the oscillatory behaviour of $`P_{ee}`$ before and after $`0.86`$ MeV cannot in general satisfy the other two requirements mentioned above. If one uses the first node to “kill 7Be”, then for 8B neutrino energies the $`P_{ee}`$ is close to one and not half as would be desirable; on the other hand, if one uses one of the higher nodes (higher $`\mathrm{\Delta }m^2`$), then pp neutrinos get reduced. The strategy generally employed is not to “shoot” for a node at the precise $`{}_{}{}^{7}Be`$ energy but rather somewhat away so that it reduces 7Be to a value above zero. This requires that one must reduce the 8B neutrinos by much more than 50%, so one can fit Chlorine data. The water data then requires an additional contribution, which, in the case of active VO, is provided by the neutral current cross section amounting to about 16% of the charged current one. Thus in a pure two-neutrino oscillation picture, VO works for oscillation to active neutrinos but does not work for active to sterile oscillation. It is here that the large extra dimensions come to the rescue as we see below, thereby keeping the active to sterile VO fit to solar neutrino data viable.

To see this note that in our model both vacuum oscillations and MSW oscillations are important. This is because the lowest mass pair of neutrinos is split by a very small mass difference, whereas the KK states have to be separated by $`<10^2`$ eV because of the limits from gravity experiments. We can then use the first node of $`P_{ee}`$ to suppress the 7Be. Going up in energy toward 8B neutrinos, the survival probability, which in the VO case would have risen to very near one, is suppressed by the small-angle MSW transitions to the different KK excitations of the bulk neutrino. This is the essence of our new way to fit the solar neutrino observations.

To carry out the fit, we studied the time evolution of the $`\nu _e`$ state using a program that evolved from one supplied by W. Haxton. The program was updated to use the solar model of BP98 and modified to do all neutrino transport within the sun numerically. For example, no adiabatic approximation was used. Changes were also necessary for oscillations into sterile neutrinos and to generalize beyond the two-neutrino model. Up to 16 neutrinos were allowed, but no more than 14 contribute for the solutions we considered.

For comparison with experimental results, tables of detector sensitivity for the Chlorine and Gallium experiments were taken from Bahcall’s web site. The Super-Kamiokande detector sensitivity was modeled using , where the percent resolution in the signal from Cerenkov light, averaged over the detector for various total electron energies, is provided. For details see .

Calculations of electron neutrino survival probability, averaged over the response of detectors, were compared with measurements. While theoretical uncertainties in the solar model and detector response were included in the computation of $`\chi ^2`$ as described in Ref. , the measurement results given here include only experimental statistical and systematic errors added in quadrature. The Chlorine survival probability, from Homestake, is $`0.332\pm 0.030`$. Gallium results for SAGE, GALLEX and GNO were combined to give a survival probability of $`0.579\pm 0.039`$. The $`5.520`$ MeV 1117 day Super-K experimental survival probability is $`0.465\pm 0.015`$. The best fits were with $`R58\mu m`$, $`mR`$ around $`0.0093`$, and $`\delta _{ee}0.84\times 10^7`$ eV, corresponding to $`\delta m^20.53\times 10^{11}`$ eV2. These parameters give average survival probabilities for Chlorine, Gallium, and water of 0.386, 0.533, and 0.460, respectively. They give a $`\nu _e`$ survival probability whose energy dependence is shown in Fig. 1. For two-neutrino oscillations, the coupling between $`\nu _e`$ and the higher mass neutrino eigenstate is given by $`sin^22\theta `$, whereas here the coupling between $`\nu _e`$ and the first KK excitation replaces $`sin^22\theta `$ by $`4m^2R^2=0.00035`$.

Vacuum oscillations between the lowest two mass eigenstates nearly eliminate electron neutrinos with energies of .63MeV/(2n+1) for n = 0, 1, 2, … . Thus Fig. 1 shows nearly zero $`\nu _e`$ survival near 0.63 MeV, partly eliminating the $`{}_{}{}^{7}Be`$ contribution at $`0.862`$ MeV, and giving a dip at the lowest neutrino energy. Increasing $`\delta _{ee}`$ moves the low energy dip to the right into Gallium’s most sensitive pp energy range, making the fit worse. Decreasing $`\delta _{ee}`$ increases Gallium, but hurts the Chlorine fit by moving the higher energy vacuum oscillation dip further to the left of the $`{}_{}{}^{7}Be`$ peak. MSW resonances start causing the third and fourth eigenstates to be significantly occupied above $`0.8`$ MeV, the fifth and sixth eigenstates above $`3.7`$ MeV, the 7’th and 8’th above $`8.6`$ MeV, and the 9’th and 10’th above $`15.2`$ MeV. Fig. 1 shows dips in survival probability just above these energy thresholds. The typical values of the survival probability within the $`{}_{}{}^{8}B`$ region ($`6`$ to $`14`$ MeV) are quite sensitive to the value of $`mR`$. As can be seen from Eq. 10, higher $`mR`$ increases $`1/Nm/m_nmR/n`$ for various $`n`$, and thereby increases $`\nu _e`$ coupling to higher mass eigenstates, strengthens MSW resonances, and lowers $`\nu _e`$ survival probability.

The expected energy dependence of the $`\nu _e`$ survival probability is compared with preliminary Super-K data in Fig. 2. The uncertainties are statistical only. The parameters used in making Fig. 2 were chosen to provide a good fit to the total rates only; they were not adjusted to fit this spectrum. But combining spectrum data with rates using the method described in Ref. , with numbers supplied by Super-K, gives $`\chi ^2=15.6`$ (“probability” 74%) for the spectrum predicted from the fit to total rates. It’s incorrect to calculate probability as if there were (number of data points) \- (number of parameters) = 17 degrees of freedom, but if we do, $`\chi ^2`$ corresponds to a “probability” of 55%.

One may also seek fits with $`\delta _{ee}`$ constrained to be very small, thereby eliminating vacuum oscillations. The best such fit had $`\chi ^2=5.5`$ (“probability” $`14\%`$). The same parameters then used with the Super-K spectrum gave $`\chi ^2=18.7`$ (“probability” $`35\%`$).

The seasonal effect was computed for a few points on the earth’s orbit. If $`r`$ is the distance between the earth and the sun,

$$\frac{r_0}{r}=1+ϵcos(\theta \theta _0),$$
(14)
where $`r_0`$ is one astronomical unit, $`ϵ=0.0167`$ is the orbital eccentricity, and $`\theta \theta _02\pi (tt_0)`$, with $`t`$ in years and $`t_0=`$ January 2, 4h 52m. Table I shows very small seasonal variation.

To understand the mass pattern used, consider $`L_e+L_\mu L_\tau `$ symmetry for neutrinos, with $`L_e=1`$ for $`\nu _B`$. The allowed mass pattern for $`(\nu _e,\nu _\mu ,\nu _\tau )`$ is given by Eq. 8 with all $`\delta `$’s (except $`\delta _{e\tau ,\mu \tau }`$) set to zero. The remaining $`\delta `$ entries arise after we turn on the symmetry breaking. For an explicit realization, we augment the standard model by the singlet charged Higgs $`\eta ^+`$, $`h^{++}`$ which are blind with respect to lepton number and $`SU(2)_L`$ triplet fields, $`\mathrm{\Delta }_{e,\mu ,\tau }`$, with $`Y=2`$ which carry two negative units of lepton numbers $`L_{e,\mu ,\tau }`$, respectively. The Lagrangian involving these fields consists of two parts: one $`_0`$ which is invariant under $`(L_e+L_\mu L_\tau )`$ number and contains terms ($`\eta L_\mu L_\tau ,\eta L_eL_\tau ,h^{++}e_R^{}\tau _R^{},h^{++}\mu _R^{}\tau _R^{})`$ and $`L_eL_e\mathrm{\Delta }_e,L_\mu L_\mu \mathrm{\Delta }_\mu ,L_\tau L_\tau \mathrm{\Delta }_\tau `$ and a soft breaking term $`_1=h^{++}(_{i=e,\mu ,\tau }M_{ii}\mathrm{\Delta }_i^2+M_{0\mu }\mathrm{\Delta }_\mu \mathrm{\Delta }_\tau +M_{0\tau }\mathrm{\Delta }_e\mathrm{\Delta }_\tau )+h.c.`$ With these couplings, the neutrino Majorana masses arise from two-loop effects (similar to the mechanism of ref.), and we have the $`L_e+L_\mu L_\tau `$ violating entries $`\delta _{ij}cm_{e_i}m_{e_j}`$. Using $`\delta _{ee}10^8`$ eV, then $`\delta _{\tau \tau }10^2`$, as would be required to understand the atmospheric neutrino data.

Finally, models with bulk sterile neutrinos lead to new contributions to low energy weak processes, which have been addressed in several papers. In the domain of astrophysics, they lead to new contributions to supernova energy loss, as well as to the energy density in the early universe, which can influence the evolution of the universe. Currently these issues are under discussion, and if neutrino oscillation data favor these models, any cosmological constraints must be addressed. Note, however, that this model is less sensitive to these issues than other fits because for VO $`\mathrm{\Delta }m^2`$ is an order of magnitude smaller, for MSW $`sin^22\theta `$ is more than an order of magnitude smaller, and there is only a single KK tower starting from a very small mass.

The work of R.N.M. is supported by a grant from the National Science Foundation No. PHY-9802551. The work of D.O.C. and S.J.Y. is supported by a grant from the Department Of Energy No. DE-FG03-91ER40618. One of us (R. N. M.) would like to thank the organizers of the “Lepton Flavor Violation Workshop” for kind hospitality and partial support.