format: text

text:
# Untitled Document

hep-ph/0002268 IASSNS–HEP–00–08

Pseudo-Dirac Solar Neutrinos

Yosef Nir

School of Natural Science, Institute for Advanced Study

Princeton, NJ 08540, USA1 Address for academic year 1999-2000

nir@ias.edu

Department of Particle Physics

Weizmann Institute of Science, Rehovot 76100, Israel

ftnir@wicc.weizmann.ac.il

Three of the viable solutions of the solar neutrino problem are consistent with close to maximal leptonic mixing: $`\mathrm{sin}^2\theta _{12}=\frac{1}{2}\left(1ϵ_{12}\right)`$ with $`\left|ϵ_{12}\right|1`$. Flavor models can naturally explain close to maximal mixing if approximate horizontal symmetries force a pseudo-Dirac structure on the neutrino mass matrix. An experimental determination of $`\left|ϵ_{12}\right|`$ and sign($`ϵ_{12}`$) can constrain the structure of the lepton mass matrices and consequently provide stringent tests of such flavor models. If both $`\left|ϵ_{12}\right|`$ and $`\mathrm{\Delta }m_{21}^2`$ are known, it may be possible to estimate the mass scale of the pseudo-Dirac neutrinos. Radiative corrections to close to maximal mixing are negligible. Subtleties related to the kinetic terms in Froggatt-Nielsen models are clarified.

2/00

1. Introduction

Three of the solutions of the solar neutrino problem require a large mixing angle \[1--4\]:

$$\begin{array}{cc}\hfill \mathrm{LMA}:& \mathrm{sin}^22\theta _{12}0.71,\mathrm{\Delta }m_{21}^2(120)\times 10^5eV^2,\hfill \\ \hfill \mathrm{LOW}:& \mathrm{sin}^22\theta _{12}0.81,\mathrm{\Delta }m_{21}^2(330)\times 10^8eV^2,\hfill \\ \hfill \mathrm{VAC}_\mathrm{L}:& \mathrm{sin}^22\theta _{12}0.71,\mathrm{\Delta }m_{21}^2(410)\times 10^{10}eV^2.\hfill \end{array}$$

Here LMA and LOW refer to matter-enhanced oscillations with a large mixing angle in the high and low $`\mathrm{\Delta }m^2`$ ranges, respectively, while VACL refers to vacuum oscillations with relatively large $`\mathrm{\Delta }m^2`$. The range for the mixing angle in (1.1) is close to maximal mixing, $`\mathrm{sin}^22\theta _{12}=1`$. This case is particularly interesting from the theoretical point of view. It follows from a simple structure of the relevant $`2\times 2`$ block in the neutrino mass matrix in the basis where the charged lepton mass matrix is diagonal:

$$M_\nu ^{(2)}=m\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right).$$

Such a structure is easily obtained in models of horizontal symmetries \[5--9\] that try to explain the observed smallness and hierarchy in the charged fermion parameters (mass ratios and mixing angles). For example, if the lepton doublets of the first two generations carry an opposite charge under a U(1) symmetry (and the relevant scalar field is neutral), then $`M_\nu ^{(2)}`$ has the structure (1.1) in the symmetry limit.

Any horizontal symmetry must be broken in Nature. An unbroken horizontal symmetry leads to either degeneracy between fermions of different generations or vanishing mixing angles (see e.g.  and references therein). In particular, the mass degeneracy implied by (1.1) must be broken to satisfy (1.1). The horizontal symmetry still has observable consequences if the breaking parameters are small. Then the low energy effective theory is subject to selection rules that are manifested in the smallness and hierarchy of the flavor parameters. In the case of close-to-maximal mixing, the small breaking leads to a small splitting between the masses of the two neutrinos and to a small deviation from maximal mixing, that is, the two Majorana neutrinos form a pseudo-Dirac neutrino:

$$\frac{\mathrm{\Delta }m_{21}^2}{m^2}1,1\mathrm{sin}^22\theta _{12}1.$$

Here $`m`$ denotes the average of $`m_1`$ and $`m_2`$. A measurement of these small effects will provide further information about the pattern of symmetry breaking and guide us in the process of selecting among the many presently viable models of horizontal symmetries. (For interesting studies of the implications of solar neutrino measurements for small entries in the neutrino mass matrix, see refs. \[11,,12\].)

There are three light active neutrinos in Nature. (In this work we assume that these are the only light neutrinos and do not consider the possibility of light sterile neutrinos. Note that the large angle solutions of the solar neutrino problem are inconsistent with a pseudo-Dirac $`\nu _e\nu _s`$ combination, such as in the model of ref. .) In the two generation framework, where there is a single mixing angle $`\theta `$, maximal mixing is defined by maximal oscillation depth in vacuum and corresponds to $`\mathrm{sin}^22\theta =1`$. In the three generation framework, what we mean by maximal mixing is that the disappearance probability is equivalent to that for maximal two neutrino mixing at the relevant mass scale . The disappearance probability for $`\nu _e`$ in vacuum, $`P_e\mathit{}`$, is given by

$$P_e\mathit{}=4|V_{e1}|^2|V_{e2}|^2\mathrm{sin}^2\mathrm{\Delta }_{\mathrm{sun}}4|V_{e3}|^2\mathrm{sin}^2\mathrm{\Delta }_{\mathrm{atm}}.$$

Here $`V_{ij}`$ are the elements of the Maki-Nakagawa-Sakata (MNS) matrix  and we use the definition $`\mathrm{\Delta }_{jk}\frac{\mathrm{\Delta }m_{jk}^2L}{4E}`$, and the following input from solar and atmospheric neutrino experiments:

$$|\mathrm{\Delta }_{\mathrm{sun}}|=|\mathrm{\Delta }_{21}||\mathrm{\Delta }_{31}||\mathrm{\Delta }_{32}|=|\mathrm{\Delta }_{\mathrm{atm}}|.$$

At the $`L/E`$-scale that is relevant to solar neutrinos, the second term in (1.1) averages out to $`2|V_{e3}|^2`$. The only oscillatory term is the first one, and our definition of maximal mixing corresponds to

$$4|V_{e1}|^2|V_{e2}|^2=1,$$

which leads to

$$|V_{e1}|^2=|V_{e2}|^2=1/2,|V_{e3}|^2=0.$$

In the standard parametrization of the $`V_{\mathrm{MNS}}`$ matrix,

$$V_{\mathrm{MNS}}=\left(\begin{array}{ccc}c_{12}c_{13}& s_{12}c_{13}& s_{13}e^{i\delta }\\ s_{12}c_{23}c_{12}s_{23}s_{13}e^{i\delta }& c_{12}c_{23}s_{12}s_{23}s_{13}e^{i\delta }& s_{23}c_{13}\\ s_{12}s_{23}c_{12}c_{23}s_{13}e^{i\delta }& c_{12}s_{23}s_{12}c_{23}s_{13}e^{i\delta }& c_{23}c_{13}\end{array}\right),$$

where $`s_{ij}\mathrm{sin}\theta _{ij}`$ and $`c_{ij}\mathrm{cos}\theta _{ij}`$, the conditions (1.1) translate into

$$s_{12}^2=1/2,s_{13}=0.$$

By close-to-maximal-mixing we refer to a situation close to (1.1) or, equivalently, to (1.1):

$$ϵ_{12}12s_{12}^21,s_{13}1.$$

Our convention here is that $`\mathrm{\Delta }m_{21}^2m_2^2m_1^2>0`$, so that $`ϵ_{12}>0`$ ($`ϵ_{12}<0`$) corresponds to a situation where the lighter (heavier) state has a larger component of $`\nu _e`$.

Solar neutrino experiments (and, more generally, any oscillation experiments) are sensitive to the mass-squared difference $`\mathrm{\Delta }m_{12}^2`$ but not to the masses themselves. On the other hand, they can be sensitive to small deviations from maximal mixing \[16--18\]. Moreover, matter oscillations (but not vacuum oscillations) are affected differently by $`ϵ_{12}>0`$ and by $`ϵ_{12}<0`$, that is, they are sensitive not only to $`\mathrm{sin}^22\theta _{12}`$ but also to $`\mathrm{sin}^2\theta _{12}`$. In other words, if the solar neutrino problem is solved by one of the large angle solutions, then experiments may provide us with a measurement of the sign and the size of the small parameter $`ϵ_{12}`$. The purpose of this work is to understand the potential lessons for flavor model building from solar neutrino measurements of $`ϵ_{12}`$.

Our interest lies in models where $`ϵ_{12}`$, $`s_{13}`$ and $`\mathrm{\Delta }m_{21}^2/m^2`$ are naturally small. We focus on models where there are no exact relations between entries of the lepton mass matrices (beyond the symmetric structure of the neutrino Majorana mass matrix). The smallness of physical parameters must then be related to the smallness of various entries in the mass matrices and not to fine tuned cancellations between various contributions. As a concerete example of such a framework we think of models of approximate Abelian horizontal symmetries, but most of our results have more general applicability.

Horizontal symmetries constrain the structure of the mass matrices $`M_\nu `$ and $`M_{\mathrm{}}`$. In section 2 we derive the dependence of the mixing angles and of $`\mathrm{\Delta }m_{12}^2/m^2`$ on the entries of the lepton mass matrices. Usually, the constraints of the horizontal symmetries apply at a high energy scale. The effects of renormalization group evolution (RGE) are analyzed in section 3. For specific high energy theories of flavor, such as the Froggatt-Nielsen mechanism , the kinetic terms are corrected in a flavor dependent way when heavy degrees of freedom are integrated out. The effects of non-canonical kinetic terms are studied in section 4. The analysis in sections 2$``$4 is carried out under the simplifying assumption of CP symmetry. Effects of phases are studied within a two generation model in section 5. We apply our results to various models of Abelian flavor symmetries in section 6. We summarize our conclusions in section 7.

2. From Interaction Basis Parameters to Physical Parameters

Flavor models and, in particular, models with horizontal symmetries, constrain the entries of the lepton mass matrices in the interaction basis. To understand the implications of experimental constraints, one needs to express the physical observables (masses and mixing angles) in terms of the interaction basis parameters.

Given the charged lepton mass matrix $`M_{\mathrm{}}`$ and the neutrino mass matrix $`M_\nu `$ in some interaction basis,

$$_M=\left(\begin{array}{ccc}\overline{e_L}& \overline{\mu _L}& \overline{\tau _L}\end{array}\right)M_{\mathrm{}}\left(\begin{array}{c}e_R\\ \mu _R\\ \tau _R\end{array}\right)+\left(\begin{array}{ccc}\nu _e^T& \nu _\mu ^T& \nu _\tau ^T\end{array}\right)M_\nu \left(\begin{array}{c}\nu _e\\ \nu _\mu \\ \nu _\tau \end{array}\right)+\mathrm{h}.\mathrm{c}.,$$

$`V_{\mathrm{MNS}}`$ can be found from the diagonalizing matrices $`V_{\mathrm{}}`$ and $`V_\nu `$:

$$V_{\mathrm{MNS}}=P_{\mathrm{}}V_{\mathrm{}}V_\nu ^{},$$

where $`P_{\mathrm{}}`$ is a diagonal phase matrix. The unitary matrices $`V_\mathrm{}L`$ and $`V_\nu `$ are found from

$$V_\mathrm{}LM_{\mathrm{}}M_{\mathrm{}}^{}V_\mathrm{}L^{}=\mathrm{diag}(m_e^2,m_\mu ^2,m_\tau ^2),V_\nu M_\nu ^{}M_\nu V_\nu ^{}=\mathrm{diag}(m_1^2,m_2^2,m_3^2).$$

Our first step is to express the physical mixing angles in terms of the parameters of the diagonalizing matrices. For simplicity, we ignore CP violation, so that the mass matrices and, consequently, the diagonalizing matrices are real. (We comment on the effects of CP violating phases in section 5.) Let us define the three unitary matrices

$$\begin{array}{cc}\hfill R_{12}(\theta _{12})& \left(\begin{array}{ccc}c_{12}& s_{12}& 0\\ s_{12}& c_{12}& 0\\ 0& 0& 1\end{array}\right),\hfill \\ \hfill R_{13}(\theta _{13})& \left(\begin{array}{ccc}c_{13}& 0& s_{13}\\ 0& 1& 0\\ s_{13}& 0& c_{13}\end{array}\right),\hfill \\ \hfill R_{23}(\theta _{23})& \left(\begin{array}{ccc}1& 0& 0\\ 0& c_{23}& s_{23}\\ 0& s_{23}& c_{23}\end{array}\right).\hfill \end{array}$$

Then, eq. (1.1) (with $`\delta `$ set to zero) can be rewritten as

$$V_{\mathrm{MNS}}=R_{23}(\theta _{23})R_{13}(\theta _{13})R_{12}(\theta _{12}).$$

We further parametrize the diagonalizing matrices as follows:

$$\begin{array}{cc}\hfill V_\nu ^{}=& R_{23}(\theta _{23}^\nu )R_{13}(\theta _{13}^\nu )R_{12}(\theta _{12}^\nu ),\hfill \\ \hfill V_{\mathrm{}}=& R_{23}(\theta _{23}^{\mathrm{}})R_{13}(\theta _{13}^{\mathrm{}})R_{12}(\theta _{12}^{\mathrm{}}).\hfill \end{array}$$

We limit ourselves to the large class of models where there are no exact relations between the entries in $`M_\nu `$ (up to the fact that it is symmetric, that is, $`(M_\nu )_{ij}=(M_\nu )_{ji}`$) and in $`M_{\mathrm{}}`$. Then the smallness of $`ϵ_{12}`$ and $`s_{13}`$ requires that the following parameters are small:

$$s_{12}^{\mathrm{}},s_{13}^{\mathrm{}},ϵ_{12}^\nu ,s_{13}^\nu 1,$$

where

$$ϵ_{12}^\nu 12(s_{12}^\nu )^2.$$

Evaluating to first order in the small parameters of (2.1), we obtain:

$$\begin{array}{cc}\hfill ϵ_{12}=& ϵ_{12}^\nu +2c_{23}^\nu s_{12}^{\mathrm{}}2s_{23}^\nu s_{13}^{\mathrm{}},\hfill \\ \hfill s_{13}=& s_{13}^\nu s_{23}^\nu s_{12}^{\mathrm{}}c_{23}^\nu s_{13}^{\mathrm{}}.\hfill \end{array}$$

We caution the reader that the sign of the terms that depend on $`s_{12}^{\mathrm{}}`$ and $`s_{13}^{\mathrm{}}`$ is ambiguous. In particular, we approximated $`\mathrm{sin}2\theta _{12}^\nu =1`$, but with the parametrization (2.1) it could equal $`1`$. A full treatment of the sign and phase dependence of the $`s_{12}^{\mathrm{}}`$ contribution to $`ϵ_{12}`$ is given in section 5.

Our next step is to express the parameters of the diagonalizing matrices in terms of the mass matrices. For the charged lepton sector, the expressions can be found in refs. \[20,,21\]. Typically, one finds $`s_{12}^{\mathrm{}}(M_{\mathrm{}})_{12}/(M_{\mathrm{}})_{22}`$ and $`s_{13}^{\mathrm{}}(M_{\mathrm{}})_{13}/(M_{\mathrm{}})_{33}`$. Here, we focus on the neutrino mass matrix with a pseudo-Dirac structure. If there are no exact relations between different entries in $`M_\nu `$, then the most general structure that is consistent with $`ϵ_{12}^\nu 1`$ and $`s_{13}^\nu 1`$ is

$$M_\nu =m\left(\begin{array}{ccc}y_{11}& Y_{12}& Y_{13}\\ Y_{12}& y_{22}& y_{23}\\ Y_{13}& y_{23}& Y_{33}\end{array}\right),$$

where

$$Y_{12}1,y_{ij}1.$$

As concerns $`Y_{13}`$ and $`Y_{33}`$, there are three different options:

$$\begin{array}{cc}\hfill (i)& Y_{13}<1,Y_{33}1,\hfill \\ \hfill (ii)& Y_{13}<1,Y_{33}y_{33}1,\hfill \\ \hfill (iii)& Y_{13}y_{13}1,Y_{33}1.\hfill \end{array}$$

(Explicit examples of models in the literature that realize these options are presented in section 7.) It is also convenient to define the matrix

$$\widehat{M}_\nu =R_{13}^T(\theta _{13}^\nu )R_{23}^T(\theta _{23}^\nu )M_\nu R_{23}(\theta _{23}^\nu )R_{13}(\theta _{13}^\nu ).$$

By definition, it is block diagonal. The requirement that $`ϵ_{12}^\nu 1`$ restricts the form of the (12) block:

$$\widehat{M}_\nu =m\left(\begin{array}{ccc}\delta _1& 1& 0\\ 1& \delta _2& 0\\ 0& 0& Y_3\end{array}\right);|\delta _1|,|\delta _2|1.$$

Both $`\mathrm{\Delta }m_{21}^2/m^2`$ and $`ϵ_{12}^\nu `$ depend on only $`\delta _1`$ and $`\delta _2`$:

$$\begin{array}{cc}\hfill \frac{\mathrm{\Delta }m_{12}^2}{m^2}=& 2|\delta _1^{}+\delta _2|,\hfill \\ \hfill ϵ_{12}^\nu =& \frac{|\delta _2|^2|\delta _1|^2}{2|\delta _1^{}+\delta _2|}.\hfill \end{array}$$

We now present, for the three cases of eq. (2.1), expressions for $`s_{23}^\nu `$, $`s_{13}^\nu `$, $`\delta _1`$ and $`\delta _2`$ to first order in the small parameters $`y_{ij}`$. For case (i), we take $`Y_{12}=1`$ and obtain:

$$\begin{array}{cc}\hfill s_{13}^\nu =& Y_{13}/Y_{33},s_{23}^\nu =(s_{13}^\nu Y_{12}+y_{23})/Y_{33},\hfill \\ \hfill \delta _1=& y_{11}Y_{13}^2/Y_{33},\delta _2=y_{22}.\hfill \end{array}$$

For case (ii), we take $`c_{23}^\nu Y_{12}s_{23}^\nu Y_{13}=1`$ and obtain:

$$\begin{array}{cc}\hfill \mathrm{tan}\theta _{23}^\nu =& Y_{13}/Y_{12},s_{13}^\nu =c_{23}^\nu s_{23}^\nu (y_{33}y_{22})((c_{23}^\nu )^2(s_{23}^\nu )^2)y_{23},\hfill \\ \hfill \delta _1=& y_{11},\delta _2=(c_{23}^\nu )^2y_{22}+(s_{23}^\nu )^2y_{33}2s_{23}^\nu c_{23}^\nu y_{23}.\hfill \end{array}$$

For case (iii), we take $`Y_{12}=1`$ and obtain:

$$\begin{array}{cc}\hfill s_{23}^\nu =& \frac{y_{23}Y_{33}+y_{13}Y_{12}}{Y_{33}^2Y_{12}^2},s_{13}^\nu =\frac{y_{13}Y_{33}+y_{23}Y_{12}}{Y_{33}^2Y_{12}^2},\hfill \\ \hfill \delta _1=& y_{11},\delta _2=y_{22}.\hfill \end{array}$$

We would like to emphasize several points related to the results derived above (some of the statements below were previously made in ref.  in the context of a specific class of textures for the Dirac and Majorana mass matrices in the seesaw model):

(a) Eq. (2.1) implies that flavor models where $`ϵ_{12}^\nu `$ gives the dominant contribution to $`ϵ_{12}`$ can be strongly constrained by a measurement of $`ϵ_{12}`$. Since $`\delta _1`$ and $`\delta _2`$ depend on different entries of $`M_\nu `$, we expect no exact cancellations in their contribution to $`ϵ_{12}^\nu `$. Consequently, one will be able to use the measured size of $`ϵ_{12}`$ to estimate the size of the larger between $`|\delta _1|`$ and $`|\delta _2|`$, and the sign of $`ϵ_{12}`$ to tell which is larger.

(b) Eq. (2.1) implies that observable deviations from maximal mixing in vacuum oscillations, $`1\mathrm{sin}^22\theta _{12}=ϵ_{12}^20`$, can strongly constrain flavor models. For models with a small $`s_{23}^\nu `$, we have

$$ϵ_{12}^2\mathrm{max}(\frac{\delta _2^2}{4},\frac{\delta _1^2}{4},4(s_{12}^{\mathrm{}})^2).$$

If vacuum oscillations show an observable deviation from maximal mixing, say, $`ϵ_{12}^20.1`$, it would be difficult to explain it with a parametrically suppressed $`\delta _2`$, $`\delta _1`$ and $`s_{12}^{\mathrm{}}`$. The accidental factor of sixteen, however, between the $`s_{12}^{\mathrm{}}`$ and the $`\delta _i`$ contributions in eq. (2.1) favors $`s_{12}^{\mathrm{}}`$ as the major source for such a large effect.

(c) Eq. (2.1) reveals interesting relations between the mass hierarchy and the mixing. The parameters $`ϵ_{12}^\nu `$ and $`\mathrm{\Delta }m_{21}^2/m^2`$ are of the same order of magnitude. Therefore, within models where $`ϵ_{12}^\nu `$ gives the dominant contribution to $`ϵ_{12}`$, one will be able to use the measured values of $`ϵ_{12}`$ and $`\mathrm{\Delta }m_{21}^2`$ to estimate the mass scale $`m`$ of the pseudo-Dirac neutrino pair. If the contributions to $`ϵ_{12}`$ related to $`s_{12}^{\mathrm{}}`$ and/or to $`s_{13}^{\mathrm{}}`$ are larger than the contribution related to $`ϵ_{12}^\nu `$, then the relation between $`ϵ_{12}`$ and $`\mathrm{\Delta }m_{21}^2/m^2`$ is lost and, in particular, $`ϵ_{12}\mathrm{\Delta }m_{21}^2/m^2`$ is possible. In any case, if there are no exact relations between entries of the lepton mass matrices, we expect

$$|ϵ_{12}|>\mathrm{\Delta }m_{21}^2/m^2.$$

This relation can be used in two ways. First, measurements of $`|ϵ_{12}|`$ and of $`\mathrm{\Delta }m_{21}^2`$ would give a lower bound on $`m`$. Second, in our framework we have $`m^2<\mathrm{\Delta }m_{\mathrm{atm}}^2`$ and therefore we expect

$$|ϵ_{12}|>\frac{\mathrm{\Delta }m_{\mathrm{sun}}^2}{\mathrm{\Delta }m_{\mathrm{atm}}^2}.$$

This constraint is particularly powerful if the LMA solution (see eq. (1.1)) is realized in Nature since then $`\mathrm{\Delta }m_{\mathrm{sun}}^2/\mathrm{\Delta }m_{\mathrm{atm}}^20.01`$.

3. Radiative Corrections

We consider the effect of radiative corrections on neutrino mass matrices which at a high energy scale $`\mathrm{\Lambda }`$ have the pseudo-Dirac structure (2.1). In particular, we ask whether at some low energy scale $`\mu `$ that is relevant to the solar neutrinos, a significant deviation from maximal mixing could be induced by renormalization group evolution (RGE). We take as our framework the minimal supersymmetric Standard Model. (Our results apply also to the Standard Model, but there the smallness of the charged lepton Yukawa couplings guarantees that the radiative corrections are negligible for our purposes.)

The important parameter for our purposes is related to the Yuakawa coupling of the tau lepton:

$$ϵ_\tau \frac{g_\tau ^2}{(4\pi )^2}(1+\mathrm{tan}^2\beta )\mathrm{ln}\frac{\mathrm{\Lambda }}{\mu }.$$

Here $`g_\tau (1+\mathrm{tan}^2\beta )^{1/2}=m_\tau /\varphi _d`$ is the tau Yukawa coupling in the supersymmetric standard model. The $`ϵ_\tau `$ parameter could be of $`𝒪(0.01)`$ for large $`\mathrm{tan}\beta `$. (Within the SM, one has to replace $`(1+\mathrm{tan}^2\beta )`$ with $`1/2`$, which gives $`ϵ_\tau 10^6`$.) Define a matrix

$$I_\tau =\mathrm{diag}(1,1,1+ϵ_\tau ).$$

We denote the neutrino mass scale at the high scale $`\mathrm{\Lambda }`$ by $`M_\nu ^{\mathrm{HE}}`$. Then, up to universal corrections and negligibly small effects of the muon and electron Yukawa couplings, the renormalized neutrino mass matrix at a scale $`\mu `$ below $`\mathrm{\Lambda }`$ is given in logarithmic approximation by \[23--31\]

$$M_\nu =I_\tau M_\nu ^{\mathrm{HE}}I_\tau .$$

In this section, the parameters that relate to $`M_\nu ^{\mathrm{HE}}`$ and to its diagonalization are denoted, as before, by $`s_{ij}^\nu `$ and $`\delta _i`$. They can be expressed in terms of the entries of $`M_\nu ^{\mathrm{HE}}`$ according to equations (2.1), (2.1) and (2.1). In other words, we have

$$M_\nu ^{\mathrm{HE}}=mR_{23}(\theta _{23})R_{13}(\theta _{13})\left(\begin{array}{ccc}\delta _1& 1& 0\\ 1& \delta _2& 0\\ 0& 0& Y_3\end{array}\right)R_{13}^T(\theta _{13})R_{23}^T(\theta _{23}).$$

The parameters that relate to $`M_\nu `$ and to its diagonalization will be denoted by $`\widehat{s}_{ij}`$ and $`\widehat{\delta }_i`$. The difference between them and the corresponding $`s_{ij}`$ and $`\delta _i`$ parameters vanishes in the limit $`ϵ_\tau 0`$. The main question that we would like to investigate is whether the differences $`\widehat{\delta }_{1,2}\delta _{1,2}`$ are of $`𝒪(ϵ_\tau )`$ or $`𝒪(ϵ_\tau )`$. In the latter case the radiative corrections can be safely neglected.

After a cumbersome but straightforward calculation, we find the following leading corrections:

$$\begin{array}{cc}\hfill \widehat{s}_{23}^\nu s_{23}^\nu =& \frac{1+Y_3^2}{1Y_3^2}(c_{23}^\nu )^2s_{23}^\nu ϵ_\tau +𝒪(s_{13}^\nu ϵ_\tau ),\hfill \\ \hfill \widehat{s}_{13}^\nu s_{13}^\nu =& \frac{2Y_3}{1Y_3^2}c_{23}^\nu s_{23}^\nu ϵ_\tau +𝒪(s_{13}^\nu ϵ_\tau ),\hfill \\ \hfill \widehat{\delta }_1\delta _1=& \frac{4Y_3^2}{Y_3^21}s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau +𝒪(ϵ_\tau ^2),\hfill \\ \hfill \widehat{\delta }_2\delta _2=& \frac{4}{Y_3^21}s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau +𝒪(ϵ_\tau ^2).\hfill \end{array}$$

From the expressions for $`\widehat{\delta }_1`$ and $`\widehat{\delta }_2`$, we obtain:

$$\begin{array}{cc}\hfill \frac{\mathrm{\Delta }m_{21}^2}{m^2}=& 2|\delta _1+\delta _2+4s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau |,\hfill \\ \hfill \widehat{ϵ}_{12}^\nu =& \frac{(\delta _2\delta _1)}{2}\frac{\delta _2+\delta _1+4s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau }{|\delta _1+\delta _2+4s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau |}+2\frac{1+Y_3^2}{1Y_3^2}\frac{(\delta _2+\delta _1)s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau }{|\delta _1+\delta _2+4s_{13}^\nu c_{23}^\nu s_{23}^\nu ϵ_\tau |}.\hfill \end{array}$$

Eq. (3.1) should be compared to eq. (2.1).

We would like to emphasize the following points concerning equations (3.1) and (3.1):

(a) The change in $`s_{23}^\nu `$ is small, $`(\widehat{s}_{23}^\nu s_{23}^\nu )/s_{23}^\nu =𝒪(ϵ_\tau )`$.

(b) The difference $`\widehat{s}_{13}^\nu s_{13}^\nu `$ is suppressed beyond the naive estimate of $`ϵ_\tau `$: for large (small) $`Y_3`$ it is further suppressed by $`1/Y_3`$ ($`Y_3`$) while for $`Y_31`$ it is suppressed by the small $`s_{23}^\nu `$. Effectively then we have $`(\widehat{s}_{13}^\nu s_{13}^\nu )/s_{13}^\nu =𝒪(ϵ_\tau )`$.

(c) Our main result is that the RGE-induced deviation from maximal mixing and mass splitting are suppressed by $`ϵ_\tau s_{13}`$. A combination of the CHOOZ result  and the SuperKamiokande results on atmospheric neutrinos  implies that $`s_{13}`$ is small \[34,,35\]. Consequently, the $`ϵ_\tau s_{13}`$ suppression factor is constrained to be below $`𝒪(10^3)`$. In the limit $`s_{13}=0`$, the leading effects are of order $`ϵ_\tau ^2`$ .

To summarize the results of this section: We find that the contribution from radiative corrections to the deviation from maximal mixing is suppressed beyond the smallness of $`ϵ_\tau `$. The leading corrections to $`ϵ_{12}^\nu `$, $`s_{13}^\nu `$ and $`\mathrm{\Delta }m_{21}^2/m^2`$ are $`𝒪[ϵ_\tau \times \mathrm{max}(s_{13},ϵ_\tau )]`$. Model independently, the size of the effect is not larger than $`𝒪(10^3)`$. (This correction could be important for the mass splitting in the VACL solution.) For the deviation from maximal mixing, this correction is too small to be observed.

4. Non-Canonical Kinetic Terms

Models with horizontal symmetries predict the structure of the mass matrices in the basis where the horizontal charges are well defined. This preferred interaction basis can, in general, be different from the basis where the kinetic terms are canonically normalized \[21,,36,,37\]. In particular, when heavy degrees of freedom related to flavor physics are integrated out, the kinetic terms for the left-handed lepton doublets $`L_i`$ ($`i=1,2,3`$) can be modified to

$$R_{ij}L_i^{}\gamma ^\mu _\mu L_j,$$

where $`R`$ is a hermitian matrix. By rescaling of the $`L_i`$ fields we can bring the diagonal entries of $`R`$ to equal unity,

$$R_{ii}=1,|R_{ij}|1(ij).$$

One can always find a hermitian matrix $`K`$ that brings the kinetic terms back to canonical normalization :

$$K^{}RK=\mathrm{diag}(1,1,1).$$

If the mass matrix in the basis where the kinetic terms are of the non-canonical form (4.1) is $`M_\nu ^{\mathrm{NC}}`$, then the true mass matrix, that is the matrix in the basis with canonical kinetic terms, is given by

$$M_\nu =K^TM_\nu ^{\mathrm{NC}}K.$$

$`K`$ has the form

$$K=\left(\begin{array}{ccc}1& k_{12}& k_{13}\\ k_{12}^{}& 1& k_{23}\\ k_{13}^{}& k_{23}^{}& 1\end{array}\right).$$

For simplicity, we again neglect CP violation and take $`R`$ and, consequently $`K`$, to be real.

We are interested in finding the effects of $`k_{ij}0`$ on the deviation from maximal mixing and on the mass splitting. Our analysis follows similar lines to our study of radiative corrections in the previous section. We take

$$M_\nu ^{\mathrm{NC}}=mR_{23}(\theta _{23}^\nu )R_{13}(\theta _{13}^\nu )\left(\begin{array}{ccc}\delta _1& 1& 0\\ 1& \delta _2& 0\\ 0& 0& Y_3\end{array}\right)R_{13}^T(\theta _{13}^\nu )R_{23}^T(\theta _{23}^\nu ).$$

The parameters that relate to the matrix $`M_\nu `$ of eq. (4.1) are denoted by $`\widehat{s}_{ij}^\nu `$ and $`\widehat{\delta }_i`$.

For the differences between $`\widehat{s}_{ij}^\nu ,\widehat{\delta }_i`$ and the corresponding $`s_{ij}^\nu ,\delta _i`$, we find:

$$\begin{array}{cc}\hfill \widehat{s}_{23}^\nu s_{23}^\nu =& \frac{c_{23}^\nu }{Y_3^21}\left[(Y_3^2+1)k_{23}((c_{23}^\nu )^2(s_{23}^\nu )^2)+2Y_3(k_{12}s_{23}^\nu +k_{13}c_{23}^\nu )\right],\hfill \\ \hfill \widehat{s}_{13}^\nu s_{13}^\nu =& \frac{1}{Y_3^21}\left[(Y_3^2+1)(k_{12}s_{23}^\nu +k_{13}c_{23}^\nu )+2Y_3k_{23}((c_{23}^\nu )^2(s_{23}^\nu )^2)\right],\hfill \\ \hfill \widehat{\delta }_{2,1}\delta _{2,1}=& 2(k_{12}c_{23}^\nu k_{13}s_{23}^\nu ).\hfill \end{array}$$

For the mass difference and deviation from maximal mixing, we obtain

$$\begin{array}{cc}\hfill \frac{\mathrm{\Delta }m_{21}^2}{m^2}=& 2|\delta _1+\delta _2+4(k_{12}c_{23}^\nu k_{13}s_{23}^\nu )|,\hfill \\ \hfill \widehat{ϵ}_{12}^\nu =& \frac{(\delta _2\delta _1)}{2}\frac{\delta _1+\delta _2+4(k_{12}c_{23}^\nu k_{13}s_{23}^\nu )}{|\delta _1+\delta _2+4(k_{12}c_{23}^\nu k_{13}s_{23}^\nu )|}.\hfill \end{array}$$

Terms of order $`s_{13}^\nu k_{ij}`$ contribute with different signs to $`\widehat{\delta }_1`$ and $`\widehat{\delta }_2`$ and modify $`\widehat{ϵ}_{12}^\nu `$ in a qualitatively different way. Quantitatively, however, these effects are negligible.

Before we analyze the consequences eqs. (4.1) and (4.1), we would like to make two comments regarding the size of $`k_{ij}`$ in Froggett-Nielsen type models:

1. In most models of horizontal symmetries, we have

$$k_{ij}<s_{ij}^{\mathrm{}}.$$

2. If two fields carry the same horizontal quantum numbers, $`H(L_i)=H(L_j)`$, we can always define these fields in such a way that $`k_{ij}=0`$.

We would like to emphasize the following points:

(a) The changes in $`s_{23}^\nu `$ and $`s_{13}^\nu `$ are of $`𝒪(k_{ij})`$. The effect can be significant for a small mixing angle. In particular, in the supersymmetric framework, if $`s_{ij}^\nu `$ vanishes because of holomorphy , we expect such zeros to be lifted by these corrections. The mixing angle is, however, still parametrically suppressed.

(b) The leading effect on $`ϵ_{12}^\nu `$ does not change its size but, if $`4(k_{12}c_{23}^\nu k_{13}s_{23}^\nu )`$ is not much smaller than max($`\delta _1,\delta _2`$), can affect its sign. (With CP violating phases, also the size of $`ϵ_{12}^\nu `$ is affected, but in the Froggatt-Nielsen class of models the parametric suppression remains the same.)

(c) In principle, the $`k_{ij}`$-related corrections could enhance $`\mathrm{\Delta }m_{21}^2/m^2`$ compared to $`ϵ_{12}^\nu `$ and therefore avoid (2.1). However, in models where the constraint (4.1) holds, (2.1) is valid.

We conclude that, in general, in models where the kinetic terms are normalized according to (4.1), sign($`ϵ_{12}`$) does not give a useful constraint.

5. An Effective Two Generation Framework

In previous sections we took all the parameters in the Lagrangian to be real. To understand some of the effects of phases, we analyze a two generation model allowing for the most general phase structure.

We parametrize the two generation mixing matrix by

$$V=\left(\begin{array}{cc}c& se^{i\beta }\\ s& ce^{i\beta }\end{array}\right),$$

where $`c\mathrm{cos}\theta _{12}`$, $`s\mathrm{sin}\theta _{12}`$ and the phase $`\beta `$ is physical but does not play a role in oscillation experiments. We parametrize the diagonalizing matrices $`V_{\mathrm{}}`$ and $`V_\nu `$ in the following way:

$$V_{\mathrm{}}=\left(\begin{array}{cc}c_{\mathrm{}}& s_{\mathrm{}}e^{i\beta _{\mathrm{}}}\\ s_{\mathrm{}}& c_{\mathrm{}}e^{i\beta _{\mathrm{}}}\end{array}\right),V_\nu =\left(\begin{array}{cc}c_\nu & s_\nu e^{i\beta _\nu }\\ s_\nu & c_\nu e^{i\beta _\nu }\end{array}\right).$$

Using (2.1), we can express the size of the mixing angle in terms of the four parameters $`s_\nu `$, $`s_{\mathrm{}}`$, $`\beta _\nu `$ and $`\beta _{\mathrm{}}`$:

$$s^2=c_{\mathrm{}}^2s_\nu ^2+s_{\mathrm{}}^2c_\nu ^22e(c_{\mathrm{}}s_{\mathrm{}}c_\nu s_\nu e^{i(\beta _{\mathrm{}}\beta _\nu )}).$$

The charged lepton mass matrix can be written as

$$M_{\mathrm{}}=\left(\begin{array}{cc}m_{11}& m_{12}\\ m_{21}& m_{22}\end{array}\right).$$

Our assumption that $`s_{\mathrm{}}`$ is small requires that a certain combination of entries is small:

$$|s_{\mathrm{}}||\delta _{\mathrm{}}|1,$$

where

$$\delta _{\mathrm{}}\frac{m_{11}m_{21}^{}+m_{12}m_{22}^{}}{|m_{21}|^2+|m_{22}|^2|m_{12}|^2|m_{11}|^2}.$$

The neutrino mass matrix is given by

$$M_\nu =m\left(\begin{array}{cc}\delta _1& 1\\ 1& \delta _2\end{array}\right),|\delta _i|1.$$

We include the effects of radiative corrections, parametrized by

$$I_\mu =\mathrm{diag}(1,1+ϵ_\mu );ϵ_\mu \frac{g_\mu ^2}{(4\pi )^2}(1+\mathrm{tan}^2\beta )\mathrm{ln}\frac{\mathrm{\Lambda }}{\mu },$$

and of non-canonical kinetic terms, parametrized by

$$K=\left(\begin{array}{cc}1& k\\ k^{}& 1\end{array}\right).$$

We find:

$$\begin{array}{cc}\hfill \widehat{\delta }_1=& \delta _1+2k^{},\hfill \\ \hfill \widehat{\delta }_2=& \delta _2+2k.\hfill \end{array}$$

We can now express the mass splitting $`\mathrm{\Delta }m^2/m^2`$ and the deviation from maximal mixing $`ϵ_{12}`$ in terms of the three parameters $`\delta _{\mathrm{}}`$ of eq. (5.1) and $`\widehat{\delta }_1`$ and $`\widehat{\delta }_2`$ of eq. (5.1):

$$\begin{array}{cc}\hfill \frac{\mathrm{\Delta }m^2}{m^2}=& 2|\widehat{\delta }_1^{}+\widehat{\delta }_2|,\hfill \\ \hfill ϵ_{12}=& \frac{|\widehat{\delta }_2|^2|\widehat{\delta }_1|^2}{2|\widehat{\delta }_1^{}+\widehat{\delta }_2|}+2e\left[\delta _{\mathrm{}}\left(\frac{|\widehat{\delta }_1^{}+\widehat{\delta }_2|}{\widehat{\delta }_1^{}+\widehat{\delta }_2}\right)\right].\hfill \end{array}$$

Eq. (5.1) allows us to make (or to re-emphasize) the following points:

1. We again observe the accidental factor of four between the $`\delta _{\mathrm{}}`$ contribution and the $`\widehat{\delta }_i`$ contribution to $`ϵ_{12}`$. A large measured value of $`ϵ_{12}`$ might be a hint then to the size of $`\delta _{\mathrm{}}`$.

2. The usefulness of an experimental determination of sign($`ϵ_{12}`$) depends on the relative size of the small parameters. If $`|\delta _1|,|\delta _2||\delta _{\mathrm{}}|,|k|`$, then sign($`ϵ_{12}`$) depends on the relative size of $`|\delta _1|`$ and $`|\delta _2|`$ which is predicted by the models and a useful constraint can be derived. On the other hand, if $`|\delta _{\mathrm{}}|`$ and/or $`|k|`$ are not smaller than both $`|\delta _1|`$ and $`|\delta _2|`$, then sign($`ϵ_{12}`$) depends on the relative phases between $`\delta _{\mathrm{}}`$ or $`k`$ and $`(\delta _1^{}+\delta _2)`$. Since generic models of approximate horizontal symmetries do not predict the phases, we cannot derive any useful constraint.

6. Abelian horizontal symmetries

The most natural application of our results is in the framework of approximate Abelian horizontal symmetries. To understand the principles of this framework, let us take the simplest example of a horizontal symmetry, $`H=U(1)`$, that is broken by a single small parameter. We denote the breaking parameter by $`\lambda `$ and assign to it a horizontal charge $`1`$. Wherever numerical values are relevant, we take $`\lambda =0.2`$ (so that it is of the order of the Cabibbo angle). Within a supersymmetric framework, the following selection rules apply:

a. Terms in the superpotential that carry an integer $`H`$-charge $`n0`$ are suppressed by $`\lambda ^n`$. Terms with $`n<0`$ vanish by holomorphy.

b. Terms in the Kähler potential that carry an integer $`H`$-charge $`n`$ are suppressed by $`\lambda ^{|n|}`$.

We are particularly interested in the leptonic Yukawa terms:

$$_Y=Y_{ij}^{\mathrm{}}L_i\overline{\mathrm{}}_j\varphi _d+\frac{Y_{ij}^\nu }{M}L_iL_j\varphi _u\varphi _u+\mathrm{h}.\mathrm{c}.,$$

where $`i=1,2,3`$ is a generation index, $`L_i`$ are lepton doublet fields, $`\overline{\mathrm{}}_j`$ are lepton charged singlet fields, and $`\varphi _u`$ and $`\varphi _d`$ are the two Higgs fields. The couplings $`Y_{ij}`$ are dimensionless Yukawa couplings and $`M`$ is a high energy scale. The Yukawa terms come from the superpotential. If the sum of the horizontal charges in a particular term is a positive integer, then the resulting mass term is suppressed as follows:

$$\begin{array}{cc}\hfill (M_{\mathrm{}})_{ij}& \varphi _d\lambda ^{H(L_i)+H(\overline{\mathrm{}}_j)+H(\varphi _d)},\hfill \\ \hfill (M_\nu )_{ij}& \frac{\varphi _u^2}{M}\lambda ^{H(L_i)+H(L_j)+2H(\varphi _u)}.\hfill \end{array}$$

Otherwise, i.e. if the sum of charges is negative or non-integer, the Yukawa coupling vanishes. We use the $``$ sign to emphasize that there is an unknown, independent, order one coefficient for each term (except for the relation $`(M_\nu )_{ij}=(M_\nu )_{ji}`$)).

To understand the possible implications of close-to-maximal mixing on theoretical model building, we imagine that future measurements will give

$$ϵ_{12}\lambda .$$

We examine the consequences of such a constraint on three classes of models in the literature. We find that two classes of models will be excluded, while in the other a unique model is singled out that is consistent with all the requirements.

6.1. Holomorphic zeros

Option (i) of eq. (2.1) has been realized in the framework of supersymmetric Abelian horizontal symmetries, where holomorphic zeros can induce a large 23 mixing together with large 23 mass hierarchy . The horizontal symmetry is $`U(1)_1\times U(1)_2`$ with breaking parameters

$$\lambda _1(1,0),\lambda _2(0,1);\lambda _1\lambda _2\lambda =0.2.$$

We impose four requirements on the model: Large 23 mixing, $`s_{23}1`$; Large hierarchy, $`m_2/m_31`$; $`\nu _1`$ and $`\nu _2`$ form a pseudo-Dirac neutrino, $`\mathrm{\Delta }m_{12}^2m^2`$; A deviation from maximal mixing given by $`ϵ_{12}\lambda `$ (this is the hypothetical constraint from solar neutrino measurements). We find that there is a single set of horizontal charge assignments to the Higgs and lepton doublets that is consistent with all four requirements:

$$\varphi _u(0,0),\varphi _d(0,0),L_1(1,0),L_2(1,1),L_3(0,0).$$

(The choice is single up to trivial shifts by hypercharge which is an exact symmetry of the model, by a Peccei-Quinn symmetry that is an accidental symmetry of the Yukawa sector, and by lepton number if it only changes the overall neutrino mass scale and can be absorbed in the parameter $`M`$, and up to trivial exchange of $`U(1)_1U(1)_2`$.) We find then a unique structure for $`M_\nu `$:

$$M_\nu \frac{\varphi _u^2}{M}\left(\begin{array}{ccc}\lambda ^2& \lambda & \lambda \\ \lambda & 0& 0\\ \lambda & 0& 1\end{array}\right).$$

This matrix is of the form (2.1) with option (i) of eq. (2.1). Therefore, eqs. (2.1) can be applied. To have $`s_{23}1`$ and large enough $`ϵ_{12}`$, together with acceptable charged lepton mass hierarchy, we can choose, for example,

$$\overline{\mathrm{}}_1(3,4),\overline{\mathrm{}}_2(3,2),\overline{\mathrm{}}_3(3,0),$$

which gives

$$M_{\mathrm{}}\varphi _d\left(\begin{array}{ccc}\lambda ^8& \lambda ^6& \lambda ^4\\ \lambda ^7& \lambda ^5& \lambda ^3\\ \lambda ^7& \lambda ^5& \lambda ^3\end{array}\right).$$

The parametric suppression of the physical parameters is then as follows:

$$m_\tau /\varphi _d\lambda ^3,m_\mu /m_\tau \lambda ^2,m_e/m_\mu \lambda ^3,$$

$$\mathrm{\Delta }m_{21}^2/\mathrm{\Delta }m_{23}^2\lambda ^3,\mathrm{\Delta }m_{12}^2/m^2\lambda ,$$

$$s_{23}1,s_{13}\lambda ,ϵ_{12}\lambda .$$

The corrections due to a non-canonical kinetic terms,

$$k_{23}\lambda ^2,k_{12}\lambda ^3,k_{13}\lambda ,$$

leave eqs. (6.1) and (6.1) unchanged.

Within the framework of Abelian horizontal symmetries, it is particularly interesting to find predictions for relations among the physical parameters that are independent of a specific choice of horizontal charges. In the quark sector, there is a single such relation , $`|V_{us}||V_{ub}/V_{cb}|`$. In the lepton sector, when singlet neutrinos play no role, there are three such relations . For the class of models where holomorphic zeros give a pseudo-Dirac structure in the 12 sector but do not affect the parameters that are related to the third generation (the model presented in this subsection belongs to this class), we have the following relations:

$$\begin{array}{cc}\hfill ϵ_{12}& s_{13}/s_{23},\hfill \\ \hfill m/m_3& s_{13}s_{23},\hfill \\ \hfill \mathrm{\Delta }m_{12}^2/m^2& s_{13}/s_{23}.\hfill \end{array}$$

The first of these relations, which involves only mixing angles, can be tested if oscillation experiments measure $`ϵ_{12}`$ and $`s_{13}`$. The last two relations can be combined to give another testable relation.

$$\mathrm{\Delta }m_{\mathrm{sun}}^2/\mathrm{\Delta }m_{\mathrm{atm}}^2s_{13}^3s_{23}.$$

6.2. $`L_eL_\mu L_\tau `$ symmetry

Option (ii) of eq. (2.1) can be realized in a particulary interesting frameowork of approximate $`L_eL_\mu L_\tau `$ symmetry \[40--45\]. The symmetry is broken by small parameters, $`\epsilon _+`$ and $`\epsilon _{}`$ of charges $`+2`$ and $`2`$, respectively . The neutrino mass matrix has the following form:

$$M_\nu \frac{\varphi _u^2}{M}\left(\begin{array}{ccc}\epsilon _{}& 1& 1\\ 1& \epsilon _+& \epsilon _+\\ 1& \epsilon _+& \epsilon _+\end{array}\right).$$

This matrix is of the form (2.1) with option (ii) of eq. (2.1). Therefore, eqs. (2.1) can be applied. We find:

$$m_{1,2}=m\left(1\pm 𝒪[\mathrm{max}(\epsilon _+,\epsilon _{})]\right),m_3=m𝒪(\epsilon _+),$$

$$s_{23}^\nu =𝒪(1),s_{13}^\nu =𝒪(\epsilon _+),ϵ_{12}^\nu =𝒪[\mathrm{max}(\epsilon _+,\epsilon _{})].$$

The charged lepton mass matrix has the form :

$$M_{\mathrm{}}\varphi _d\left(\begin{array}{ccc}\lambda _e& \lambda _\mu \epsilon _{}& \lambda _\tau \epsilon _{}\\ \lambda _e\epsilon _+& \lambda _\mu & \lambda _\tau \\ \lambda _e\epsilon _+& \lambda _\mu & \lambda _\tau \end{array}\right),$$

where the $`\lambda _i`$ allow for a generic approximate symmetry that acts on the SU(2)-singlet charged leptons. Such a symmetry, however, does not affect the relevant diagonalizing angles:

$$s_{23}^{\mathrm{}}=𝒪(1),s_{13}^{\mathrm{}}=𝒪(\epsilon _{}),s_{12}^{\mathrm{}}=𝒪(\epsilon _{}).$$

Eqs. (6.1) and (6.1) lead to the following estimates of the physical mixing angles:

$$s_{23}=𝒪(1),s_{13}=𝒪[\mathrm{max}(\epsilon _+,\epsilon _{})],ϵ_{12}=𝒪[\mathrm{max}(\epsilon _+,\epsilon _{})].$$

We can also estimate the corrections due to non-canonical kinetic terms:

$$k_{23}=0,k_{12},k_{13}=𝒪[\mathrm{max}(\epsilon _+,\epsilon _{})].$$

This leaves the parametric suppression of the physical parameters unchanged.

From eqs. (6.1) and (6.1) we obtain:

$$ϵ_{12}=𝒪(\mathrm{\Delta }m_{\mathrm{sun}}^2/\mathrm{\Delta }m_{\mathrm{atm}}^2).$$

Measurements of $`\mathrm{\Delta }m_{ij}^2`$ and of $`ϵ_{12}`$ can then lead to the exclusion of this model . For example, if $`\mathrm{\Delta }m_{\mathrm{sun}}^2/\mathrm{\Delta }m_{\mathrm{atm}}^210^2`$ and $`ϵ_{12}0.1`$ are established, the model will be excluded.

6.3. Models with two breaking parameters

Option (iii) of eq. (2.1), that is hierarchy of mass splittings without hierarchy of masses, has been realized in the framework of non-anomalous horizontal $`U(1)_H`$ symmetry . The symmetry is broken by two small parameters of opposite charges and equal magnitudes:

$$H(\lambda )=+1,H(\overline{\lambda })=1;\lambda =\overline{\lambda }0.2.$$

Then, the following selection rule applies: terms in the superpotential or in the Kahler potential that carry an (integer) $`H`$-charge $`n`$ are suppressed by $`\lambda ^{|n|}`$. The three neutrino masses are of the same order of magnitude, but the mass splitting between $`\nu _1`$ and $`\nu _2`$ is small if we have

$$\begin{array}{cc}\hfill |H(L_1)+H(L_2)|=& 2|H(L_3)|,\hfill \\ \hfill |H(L_1)+H(L_2)|<& 2|H(L_1)|,2|H(L_2)|.\hfill \end{array}$$

From eq. (2.1) we learn that

$$ϵ_{12}^\nu \mathrm{max}(\lambda ^{2|H(L_1)||H(L_1)+H(L_2)|},\lambda ^{2|H(L_2)||H(L_1)+H(L_2)|}).$$

A typical contribution to $`s_{12}^{\mathrm{}}`$ is given by

$$s_{12}^{\mathrm{}}\lambda ^{|H(L_1)+H(\overline{\mathrm{}}_2)||H(L_2)+H(\overline{\mathrm{}}_2)|}.$$

The important point here is that the first condition in eq. (6.1) requires that $`H(L_1)`$ and $`H(L_2)`$ are either both even or both odd. Eqs. (6.1) and (6.1) give then an upper bound on $`ϵ_{12}`$,

$$ϵ_{12}<\lambda ^2.$$

We conclude that if experiments find $`ϵ_{12}\lambda `$, this type of models will be strongly disfavored.

6.4. Alignment

We would like to make a comment on a particular class of supersymmetric models, where there is no degeneracy among the sleptons and the only mechanism to suppress the supersymmetric contributions to lepton flavor changing decays is alignment \[47,,21,,39\], that is small mixing angles in the neutralino-lepton-slepton couplings. In such models, there is a strong constraint on $`s_{12}^{\mathrm{}}`$ (see e.g. ):

$$\frac{B(\mu e\gamma )}{1.2\times 10^{11}}\left(\frac{s_{12}^{\mathrm{}}}{2\times 10^3}\right)^2\left(\frac{100GeV}{m(\stackrel{~}{\mathrm{}})}\right)^4<1,$$

where $`m(\stackrel{~}{\mathrm{}})`$ is the average slepton mass. In these models it is then particularly difficult to explain a large deviation from maximal mixing. If the dominant source of deviation from maximal mixing is $`s_{12}^{\mathrm{}}`$, we have

$$ϵ_{12}2s_{12}^{\mathrm{}}<4\times 10^3\left(\frac{m(\stackrel{~}{\mathrm{}})}{100GeV}\right)^2.$$

7. Conclusions

If the solar neutrino problem is solved by a large mixing angle solution, and if the mixing is established to be close to maximal but not precisely maximal, then interesting constraints for theoretical model building would arise. Specifically, experiments may measure the size and the sign of the small parameter $`ϵ_{12}`$ defined by

$$\mathrm{sin}^2\theta _{12}\frac{1}{2}(1ϵ_{12}).$$

Flavor models can account for a small $`ϵ_{12}`$ by forcing a pseudo-Dirac structure on the neutrino mass matrix through an approximate horizontal symmetry,

$$M_\nu ^{(2)}m\left(\begin{array}{cc}\delta _1& 1\\ 1& \delta _2\end{array}\right),|\delta _1|,|\delta _2|1.$$

We focus on models where there are no exact relations between different entries of the lepton mass matrices (except for $`(M_\nu )_{ij}=(M_\nu )_{ji}`$). Our main points are the following:

1. The most powerful constraints would arise if $`\delta _1`$ and/or $`\delta _2`$ are the dominant sources of $`ϵ_{12}`$. Then the size of $`|ϵ_{12}|`$ gives the size of the larger between $`|\delta _1|`$ and $`|\delta _2|`$ while the sign of $`ϵ_{12}`$ determines which of the two is larger. Moreover, the mass scale of the solar neutrinos (and not only their mass-squared splitting) can be estimated, $`m^2\mathrm{\Delta }m_{21}^2/|ϵ_{12}|`$.

2. If the dominant source of $`ϵ_{12}`$ is a small angle in the diagonalizing matrix for the charged lepton mass matrix, $`s_{12}^{\mathrm{}}`$, then $`|ϵ_{12}|`$ constrains the size of $`s_{12}^{\mathrm{}}`$ but sign($`ϵ_{12}`$) is unlikely to test the theoretical models. The order of magnitude relation between $`|ϵ_{12}|`$ and $`\mathrm{\Delta }m_{12}^2/m^2`$ is lost, but there is still a useful inequality, $`|ϵ_{12}|>\mathrm{\Delta }m_{21}^2/m`$.

3. Radiative corrections do not play a significant role in $`ϵ_{12}`$ and in $`s_{13}`$. They are supppressed by the tau Yukawa coupling, by a loop factor and by $`s_{13}`$. Consequently, their effect is below the level of $`10^3`$.

4. In models of horizontal symmetries where the kinetic terms are not canonically normalized, sign($`ϵ_{12}`$) depends on the kinetic terms as well and is unlikely to test the models.

It remains to be seen whether future developments in solar neutrino experiments would make a convincing case for the intriguing scenario of pseudo-Dirac neutrinos .

Acknowledgments

I thank John Bahcall, Plamen Krastev and Alexei Smirnov for useful discussions and correspondence. Partial support to this work was provided by the Department of Energy under contract No. DE–FG02–90ER40542, by the Ambrose Monell Foundation, by AMIAS (Association of Members of the Institute for Advanced Study), by the Israel Science Foundation founded by the Israel Academy of Sciences and Humanities, and by the Minerva Foundation (Munich).

References

relax J.N. Bahcall, P.I. Krastev and A.Yu. Smirnov, Phys. Rev. D60 (1999) 093001, hep-ph/9905220; Phys. Lett. B477 (2000) 401, hep-ph/9911248. relax M.C. Gonzalez-Garcia, P.C. de Holanda, C. Pena-Garay and J.W.F. Valle, Nucl. Phys. B573 (2000) 3, hep-ph/9906469. relax G.L. Fogli, E. Lisi, D. Montanino and A. Palazzo, Phys. Rev. D62 (2000) 013002, hep-ph/9912231. relax C. Giunti, M.C. Gonzalez-Garcia and C. Pena-Garay, Phys. Rev. D62 (2000) 013005, hep-ph/0001101. relax S.T. Petcov, Phys. Lett. B110 (1982) 245. relax C.N. Leung and S.T. Petcov, Phys. Lett. B125 (1983) 461. relax G.C. Branco, W. Grimus and L. Lavoura, Nucl. Phys. B312 (1989) 492. relax S.T. Petcov and A.Yu. Smirnov, Phys. Lett. B322 (1994) 109, hep-ph/9311204. relax P. Binetruy, S. Lavignac, S. Petcov and P. Ramond, Nucl. Phys. B496 (1997) 3, hep-ph/9610481. relax M. Leurer, Y. Nir and N. Seiberg, Nucl. Phys. B398 (1993) 319, hep-ph/9212278. relax E.Kh. Akhmedov, Phys. Lett. B467 (1999) 95, hep-ph/9909217. relax E.Kh. Akhmedov, G.C. Branco and M.N. Rebelo, Phys. Rev. Lett. 84 (2000) 3535, hep-ph/9912205. relax J.W. Valle and M. Singer, Phys. Rev. D28 (1983) 540. relax V. Barger, S. Pakvasa, T.J. Weiler and K. Whisnant, Phys. Lett. B437 (1998) 107, hep-ph/9806387. relax Z. Maki, M. Nakagawa and S. Sakata, Prog. Theo. Phys. 28 (1962) 870. relax A.H. Guth, L. Randall and M. Serna, JHEP 9908 (1999) 018, hep-ph/9903464. relax A. de Gouvea, A. Friedland and H. Murayama, hep-ph/9910286; hep-ph/0002064. relax A. Friedland, hep-ph/0002063. relax C.D. Froggatt and H.B. Nielsen, Nucl. Phys. B147 (1979) 277. relax L.J. Hall and A. Rasin, Phys. Lett. B315 (1993) 164. relax M. Leurer, Y. Nir and N. Seiberg, Nucl. Phys. B420 (1994) 468, hep-ph/9410320. relax G. Dutta and A.S. Joshipura, Phys. Rev. D51 (1994) 3838, hep-ph/9405291. relax P.H. Chankowski and Z. Pluciennik, Phys. Lett. B316 (1993) 312, hep-ph/9306333. relax K.S. Babu, C.N. Leung and J. Pantaleone, Phys. Lett. B319 (1993) 191, hep-ph/9309223. relax J. Ellis and S. Lola, Phys. Lett. B458 (1999) 310, hep-ph/9904279. relax J.A. Casas, J.R. Espinosa, A. Ibarra and I. Navarro, Nucl. Phys. B556 (1999) 3, hep-ph/9904395; Nucl. Phys. B569 (2000) 82, hep-ph/9905381; JHEP 9909 (1999) 015, hep-ph/9906281. relax R. Barbieri, G.G. Ross and A. Strumia, JHEP 9910 (1999) 020, hep-ph/9906470. relax E. Ma, J. Phys. G25 (1999) L97, hep-ph/9907400 . relax N. Haba, Y. Matsui, N. Okamura and M. Sugiura, Prog. Theor. Phys. 103 (2000) 145, hep-ph/9908429. relax P.H. Chankowski, W. Krolikowski and S. Pokorski, Phys. Lett. B473 (2000) 109, hep-ph/9910231. relax K.R.S. Balaji, A.S. Dighe, R.N. Mohapatra and M.K. Parida, Phys. Rev. Lett. 84 (2000) 5034, hep-ph/0001310. relax CHOOZ collaboration, M. Apollonio et al., Phys. Lett. B420 (1998) 397, hep-ex/9711002; Phys. Lett. B466 (1999) 415, hep-ex/9907037. relax Super-Kamiokande Collaboration, Y. Fukuda et al., Phys. Rev. Lett. 81 (1998) 1562, hep-ex/9807003. relax V. Barger, T.J. Weiler and K. Whishnant, Phys. Lett. B440 (1998) 1, hep-ph/9807319. relax G.L. Fogli, E. Lisi, A. Maronne and G. Scioscia, Phys. Rev. D59 (1999) 033001, hep-ph/9808205. relax E. Dudas, S. Pokorski and C.A. Savoy, Phys. Lett. B356 (1995) 45, hep-ph/9504292. relax G. Eyal and Y. Nir, Nucl. Phys. B528 (1998) 21, hep-ph/9801411. relax Y. Grossman, Y. Nir and Y. Shadmi, JHEP 9810 (1998) 007, hep-ph/9808355. relax Y. Grossman and Y. Nir, Nucl. Phys. B448 (1994) 30, hep-ph/9502418. relax R. Barbieri, L.J. Hall, D. Smith, A. Strumia and N. Weiner, JHEP 9812 (1998) 017, hep-ph/9807235. relax A.S. Joshipura and S.D. Rindani, hep-ph/9811252. relax P.H. Frampton and S.L. Glashow, Phys. Lett. B461 (1999) 95, hep-ph/9906375. relax R.N. Mohapatra, A. Perez-Lorenzana and C.A. de S. Pires, Phys. Lett. B474 (2000) 355, hep-ph/9911395. relax K. Cheung and O.C.W. Kong, Phys. Rev. D61 (2000) 113012, hep-ph/9912238. relax Q. Shafi and Z. Tavartkiladze, Phys. Lett. B482 (2000) 145, hep-ph/0002150. relax Y. Nir and Y. Shadmi, JHEP 9905 (1999) 023, hep-ph/9902293. relax Y. Nir and N. Seiberg, Phys. Lett. B309 (1993) 337, hep-ph/9304307. relax J. Feng, Y. Nir and Y. Shadmi, Phys. Rev. D61 (2000) 113005, hep-ph/9911370. relax M.C. Gonzalez-Garcia, Y. Nir, C. Pena-Garay and A. Smirnov, to appear.


text:
# Stability analysis of the limit-from𝐷-dimensional nonlinear Schrödinger equation with trap and two- and three-body interactions

## Abstract

Considering the static solutions of the $`D`$dimensional nonlinear Schrödinger equation with trap and attractive two-body interactions, the existence of stable solutions is limited to a maximum critical number of particles, when $`D2`$. In case $`D=2`$, we compare the variational approach with the exact numerical calculations. We show that, the addition of a positive three-body interaction allows stable solutions beyond the critical number. In this case, we also introduce a dynamical analysis of the conditions for the collapse.
PACS: 03.75.Fi; 47.20.Ky; 02.30.Jr; 31.75.Pf
Keywords: Nonlinear Schrödinger Equation; trapped two and three-body atomic systems; multidimensional systems

Recent experiments on Bose Einstein Condensation (BEC)  have brought great attention to its theoretical formulation. Atomic traps are effectively described by the Ginzburg-Pitaevskii-Gross (GPG) formulation of the nonlinear Schrödinger equation (NLSE) , which includes two-body interaction. When the atoms have negative two-body scattering lengths, a formula for the critical maximum number of atoms was presented in ref. . In ref. , the formulation was extended in order to include the effective potential originated from the three-body interaction. In this case, in three-dimensions, it was shown that a kind of first order phase-transition occurs. In this connection, as also considered in the motivations given in , it is relevant to observe that recently it was reported the possibility of altering continuously the two-body scattering length, from positive to negative values, by means of an external magnetic field . Within such perspective, the two-body binding energy can be close to zero, and one can approach the so-called Efimov limit, which corresponds to an increasing number of three-body bound states . Near this limit, nontrivial consequences can occur in the dynamics of the condensate, such that one should also consider three-body effects in the effective nonlinear potential.

In the present work, we study the critical number of atoms in arbitrary $`D`$dimensions, using a variational procedure; and also by an exact numerical approach in the case of dimension $`D=2`$. The $`D`$dimensional NLSE, with attractive two-body interactions, was previously analyzed in models of plasma and light waves in nonlinear media . The collapse conditions, in this case, were investigated without  and with  the harmonic potential term. In case of $`D=`$3, it was shown that a repulsive nonlinear three-body interaction term can extend considerably the critical limit for the existence of stable solutions .

Motivated by the observed high interest in stable solutions for arbitrary $`D`$, we look for variational solutions in a few significant cases ($`D=`$1,2,4 and 5) not previously considered, when a three-body interaction term, parametrized by $`\lambda _3`$, is added to the effective non-linear interaction that contains a two-body attractive term. Our analysis also shows that, as in case of $`D=3`$, a kind of first-order phase-transition can occur when $`D4`$, for certain cases of $`\lambda _30`$. In the present paper, we have also considered the approach given in , in order to study the stability conditions in the case of arbitrary $`D`$, when the non-linear interaction contains two (attractive) and three-body terms.

In order to obtain an analytical approach and verify the validity of the variational Ritz method, we consider in detail the case of $`D=2`$, with and without the three-body term, comparing the variational results with exact numerical calculations for some relevant physical observables. In this case, we also discuss how the method given in  can be extended in order to approach analytically the exact value for the total energy.

By extending the GPG formalism from three to $`D`$ dimensions, including two  and three-body interactions in the effective non-linear potential , we obtain

$`i\mathrm{}{\displaystyle \frac{d\psi }{dt}}=\left[{\displaystyle \frac{\mathrm{}^2}{2m}}^2+{\displaystyle \frac{m\omega ^2r^2}{2}}+\lambda _2|\psi |^2+\lambda _3|\psi |^4\right]\psi ,`$ (1)
where $`\psi \psi (\stackrel{}{r},t)`$ is the wave-function normalized to the number of atoms $`N`$, $`\omega `$ is the frequency of the trap harmonic potential and $`m`$ is the mass of the atom. $`\lambda _2`$ and $`\lambda _3`$ are, respectively, the strength of the two- and three-body effective interaction, given in a $`D`$dimensional space. $`r|\stackrel{}{r}|`$ is the hyperradius, such that $`\stackrel{}{r}_{i=1}^Dr_i\widehat{e}_i`$ and $`_{i=1}^D\widehat{e}_i\frac{}{r_i}`$ ($`\widehat{e}_i`$ is the unit vector, with $`i=1,2,\mathrm{}D`$).

The stationary solutions for the chemical potential $`\mu `$ are given by

$$i\mathrm{}\frac{d\psi }{dt}=\mu \psi .$$
(2)
Considering the general solution of eq.(1), $`i\mathrm{}{\displaystyle \frac{d\psi }{dt}}={\displaystyle \frac{\delta }{\delta \psi ^{}}}`$, one can obtain the total energy $`E`$:

$`E`$ $`=`$ $`{\displaystyle d^D\stackrel{}{r}},\mathrm{with}`$ (3)
$``$ $``$ $`{\displaystyle \frac{\mathrm{}^2}{2m}}\left|\psi \right|^2+{\displaystyle \frac{m\omega ^2r^2}{2}}\left|\psi \right|^2+{\displaystyle \frac{\lambda _2}{2}}|\psi |^4+{\displaystyle \frac{\lambda _3}{3}}|\psi |^6.`$ (4)
Here we consider only attractive two-body interaction, which is more interesting in the case of trapped atoms. For $`D=`$3, $`\lambda _24\pi \mathrm{}^2|a|/m`$, where $`a`$ is the two-body scattering length and $`m`$ is the mass of the atom. In the case of arbitrary $`D`$, $`\lambda _2`$ has dimension of energy times $`L^D`$, where $`L`$ is a length scale in such space. However, a convenient redefinition of the wave-function in terms of dimensionless variables will absorb this constant, as will be shown.

Our study will be concentrated on the ground state for a spherically symmetric potential. We first consider the case of $`\lambda _3=0`$, using a variational procedure, with a trial Gaussian wave-function for $`\psi (\stackrel{}{r})`$, normalized to $`N`$, given by

$$\psi _{var}(\stackrel{}{r})=\sqrt{N}\left(\frac{1}{\pi \alpha ^2}\frac{m\omega }{\mathrm{}}\right)^{D/4}\mathrm{exp}\left[\frac{r^2}{2\alpha ^2}\left(\frac{m\omega }{\mathrm{}}\right)\right],$$
(5)
where $`\alpha `$ is a dimensionless variational parameter. From eq. (3), the corresponding expression for the total variational energy can be expressed as

$`E_{var}`$ $`=`$ $`\mathrm{}\omega {\displaystyle \frac{N}{\nu }}_{var},`$ (6)
$`_{var}`$ $``$ $`\nu \left({\displaystyle \frac{D}{4\alpha ^2}}+{\displaystyle \frac{D\alpha ^2}{4}}\right){\displaystyle \frac{\nu ^2\mathrm{\Omega }_D}{4(2\pi )^{D/2}\alpha ^D}}+{\displaystyle \frac{G_3}{6\pi ^D}}{\displaystyle \frac{\nu ^3\mathrm{\Omega }_D^2}{3^{D/2}\alpha ^{2D}}},`$ (7)
where $`\mathrm{\Omega }_D`$ is the solid angle in $`D`$ dimensions,

$`\mathrm{\Omega }_D{\displaystyle \frac{2\pi ^{D/2}}{\mathrm{\Gamma }(D/2)}},G_3{\displaystyle \frac{\lambda _3}{2(\lambda _2)^2}}\mathrm{}\omega ,`$ (9)
$`\mathrm{and}\nu {\displaystyle \frac{N}{\mathrm{\Omega }_D}}{\displaystyle \frac{2\lambda _2}{\mathrm{}\omega }}\left({\displaystyle \frac{m\omega }{\mathrm{}}}\right)^{D/2}.`$ (10)
By using dimensionless variables, $`\stackrel{}{x}\sqrt{m\omega /\mathrm{}}\stackrel{}{r}`$, we redefine the wave-function $`\psi `$ as

$$\varphi (\stackrel{}{x})\sqrt{\frac{2|\lambda _2|}{\mathrm{}\omega }}\psi (\stackrel{}{r}),$$
(11)
such that

$$|\varphi (\stackrel{}{x})|^2d^D\stackrel{}{x}=N\left[\frac{2|\lambda _2|}{\mathrm{}\omega }\right]\left(\frac{m\omega }{\mathrm{}}\right)^{D/2}=\nu \mathrm{\Omega }_D.$$
(12)
The dimensionless equation corresponding to eq. (1), can be rewritten as

$$\left[\left(\underset{1}{\overset{D}{}}\frac{d^2}{dx_i^2}+x_i^2\right)|\varphi |^2+G_3|\varphi |^42\beta \right]\varphi =0,$$
(13)
where $`\beta \mu /(\mathrm{}\omega )`$ is the dimensionless chemical potential. From eqs. (11) and (5), the trial wave-function can be written as

$$\varphi _{var}(x)\sqrt{\nu \mathrm{\Omega }_D}\left(\frac{1}{\pi \alpha ^2}\right)^{D/4}\mathrm{exp}\left(\frac{x^2}{2\alpha ^2}\right),$$
(14)
The variational results, obtained by using the above expressions can be extended analytically to non-integer values of the dimension $`D`$. Minimization of the energy \[eq. (6)\], with respect to $`\alpha ^2`$, is done numerically by sweeping over $`\alpha ^2`$ values. The results for the energy and the chemical potential are shown in Fig. 1. For each value of $`D`$, one can observe a critical number of atoms, $`N_c`$, related to the critical parameter $`\nu _c`$, only when $`D2`$. This critical limit corresponds to the cusps in the upper plot of Fig.1 and is also observed using exact numerical calculation for $`D=`$3. It is also interesting to note that for $`D>2`$ there are two branches of solutions for $`_{var}`$ and $`\beta `$, one stable and the other unstable. In the energy, the lower branch corresponds to stable solutions (minima), while the upper one gives unstable solutions (maxima).

The case with $`D=2`$ is particularly interesting, as no unstable solutions exist and there are stable solutions only for $`\nu <2`$, such that $`\nu _c=2`$. For $`D=2`$, the minimization of eq. (6) in respect to $`\alpha ^2`$ leads to

$$_{var}=\nu \sqrt{1\frac{\nu }{2}}.$$
(15)
The behavior of $`\nu `$, and the corresponding critical limits, as one alters the dimension $`D`$, has other curious particular results. For example, the critical limit $`\nu _c`$ has a minimum for $`D=3`$ ($`\nu _c^{(D)}\nu _c^{(3)}`$ for all $`D`$).

In conclusion of this part of our work, considering arbitrary $`D`$ with $`\lambda _3=0`$, there are no stable solutions for eq. (1), if the wave-function $`\varphi (x)`$, given by eq. (12), is normalized to $`\nu >\nu _c`$. Fig. 1 shows that this restriction is strongest for $`D=3`$: $`\nu _c`$ is a minimum when compared with $`\nu _c`$ for $`D3`$. This is a relevant result, considering that $`\nu `$ is directly proportional to the number of atoms. Also, it is observed that $`\nu _c`$ increases very fast for $`D>3`$.

Next, we also solve equation (13) exactly employing the shooting and Runge-Kutta methods, and compare the results with the corresponding variational ones. In this case, we consider only the particular interesting case of $`D=2`$, with $`\lambda _3=0`$. The results are shown in Fig. 2, for the chemical potential $`\beta `$, the total energy $``$, mean-square-radius $`x^2`$, and the central density $`|\varphi (0)|^2`$. In order to numerically solve eq. (13), in the $`s`$wave, we first write it in terms of the single variable $`x\sqrt{(x_1^2+x_2^2)}`$ and consider the following boundary conditions: $`\varphi ^{}(0)=0`$ (where $``$ stands for the derivative with respect to $`x`$) and $`\varphi (x)`$ $`C\mathrm{exp}(x^2/2+[\beta 1]\mathrm{ln}(x))`$ when $`x\mathrm{}`$, where $`C`$ is a constant to be determined. As observed in Fig.2, the critical limit $`\nu _c=2`$ obtained analytically using the variational approach should be compared with $`\nu _c=1.862`$, obtained by exact numerical calculation. This critical limit was first obtained by Weinstein , in a non-linear approach with two-body term, without the trapping potential. The coincidence of the value with our exact calculation is due to the fact that at the critical limit the mean square radius goes to zero.

We have also compared the results obtained by the variational approach with the exact numerical one, in the case of $`D=2`$, for several values of the three-body interaction term (positive and negative), as shown in Figs. 3 and 4. In Fig. 3 we have the exact numerical approach and in Fig. 4 we have the corresponding variational results. By comparing the results we have for $`D=2`$ (shown in Figs. 3 and 4) with the ones obtained in ref.  for $`D=3`$, we should observe that no first order phase-transition exists in two dimensions. As observed in refs. , for $`D=3`$, a first-order phase-transition can occur in trapped condensed states with negative two-body scattering length, when a repulsive three-body (quintic) term is added in the Hamiltonian. As shown in Figs. 3 and 4, with $`G_3`$ positive the range of stability for the number of atoms $`N`$ can be increased indefinitely; with $`G_3`$ negative this range is reduced.

We can analyze the collapse conditions using “the virial theorem” approach . The mean square radius, $`r^2`$, of a $`D`$dimensional condensate, is given by

$`{\displaystyle \frac{d^2r^2}{dt^2}}+4\omega ^2r^2=`$ (16)
$`{\displaystyle \frac{1}{m}}\left[4H+\lambda _2(D2)|\psi |^2+{\displaystyle \frac{4\lambda _3}{3}}(D1)|\psi |^4\right],`$ (17)
where

$$𝒪\frac{1}{N}d^D\stackrel{}{r}\psi ^{}(\stackrel{}{r},t)𝒪\psi (\stackrel{}{r},t)$$
(18)
and $`H=E/N`$. When $`\lambda _3=0`$ we obtain the equation derived in .

We can also write the eq. (16) in dimensionless units, as it was done in eqs. (11-13):

$`{\displaystyle \frac{d^2x^2}{d\tau ^2}}+4x^2={\displaystyle \frac{4}{\nu }}+2f(\tau ),`$ (19)
where

$`f(\tau ){\displaystyle \frac{\lambda _2}{|\lambda _2|}}{\displaystyle \frac{D2}{4}}|\varphi |^2+G_3{\displaystyle \frac{D1}{3}}|\varphi |^4.`$ (20)
Using the initial conditions for $`x^2`$ and $`dx^2/d\tau `$, where, for simplicity, we assume $`dx^2/d\tau =0`$, the solution of eq. (19) is given by

$`x^2`$ $`=`$ $`{\displaystyle \frac{}{\nu }}+\left[x^2|_0{\displaystyle \frac{}{\nu }}\right]\mathrm{cos}(2\tau )`$ (21)
$`+`$ $`{\displaystyle _0^\tau }f(\tau ^{})\mathrm{sin}(2(\tau \tau ^{}))𝑑\tau ^{}.`$ (22)
The stability regions and the estimates for the collapse time can be obtained from the analysis of this solution, like as performed for the case $`\lambda _3=0`$ in . Let us analyze the dynamics when $`D=2`$. In this case, $`\lambda _2`$ does not appear explicitly in $`x^2`$ and $`f(\tau )`$ also does not depend on this parameter:

1. For a positive $`G_3`$, negative $`\lambda _2`$ and $`>0`$ we observe that $`x^2`$ cannot be zero and the condensate is stable. The mean square radius of the condensate oscillates in time around a finite value. This is confirmed by the numerical simulations (see Figs. 3 and 4).
2. For a negative $`G_3`$, positive $`\lambda _2`$ an analysis of stability like the one performed in ref.  shows that
a) When the total energy $`<0`$, the condensate is unstable and the wavefields collapse in a finite time at any initial conditions;
b) When $`>0`$, as the function $`f(\tau )`$ is negative, the contribution of the integral term for $`\tau <\pi `$ is negative. Then, we found the collapse condition as

$$x^2|_02\frac{}{\nu }.$$
(23)

The same kind of analysis, for $`D>2`$, is much involved in the present approach, as the sign of the function $`f(\tau )`$ is not fixed at opposite signs for the parameters $`\lambda _2`$ and $`\lambda _3`$.

Some information about the dynamics of the collapse can also be obtained by using the techniques based on integral inequalities . For instance, when $`D=2`$, we can estimate the three-body term contribution in $`E`$, following the procedure given in

$$d^2\stackrel{}{r}|\psi |^6C_2\left(d^2\stackrel{}{r}\frac{|\psi |^2}{2m}\right)^2\left(d^2\stackrel{}{r}|\psi |^2\right)=C_2K^2N,$$
(24)
where $`K`$ is the kinetic energy and $`C_2`$ is defined from the minimization of the functional

$$𝒥=\frac{\left(d^2\stackrel{}{r}|\psi |^2\right)^2\left(d^2\stackrel{}{r}|\psi |^2\right)}{d^2\stackrel{}{r}|\psi |^6}.$$
(25)
Combining with the corresponding estimate for $`d^2\stackrel{}{r}|\psi |^4`$, we obtain $`E>E(K)`$, where

$$E(K)=K+\frac{\omega ^2N^2}{4K}+\frac{\lambda _2}{2}C_1NK+\frac{\lambda _3}{3}C_2K^2N.$$
(26)
When $`\lambda _3=0`$ we get the equation derived in . Equation (26) should be compared with the corresponding variational expression (6), where the kinetic energy is given by $`K=N\mathrm{}\omega /(2\alpha ^2)`$ and $`\alpha `$ is the width of the cloud. As we see, the expression for the energy (26) is very similar to the obtained by the variational approach. However, (26) is valid for arbitrary time and describes the nonstationary dynamics. By using the variational expression (upper limit) for the ground-state, and the right-hand-side of eq. (26) (lower limit), we can approach analytically the exact solution for the total energy

$$E(K)<E<E_{var}.$$
(27)
For a more deep insight to the problem of stability, we need to obtain the values of the constants $`C_1`$ and $`C_2`$. This problem requires a generalization of the method suggested by Weinstein in , to be considered in a future work.

We should observe that exact numerical results, when $`G_3=0`$, have already been considered in refs.  (for $`D=1`$ and $`D=3`$), in  (for $`D=`$ 3), and in  (for $`D=`$2). In , for $`D=3`$, it was also considered the case with $`G_30`$, and shown a kind of first-order phase-transition in the condensate. In the present work, we have extended the variational formalism, in case $`G_30`$, for an arbitrary $`D`$dimension. In the following Figs. 5, 6 and 7, we present our results for the chemical potential as a function of $`\nu `$, for a set of given values of $`G_3`$, in case of $`D=`$ 1, 4 and 5. As one can observe in Fig. 5, even in case of $`D=1`$ one can reach a critical maximum limit for $`\nu `$, when $`G_3`$ is enough negative. For $`D=`$ 4 and 5 (Figs. 6 and 7), we observe similar picture of first-order phase-transition occurring for some specific values of $`G_3`$.

In conclusion, in the present work, we first studied the stability and the critical number of atoms in arbitrary $`D`$dimensions using a variational procedure, for the case we have two-body (attractive) and three-body contributions. This part extends a previous analysis done in refs. . Next, we considered in more detail the case $`D=2`$.. We compared the variational results with exact numerical calculations for the chemical potential, total energy, mean-square-radius and density. Finally, we extended numerically the approach for $`D=2`$, including an effective three-body interaction term. We studied the sensibility of the critical numbers with respect to corrections in the non-linear interaction. The effective interaction considered in the equation contains a trapped harmonic interaction, and two nonlinear terms, proportional to the density $`|\psi |^2`$ (due to first-order two-body interaction) and to $`|\psi |^4`$ (due to first-order three-body interaction). We also verified, by a variational procedure, that a critical number of particles exists only for $`D2`$, when the nonlinear term of the NLSE contains just the cubic term. In case of $`D=1`$, a critical maximum number of atoms can exist with the addition of a negative quintic term (three-body) in the NLSE. In all cases where the number of atoms is limited, we observed that the addition of a positive $`|\psi |^4`$ allows stable solutions beyond the critical number. We also introduced an analysis of the collapse conditions, using “the virial theorem” approach given in . The dynamics of the collapse was discussed in terms of the techniques developed in . In particular, we showed how the exact energy can be approached in the case of $`D=2`$ with two and three-body term contributions.

Acknowledgments We are grateful to Jordan M. Gerton for the suggestions and careful reading of the manuscript. This work was partially supported by Fundação de Amparo à Pesquisa do Estado de São Paulo and Conselho Nacional de Desenvolvimento Científico e Tecnológico.


text:
# 1 Dark Matter

## 1 Dark Matter

The problem of missing or Dark Matter, namely that there is insufficient material in the form of stars to hold galaxies and clusters together, has been known since the pioneering work of Bessel, Zwicky and most recently Rubin .
The existence of non-luminous Dark Matter was first inferred in 1984 by Fredrich Bessel from gravitational effects on positional measurements of Sirius and Procyon. In 1933, Zwicky concluded that the velocity dispersion in Rich Clusters of galaxies required 10 to 100 times more mass to keep them bound than could be accounted for by luminous galaxies themselves.
Finally, Trimble  noted that the majority of galactic rotation curves, at large radii, remain flat or even rise well outside the radius of the luminous astronomical object.
The missing Dark Matter has been traditionally explained in terms of Dark Matter Halo’s , although none of the Dark Matter Halo models have been very successful in explaining the experimental data .
This paper will describe the missing matter (Dark Energy) in terms of a Cosmological Constant which leads to a constant energy density.
The experimental determination of galactic velocity rotation curves (VRC) has been one of the most important approaches used to estimate the ”local” mass (energy) density of the Universe. Several sets of data from VRC’s will be analysed and the contribution due to the Cosmogical Constant determined.

## 2 Constaints on the value of the Cosmological Constant

It is interesting to estimate the allowed range of values for the Cosmological Constant within the constraints of General Relativity and observational astronomy, (for a comprehensive review, see Bahcall et.al. ).
Starting from a General Relativity point of view, the Friedman energy equation is given by:

$$1=\frac{8\pi G_N}{3}\frac{\rho _{matter}}{H^2}\frac{kc^2}{R^2H^2}+\frac{c^2\mathrm{\Lambda }}{3H^2},$$
(1)
where the Hubble Constant is denoted by $`H`$, the curvature term by $`k`$ and $`G_N`$ denotes the Newton gravitational constant. Eq.(1) can be rewritten as

$$1=\mathrm{\Omega }_m+\mathrm{\Omega }_k+\mathrm{\Omega }_\mathrm{\Lambda }$$
(2)
Here the relative contributions to the energy density of the universe are given by, the mass, curvature and Cosmological Constant.
If we assume that the curvature contribution is small:

$$1=\mathrm{\Omega }_{Matter}+\mathrm{\Omega }_\mathrm{\Lambda }$$
(3)
In order to satisfy equation (3), it was surprising to discover that only a narrow range of values for the observed Cosmological Parameters were allowed. A ”reasonable” set of parameters consistent with observation are:

$$H_O=100Kms^1Mpc^1,\rho _{matter}=5\times 10^{30}gcm^3,\frac{\mathrm{\Omega }_\mathrm{\Lambda }}{\mathrm{\Omega }_{matter}}=4.3$$
(4)
and $`\mathrm{\Omega }_{Matter}+\mathrm{\Omega }_\mathrm{\Lambda }=1.4`$ (here we assume a small value for the curvature $`0.4`$. (For an authoritative review of the matter/energy sources of the universe, see Turner ).
It was found that observational constraints placed upon the range of values for the cosmological parameters lead to a surprisingly narrow range of possible values for the Cosmological Constant, the range being given by:

$$10^{56}<|\mathrm{\Lambda }|<5\times 10^{55}cm^2.$$
(5)
## 3 Experimental Results

It was shown , within the Weak Field Approximation, that the Cosmological Constant at large radii could be determined from galactic velocity rotation curves. This contribution is given by:

$$v_\mathrm{\Lambda }^2(r)=v_{obs}^2(r)v_{mass}^2(r),\mathrm{leading}\mathrm{to},$$
(6)
$$v_\mathrm{\Lambda }^2/r=\frac{c^2\mathrm{\Lambda }r}{3},\mathrm{at}\mathrm{large}\mathrm{r}$$
(7)
The results obtained by this analysis are shown in Table 1.
The experimental values obtained for the Cosmological Constant fall within the range determined from General Relativity and observational constraints. While the initial results are promising, a thorough and systematic analysis of galactic rotation curves needs to be undertaken in order to confirm the trend.
Previous results  reported for the value of the Cosmological Constant were 100 to 1000 times the ”allowed value”. This systematic error arose for two main reasons: the first by not taking the gradient of the curves at sufficiently large radii and the second by the lack of access to ”real” experimental data leading to crude data analysis.
The results presented in this paper suffer from the second problem, i.e. all the gradients were obtained from the data in the published literature and not from raw experimental data i.e. M33 Corbelli $`\&`$ Salucci , NGC 3198  and all others from .
However experience has taught us that a cursory look at rotation curves will determine which galaxies are candidates for explanation by a Cosmological Constant and which are not. Galaxies where the velocity rotation curve remains flat or rises at large radii, are immediate candidates. NGC 3198 is a good example, whereas others such as M33  has clearly not relaxed to the Cosmological background, even at many times the galactic radii. A full explanation for M33 has to be sought in a different direction.
Finally, a simple calculation of the effective mass density due to the Cosmological Constant in NGC 3198,

$$\rho _{eff}=\frac{c^2\mathrm{\Lambda }}{4\pi G_N}$$
(8)
leads to a value of $`5.4\times 10^{29}gcm^3`$ which is comparable to the HI mass density  at the outer disk of NGC 3198 galaxy. This is further confirmation that the Cosmological Constant effect can be seen at galactic scale lengths.

## 4 Accelerating or Decelerating Universe?

Recently there has been great interest in the Type Ia Supernovae results of Perlmutter et al  which suggest that the universe is accelerating.
In this section we will show that the Weak Field Approximation coupled with galactic velocity rotation curve data inevitably lead to a negative Cosmological Constant.
The equation for the VRC is given  by 333In Ref. eq.(15), there was a typographical sign error for one of the terms and also the negative pressure effect associated with $`\mathrm{\Lambda }`$ was not fully appreciated.

$$\frac{v^2}{r}=\frac{Gm}{r^2}+\frac{c^2\mathrm{\Lambda }}{3}r$$
(9)
We note that Eq.(9) is only strictly true for small and large radii, however it will serve to illustrate our arguments.
Using the Newtonian limit of Einstein field equations we derived equation (9). It is important to realize that the Cosmological Constant obeys the equation of state given by,

$$P_\mathrm{\Lambda }=c^2\rho _\mathrm{\Lambda },$$
(10)
where $`P_\mathrm{\Lambda }`$ is the pressure term due to $`\mathrm{\Lambda }`$. Taking the Newtonial limit in the absence of matter, $`T_{\mu \nu }=0`$, the differential equation for the static Newtonian potential becomes

$$^2\mathrm{\Phi }=c^2\mathrm{\Lambda }$$
(11)
leading to,

$$\rho _{eff}=\rho _\mathrm{\Lambda }+\frac{3P_\mathrm{\Lambda }}{c^2}=2\rho _\mathrm{\Lambda }$$
(12)
If we arbitrary set $`\mathrm{\Phi }=0`$ at the origin, then in spherical coordinates (11) has the solution $`\mathrm{\Phi }=\frac{c^2\mathrm{\Lambda }}{6}r^2`$ . Thus, the Cosmological Constant leads to the following correction to the Newtonial potential

$$\mathrm{\Phi }=\frac{Gm}{r}\frac{c^2\mathrm{\Lambda }}{6}r^2$$
(13)
At small galactic radii the velocity versus radius contribution is well known and follows Newtonian physics. For large radii a negative Cosmological Constant gives a positive contribution to the VRC which is what is actually observed. On the other hand the effect of a positive Cosmological Constant would be to lower the rotation curve below that due to matter alone.
The above simple argument, based on observational astronomy, allows only a negative Cosmological Constant as a possible explanation for the galactic velocity rotation curve data. This is in clear disagreement with the Type Ia supernovae results . However, given the uncertainties in the determination of the deceleration parameter, $`q_0`$, derived from supernovae data  the approach outlined above has certain merits worth consideration.
In summary these are , the Cosmological Constant is determined from $`direct`$ measurement unlike the Supernovae results, the experimentally determined value is the correct order of magnitude as that required from cosmological constraints, and finally a negative Cosmological Constant is consistent, and indeed a natural physical explanation , for the observed galactic velocity rotation curve data.
Finally, observations of global clusters of stars constrain the age of the universe and consequently place an observational limit on a negative Cosmological Constant  of ,

$$|\mathrm{\Lambda }|2.2\times 10^{56}cm^2.$$
(14)
Note, the Cosmological Constant derived from global cluster constraints is in agreement with the experimentally determined value derived from galactic velocity rotation curve data.

### 4.1 Experimental Tests-Dark Matter Halo vs Cosmological Constant

It would be of some interest if it was possible to experimentally distinguish between the contribution of Dark Matter Halo’s and Dark Energy (Cosmological Constant) to galactic rotation curves.
We know that Dark Matter predicts a variation of mass at large radii given by ,

$$M_{DM}(r)r$$
(15)
while for Dark Energy due to a Cosmological Constant,

$$M_\mathrm{\Lambda }(r)r^3[\rho _\mathrm{\Lambda }+(3P_\mathrm{\Lambda }/c^2)].$$
(16)
With these different types of predictive variations it should be possible to design experimental tests to distinguish between the two phenomena.

## 5 Quark Hadron Phase transition

In this section which is of more speculative nature, working within the Extended Large Number Hypothesis, and using the experimentally determined Cosmological Constant, we will demonstrate how the energy density for the Quark - Hadron can be estimated.
However, it is useful to put into context the significance of the Cosmological Constant for many seemingly disparate branches of Physics. The figure 1 below shows the Cosmological Constant at the epicentre of Physics.
The diagram demonstrates a dichotomy whereby several branches of Physics need a non-zero Cosmological Constant in order to explain key physical phenomena, whilst in others a non-zero value presents a fundamental problem.

It is also noted here that while fundamental theories of Particle Physics such as the Standard Model, Quantum Field Theory and String Theory have many major predictive successes they all have problem with a high vacuum energy density. On the other hand while the Extended LNH is formulated from a naive theory  it appears to correctly predict the correct vacuum energy density and other cosmological parameters. The Extended LNH relates the value of the Cosmological Constant to the effective mass given by:

$$|\mathrm{\Lambda }|=\frac{G_N^2m_{eff}^6}{h^4}=\frac{c^6L_s^4}{h^6}m_{eff}^6$$
(17)
Matthews  pointed out that when using the Extended LNH to determine today’s cosmological parameters, the mass of the proton originally suggested by Dirac should be replaced by the energy density of the last phase transition of the Universe : Quark - Hadron.
Note that in equation (17) there are no free parameters, $`L_s`$ is normalised to the gravitational constant and corresponds to the fundamental length of String Theory.
Using equation (17) and the Cosmological Constant determined from NGC 5033, the effective mass is given by

$$m_{Effective}=332MeV$$
(18)
We will associate this value with the Quark - Hadron phase transition energy. (Other experimentally determined Cosmological Constant data give $`m_{QH}`$ in the range 295 - 410 MeV). The experimentally determined value within the LNH predicts the correct order of magnitude for the phase transition.
The above result poses the question that it might be possible to gain insights on the quantum mechanical origin of the Universe, as Dirac  suggested, from direct observation of the present day Universe.
Finally, does the Cosmological Constant provide the key to the integration of the various Physics disciplines as Figure 1 suggests?

## 6 Discussion

Analysis of the galactic rotation curves show that the missing Galactic Dark Matter can be explained in terms of a Cosmological Constant.
This contribution can be considered a prime candidate for the ”Dark Energy” which is smoothly distributed throughout space, and contributes approximately $`70\%`$ to the mass/energy of the Universe .
However, in order to support this thesis for the Cosmological Constant, thorough and systematic analysis of galactic velocity rotation data needs to take place.
It was shown how, within the Weak Field Approximation, that VRC data inevitably lead to a negative value for the Cosmological Constant in direct disagreement with the type Ia Supernovae data. Nevertheless, given the uncertainties in determining the deceleration parameter $`q_0`$ , from the redshift-magnitude Hubble diagram using Type Ia supernovae as standard candles, we believe our approach is worth further consideration.
The experimental values determined for the Cosmological Constant are shown to lie within an acceptable range. These values, used within the Extended Large Number Hypothesis, predict values for the Quark-Hadron phase transition energy in the range 295-410 MeV.
It would be remarkable, if proved correct, that the Cosmological Constant could be directly determined from the analysis of galactic velocity rotation curves.
Equally remarkable, if proved correct, is the idea that astronomical observations can shed light on the last quantum mechanical phase transition of the Universe, namely the Quark-Hadron.

## 7 Acknowledgements

We would like to thank Paolo Salucci for invaluable discussions on Cosmological and Astronomical aspects related to this work and Alexander Love for useful comments on the manuscript. We also thank J. Hargreaves and D. Bailin for suggestions and useful discussions.
We also would like to thank Deja Whitehouse for proof reading this document.
George Kraniotis was supported for this work by PPARC.


text:
# Evidence of first-order transition between vortex glass and Bragg glass phases in high-𝑇_c superconductors with point pins: Monte Carlo simulations

\[

## Abstract

Phase transition between the vortex glass and the Bragg glass phases in high-$`T_\mathrm{c}`$ superconductors in $`\stackrel{}{B}\stackrel{}{c}`$ is studied by Monte Carlo simulations in the presence of point pins. A finite latent heat and a $`\delta `$-function peak of the specific heat are observed, which clearly indicates that this is a thermodynamic first-order phase transition. Values of the entropy jump and the Lindemann number are consistent with those of melting transitions. A large jump of the inter-layer phase difference is consistent with the recent Josephson plasma resonance experiment of Bi2Sr2CaCu2O8+y by Gaifullin et al.

\] Vortex states in high-$`T_\mathrm{c}`$ superconductors have been intensively studied experimentally and theoretically . Because of large fluctuations owing to high transition temperature and strong anisotropy, the flux-line lattice (FLL) melts at much lower temperatures than those predicted by Abrikosov’s mean-field theory. The FLL melting is a thermodynamic first-order phase transition. In pure systems, the melting line stretches up to a high magnetic field as large as $`H_{\mathrm{c}_2}`$. However, all experiments show that first-order melting lines terminate at much lower magnetic fields . Complicated phase diagrams are obtained experimentally in Bi2Sr2CaCu2O8+y (BSCCO)  and YBa2Cu3O7-δ (YBCO) , and it is believed that effects of impurities are essential in real materials. For example, vacancies of oxygen atoms, which play the role of point pins to flux lines, cannot be excluded completely even in crystals of highest quality.

Fisher et al.  studied the Ginzburg-Landau model in a random potential, and proposed the so-called vortex glass (VG) phase before experimental studies. Giamarchi and Doussal  pointed out that the Bragg glass (BG) phase can exist in weak fields. In this phase, correlation functions decay in power laws , and the structure factor shows a triangular Bragg pattern. Accordingly, a phase transition between these two glass phases may be observed when the magnetic field is sweeped across the phase boundary (see Fig. 1). The existence of the VG–BG phase transition in vortex systems with point pins shows a sharp contrast to the phase diagram of pure systems. Presuming a first-order phase transition, the VG–BG phase boundary was evaluated phenomenologically . However, physical properties around this phase boundary have not been clarified in experiments and numerical calculations until recently. The shape of the phase boundary seems to depend on observed quantities, and some experiments even suggest a crossover rather than a phase transition. Although some simulations  gave phase diagrams similar to that of Giamarchi and Doussal, numerical accuracy of these studies was not good enough to distinguish phase transitions and crossovers. The stability of the VG phase was studied by Kawamura  including the screening effect.

Quite recently, Gaifullin et al. observed  a large jump of the inter-layer phase difference on the VG–BG phase boundary of BSCCO by the Josephson plasma resonance experiment. They claimed that their observation is the evidence of a first-order phase transition. In the present Letter, we show more direct evidence of the first-order phase transition on the VG–BG boundary by large-scale Monte Carlo simulations. That is, a finite latent heat and a $`\delta `$-function peak of the specific heat are observed. Sharp jumps of the inter-layer phase difference and the averaged fluctuations of flux lines are also obtained.

In order to clarify vortex states and phase transitions of high-$`T_\mathrm{c}`$ superconductors in the presence of point pins, we start from the three-dimensional anisotropic, frustrated XY model on a simple cubic lattice . Effects of point pins are introduced into the model by randomly-distributed weakly-coupled plaquettes in the $`ab`$ plane. Since a vortex sitting on a plaquette costs an energy proportional to the couplings surrounding it, flux lines tend to penetrate plaquettes with weaker couplings in order to reduce such loss of energies. The Hamiltonian of our model is given by

$`H`$ $`=`$ $`{\displaystyle \underset{i,jab\mathrm{plane}}{}}J_{ij}\mathrm{cos}\left(\varphi _i\varphi _jA_{ij}\right)`$ (2)
$`{\displaystyle \frac{J}{\mathrm{\Gamma }^2}}{\displaystyle \underset{m,nc\mathrm{axis}}{}}\mathrm{cos}\left(\varphi _m\varphi _n\right),`$
$`A_{ij}`$ $`=`$ $`{\displaystyle \frac{2\pi }{\mathrm{\Phi }_0}}{\displaystyle _i^j}𝐀^{(2)}d𝐫^{(2)},`$ (3)
with the periodic boundary condition along all the directions. Couplings in the $`ab`$ plane are given by $`J_{ij}=bJ`$ ($`0<b<1`$) on the weakly-coupled plaquettes, and $`J_{ij}=J`$ otherwise. The density and the strength of point pins are controlled by the probability of weakly-coupled plaquettes, $`p`$, and the parameter $`b`$, respectively (see Fig. 2). The pinning energy is of order of $`(1b)J`$. A uniform magnetic field is applied along the $`c`$ axis, and its strength is proportional to the averaged number of flux lines per plaquette, $`f`$. Here we concentrate on the model with $`L_x=L_y=50`$ and $`L_c=40`$. This system size is large enough to describe the melting transition in the pure system ($`b=1`$) .

In our model, we have four adjustable parameters: the anisotropy constant $`\mathrm{\Gamma }`$, the density of flux lines $`f`$, the density of point pins $`p`$, and the strength of point pinning $`b`$. In order to investigate the VG–BG transition, we vary $`b`$, while fix the temperature at $`T=0.06J/k_\mathrm{B}`$ and other parameters at $`\mathrm{\Gamma }=20`$, $`f=1/25`$ and $`p=0.003`$. In other words, material parameters of the bulk system and the number and positions of point pins are not changed during the simulations. As will be shown later, this temperature is low enough for the study of the VG–BG phase boundary. Typical Monte Carlo steps (MCS) with the Metropolis algorithm are $`34\times 10^7`$ MCS for equilibration, and $`0.51\times 10^7`$ MCS for measurement.

First, the $`b`$ dependence of the internal energy $`e`$ and the specific heat $`C`$ per flux line per $`ab`$ plane is displayed in Fig. 3. Clearly, the internal energy shows a sharp jump at the transition point $`b^{}=0.895\pm 0.005`$ with a latent heat per flux line per $`ab`$ plane $`Q2.3\times 10^2J`$, and a $`\delta `$-function peak of the specific heat also occurs at the same parameter. These two facts indicate that the VG–BG transition is a thermodynamic first-order phase transition. From this latent heat, the entropy jump at $`b^{}`$ is estimated as

$$\mathrm{\Delta }S=Q/T0.38k_\mathrm{B},$$
(4)
which is comparable to the experimental value in the melting transition of YBCO, $`\mathrm{\Delta }S0.5k_\mathrm{B}`$ .

Second, the $`b`$ dependence of the inter-layer phase difference, $`\mathrm{cos}(\varphi _n\varphi _{n+1})`$, is plotted in Fig. 4. This quantity is related to the Josephson energy per phase variable $`e_\mathrm{J}`$ and the anisotropy constant $`\mathrm{\Gamma }`$ as

$$\mathrm{cos}(\varphi _n\varphi _{n+1})=e_\mathrm{J}\mathrm{\Gamma }^2/J,$$
(5)
and a small change of $`e_\mathrm{J}`$ is magnified in this quantity in extremely anisotropic systems. This quantity also jumps sharply at $`b^{}=0.895\pm 0.005`$, and the value of the jump at $`b^{}`$, $`\mathrm{\Delta }_{\mathrm{PD}}0.12`$, is as large as the experimental value, $`\mathrm{\Delta }_{\mathrm{PD}}0.2`$ . Moreover, the ratio of the jump of the Josephson energy to the latent heat is given by $`\mathrm{\Delta }e_\mathrm{J}/(Qf)0.34`$, which means that the latent heat is equally distributed to all the directions in the VG–BG phase transition in extremely anisotropic systems.

Third, the Lindemann number is evaluated directly . The deviation $`u`$ of a flux line is measured in each $`ab`$ plane from the projection of the mass center of the flux line, and averaged over all the flux lines and the $`ab`$ planes. Then, the Lindemann number $`c_\mathrm{L}`$ is given by

$$c_\mathrm{L}=\underset{bb^{}+0}{lim}u^2^{1/2}/a_0,$$
(6)
with the lattice constant of the triangular FLL, $`a_0=(2/\sqrt{3})^{1/2}/f^{1/2}`$. The $`b`$ dependence of $`u^2^{1/2}/a_0`$ is shown in Fig. 5, and we have $`c_\mathrm{L}0.28`$. This value is almost equal to the one obtained in the FLL melting of pure systems .

Finally, we go into some details of the present simulations. The system with $`b=0.90`$ is calculated at first. Simulations are started from a random configuration at a very high temperature, and the system is gradually cooled down to $`T=0.06J/k_\mathrm{B}`$. During the cooling process, the first-order melting transition characterized by a discontinuous appearance of the helicity modulus along the $`c`$ axis, $`\mathrm{{\rm Y}}_c`$ , takes place at $`T_\mathrm{m}0.079J/k_\mathrm{B}`$, which corresponds to the vortex liquid (VL)–BG phase transition. Then, the strength of point pinning $`b`$ is varied. Since the quantity $`\mathrm{{\rm Y}}_c`$ is proportional to the superfluid density, the region with finite $`\mathrm{{\rm Y}}_c`$ is superconducting. This quantity is nonvanishing for all the values of $`b`$ shown in Figs. 35 at $`T=0.06J/k_\mathrm{B}`$, and therefore the phase transition investigated in the present Letter is not the VL–BG one, but the VG–BG one. Equilibration in systems with point pins is much slower than that in pure systems, and only one sample can be taken for calculations at present. Nevertheless, the results obtained in the present Letter are quite clear-cut and consistent with experiments. Thus, the small number of random sampling does not seem serious. Since positions of point pins are independent in each $`ab`$ plane, the number of $`ab`$ planes, $`L_c=40`$, would be large enough for averaging effects of point pins.

Although we have concentrated on the VG–BG transition for a single density of point pins $`p`$ in the present Letter, we have also investigated the VL–BG and VL–VG transitions for various $`p`$, and obtained the overall phase diagram in the $`p`$$`T`$ plane. The structure of the $`p`$$`T`$ phase diagram is similar to that of the $`B`$$`T`$ phase diagram. Experimentally, the increase of $`p`$ corresponds to the repeated irradiation of electrons or protons, and our $`p`$$`T`$ phase diagram is consistent with recent experiments . Details of this study will be reported elsewhere .

In conclusion, the first thermodynamic evidence of the first-order transition between the vortex glass (VG) and the Bragg glass (BG) phases has been obtained in high-$`T_\mathrm{c}`$ superconductors in the presence of point pins. A finite latent heat and a $`\delta `$-function peak of the specific heat are observed by large-scale Monte Carlo simulations of the three-dimensional anisotropic, frustrated XY model with randomly-distributed weakly-coupled plaquettes. The entropy jump derived from the latent heat is nearly equal to those in the melting transition of YBCO. The Lindemann number evaluated from fluctuations of flux lines, $`c_\mathrm{L}0.28`$, is reasonable for the first-order phase transition. The inter-layer phase difference also shows a sharp jump on the VG–BG phase boundary. This property is consistent with the Josephson plasma resonance experiment of BSCCO by Gaifullin et al.

The present authors would like to thank Prof. Y. Matsuda for communications. Numerical calculations were performed on Numerical Materials Simulator (NEC SX-4) at National Research Institute for Metals, Japan.