1 INTRODUCTION

The predictions of the Standard Model (SM) are effectively confirmed by a vast number of experiments. The first discrepancy with SM was discovered in 1998 as a result of observations of atmospheric neutrino oscillations [1], proving the existence of finite neutrino masses. By adding the right neutrinos to SM while maintaining the renormalizability of the model, it is possible to generate Dirac neutrino masses. When the non-renormalizable \(d=5\) interaction of left leptons with the Higgs boson is turned on, masses of Majorana neutrinos appear [2]. A plausible hierarchy of masses is suggested by the seesaw mechanism involving both Dirac and Majorana neutrinos [3]. The mechanism of neutrino mass generation by the spontaneous chiral symmetry breaking is also discussed [4]. Only in isolated cases the theory encounters challenges in adjacent SM sectors, the interpretation of measurements of the muon’s anomalous magnetic moment [5] or the decays of beauty mesons [6].

In the searches for physics beyond SM, neutrinos occupy a unique place. The mass matrices of Majorana and Dirac neutrinos have 9 and 7 free parameters, respectively, already in the simplest case of three generations [7, 8]. This is roughly half of the total number of SM free parameters. The elements of neutrino mass matrix may contain a great deal of information on potential SM extensions.

Cross-sections of processes involving neutrinos, probabilities of oscillations and decay rates are usually expressed in terms of the Pontecorvo–Maki–Nakagawa–Sakata (PMNS) mixing matrix and the diagonal neutrino masses. Significant progress in the experimental determination of the elements of the mixing matrix and obtaining restrictions on neutrino masses has been achieved over the past decade due to the analysis of neutrino oscillations, including data from atmospheric experiments, accelerator data on the modes of appearance and disappearance of neutrinos, and data from reactor experiments [9, 10]. The elements of the mass matrix, in principle, are accessible to experimental determination in terms of the mixing matrix and the diagonal neutrino masses, but contain uncertainties associated with arbitrariness of the absolute mass scale. An indirect approach to the determination of the mass matrix can lead to loss of information at intermediate stages of statistical analysis. In this paper, we show that the processes associated with neutrinos are explicitly described in terms of the neutrino mass matrix, without using a mixing matrix. The existence of such representations makes it possible to fit the elements of neutrino mass matrix directly from experimental data, which should provide an advantage when conducting statistical analysis and, potentially, distinguish the Majorana mass matrix from the Dirac one.

Another purpose of this paper is to point out effective techniques of matrix algebra that significantly simplify the work with neutrinos of three or more flavors. The paper uses extensively Sylvester’s theorem on the representation of matrix functions [11], Frobenius covariants, which determine projection operators in terms of the matrix under consideration, and trace identities [12, 13], which provide explicit solutions to the Faddeev–LeVerrier algorithm [14, 15]. The described techniques allow to effectively construct the neutrino dispersion laws. The cases of \(n\leq 4\) flavors of Majorana neutrinos in a matter of arbitrary composition and \(n\leq 4\) flavors of Dirac neutrinos in a baryonic matter are solved in quadratures.

The outline of the paper is as follows. Neutrino oscillations in a vacuum are discussed in the next section. Neutrino oscillations in the medium are covered in Section 3. Majorana neutrinos and Dirac neutrinos are considered separately. The neutrino dispersion laws are derived, and it is demonstrated that the mass matrix explicitly describes the physical amplitudes of neutrino propagation in all cases considered. The neutrino dispersion law in the medium takes an unusual form: the minimum neutrino energy is achieved at a finite momentum. Other processes covered briefly in Section 4 include tritium \(\beta\) decay, electron capture, double-\(\beta\) decays, and quasi-elastic neutrino scattering. In all circumstances, probabilities are explicitly expressed in terms of the neutrino mass matrix rather than the mixing matrix.

2 NEUTRINO OSCILLATIONS IN VACUUM

2.1 Majorana Neutrino

Let us first consider the problem of oscillations of Majorana neutrinos. The Lagrangian of free Majorana neutrinos is written as

$$\mathcal{L}_{M}=\frac{1}{2}\bar{\nu}_{\alpha}i\hat{\nabla}\nu_{\alpha}-\frac{1}{2}\bar{\nu}_{\alpha}\left(M_{\alpha\beta}+i\gamma_{5}N_{\alpha\beta}\right)\nu_{\beta}.$$
(1)

Here \(\nu_{\alpha}\) are the real four-component spinors of Majorana neutrinos, the gamma matrices are defined in the Majorana representation where \(\gamma^{*}_{\mu}=-\gamma_{\mu}\), \(\gamma^{T}_{\mu}=-\gamma^{\mu}\), \(\gamma^{*}_{5}=\gamma^{T}_{5}=-\gamma_{5}\), \(\gamma^{\dagger}_{5}=\gamma_{5}\) and \(\nu_{\alpha}^{c}=\nu_{\alpha}^{*}=\nu_{\alpha}\), \(\bar{\nu}_{\alpha}^{c}={\nu}_{\alpha}^{T}\gamma_{0}\). The left and right neutrinos are defined by \(\nu_{\alpha L}=\Pi_{-}\nu_{\alpha}\) and \(\nu_{\alpha R}=\Pi_{+}\nu_{\alpha}\), respectively, where \(\Pi_{\pm}=(1\pm\gamma_{5})/{2}\) are the projection operators onto states with the definite chirality. The transformation properties of the spinor bilinear forms determine the mass term symmetry. We denote by \(\Gamma_{A}\) the matrices \(1,i\gamma_{5},\gamma^{\mu},\gamma_{5}\gamma^{\mu},\sigma^{\mu\nu}\), defined in the same way as in [16]. Taking into account \(\gamma_{0}\Gamma_{A}^{\dagger}\gamma_{0}=\Gamma_{A}\) and the fermionic character of the fields, for both Dirac and Majorana neutrinos

$$\bar{\nu}_{\alpha}\Gamma\nu_{\beta}=-\bar{\nu}_{\beta}^{c}\gamma^{0}(\gamma^{0}\Gamma_{A})^{T}\nu_{\alpha}^{c}=\eta\bar{\nu}_{\beta}^{c}\Gamma_{A}\nu_{\alpha}^{c}.$$
(2)

The sign \(\eta\) is read from the equation (cf. Appendix B in [17])

$$-\gamma^{0}\left(\gamma^{0}\left(\begin{matrix}1\\ i\gamma_{5}\\ \gamma^{\mu}\\ \gamma_{5}\gamma^{\mu}\\ \sigma^{\mu\nu}\end{matrix}\right)\right)^{T}=\left(\begin{matrix}1\\ i\gamma_{5}\\ -\gamma^{\mu}\\ \gamma_{5}\gamma^{\mu}\\ -\sigma^{\mu\nu}\end{matrix}\right).$$
(3)

The conditions \(\bar{\nu}_{\alpha}\nu_{\beta}=\bar{\nu}_{\beta}\nu_{\alpha}\) and \(\bar{\nu}_{\alpha}i\gamma_{5}\nu_{\beta}=\bar{\nu}_{\beta}i\gamma_{5}\nu_{\alpha}\) valid for Majorana fields imply that the mass matrices \(M\) and \(N\) are symmetric. The hermiticity of the Lagrangian ensures the hermiticity of \(\mathbb{M}\) and \(\mathbb{N}\). Since they are symmetric, they are real. We define for \(n\) flavors a \(4n\times 4n\) real matrix

$$\hat{\mathbb{M}}=\mathbb{M}+i\gamma_{5}\mathbb{N}$$
(4)

and \(n\times n\) complex symmetric matrices

$$\mathbb{M}_{\pm}=\mathbb{M}\pm i\mathbb{N}.$$
(5)

One can verify that \(\hat{\mathbb{M}}^{\dagger}=\hat{\mathbb{M}}^{T}=\gamma_{0}\hat{\mathbb{M}}\gamma_{0}=\mathbb{M}-i\gamma_{5}\mathbb{N}\) and \(\mathbb{M}_{+}^{\dagger}=\mathbb{M}^{*}_{+}=\mathbb{M}_{-}\). We also use notations \(\mathbb{M}_{R}=\mathbb{M}_{+}\), \(\mathbb{M}_{L}=\mathbb{M}_{-}\) and write \(\mathbb{M}_{\pm 1}\) for \(\mathbb{M}_{\pm}\) and \(\Pi_{\pm 1}\) for \(\Pi_{\pm}\). The matrices \(\mathbb{M}\) and \(\mathbb{N}\) are specified by \(2\times(n(n+1)/2)\) real parameters. Phase rotations of Majorana spinors \(\nu_{\alpha}\to\exp(i\gamma_{5}\phi_{\alpha})\nu_{\alpha}\) keep the kinetic term in the Lagrangian and the weak neutral current \(j_{\mu}^{\textrm{NC}}\) unchanged but modifies the mass matrix and introduces the phase multipliers into the weak charged current \(j_{\mu}^{\textrm{CC}}\). The phase multipliers in \(j_{\mu}^{\textrm{CC}}\) are absorbed further by the charged leptons. The number of remaining free parameters in the mass matrix equals \(n^{2}\).

In the representation of two-component Weyl spinors \(\nu_{\alpha L}\), the Lagrangian (1) has the form

$$\mathcal{L}_{M}=\frac{1}{2}\bar{\nu}_{\alpha L}i\hat{\nabla}\nu_{\alpha L}-\frac{1}{2}\bar{\nu}_{\alpha L}M_{\alpha\beta}^{M}(\nu_{\beta L})^{c}+\textrm{h.c.},$$
(6)

where \(\nu_{\alpha R}=(\nu_{\alpha L})^{c}\) and \(M^{M}_{\alpha\beta}=M_{\alpha\beta}+iN_{\alpha\beta}\) is a symmetric complex \(n\times n\) matrix \(\mathbb{M}_{R}\) defined above. It is useful to identify phases of the mass matrix elements: \(M^{M}_{\alpha\beta}=|M^{M}_{\alpha\beta}|e^{i\phi^{\prime}_{\alpha\beta}}.\) The chiral rotations \(\nu_{\alpha}\to\exp(-i\gamma_{5}\phi^{\prime}_{\alpha\alpha}/2)\nu_{\alpha}\) enable the diagonal elements of \(\mathbb{N}\) to be set to zero and the diagonal elements of \(\mathbb{M}_{R}\) to be considered as real and positive. For three generations

$$\mathbb{M}_{R}=\left|\left|\begin{matrix}M^{M}_{ee}&|M^{M}_{e\mu}|e^{i\phi_{e\mu}}&|M^{M}_{e\tau}|e^{i\phi_{e\tau}}\\ |M^{M}_{e\mu}|e^{i\phi_{e\mu}}&M^{M}_{\mu\mu}&|M^{M}_{\mu\tau}|e^{i\phi_{\mu\tau}}\\ |M^{M}_{e\tau}|e^{i\phi_{e\tau}}&|M^{M}_{\mu\tau}|e^{i\phi_{\mu\tau}}&M^{M}_{\tau\tau}\end{matrix}\right|\right|,$$
(7)

where \(\phi_{\alpha\beta}=\phi^{\prime}_{\alpha\beta}-(\phi^{\prime}_{\alpha\alpha}+\phi^{\prime}_{\beta\beta})/2\).

The free Lagrangian (1) is invariant with respect to the charge conjugation: \(\nu_{\alpha}\rightarrow\nu_{\alpha}^{c}=\nu_{\alpha}^{\dagger}\), however, due to the special type of the mass matrix, it is non-invariant with respect to \(P\)-parity transformations \(\nu_{\alpha}\rightarrow\mathcal{P}\nu_{\alpha}\mathcal{P}^{-1}=i\gamma_{0}\nu_{\alpha}\), and time reversal \(\nu_{\alpha}\rightarrow\mathcal{T}\nu_{\alpha}\mathcal{T}^{-1}=\gamma_{0}\gamma_{5}\nu_{\alpha}\), When weak interactions are turned on, only the combined CPT symmetry survives. \(CP\) invariance would call for \(\mathbb{N}=0\) and forbid phase rotations with the \(\gamma_{5}\) matrix, which would leave \(n(n+1)/2\) free parameters in the mass matrix.

The Feynman propagator has the form

$$\hat{\mathbb{S}}_{F}(x-y)=-i\langle 0|T{\nu}(x)\otimes\bar{\nu}(y)|0\rangle$$
$${}=\int\frac{d^{4}p}{(2\pi)^{4}}e^{-ip(x-y)}\hat{\mathbb{S}}_{F}(p),$$
(8)

where

$$\hat{\mathbb{S}}_{F}(p)=\frac{1}{\hat{p}-\hat{\mathbb{M}}}$$
(9)

is determined by the quadratic form of the Lagrangian. The matrix \(\hat{p}-\hat{\mathbb{M}}\) has the dimension \(4n\). The energy eigenvalues are grouped by pairs \(\pm\sqrt{m_{i}^{2}+\mathbf{p}^{2}}\) for different masses \(m_{i}\) (\(i=1,..,n\)). The additional doubling is due to the spin degeneracy. Using the identity \(A^{-1}=B(AB)^{-1}=(BA)^{-1}B\) for \(A=\hat{p}-\hat{\mathbb{M}}\) and \(B=\hat{p}+\hat{\mathbb{M}}^{\dagger}\) and \(\hat{p}\hat{\mathbb{M}}=\hat{\mathbb{M}}^{\dagger}\hat{p}\), we obtain

$$\hat{\mathbb{S}}_{F}(p)=(\hat{p}+\hat{\mathbb{M}}^{\dagger})\frac{1}{p^{2}-\hat{\mathbb{M}}\hat{\mathbb{M}}^{\dagger}}$$
$${}=\frac{1}{p^{2}-\hat{\mathbb{M}}^{\dagger}\hat{\mathbb{M}}}(\hat{p}+\hat{\mathbb{M}}^{\dagger}).$$
(10)

The chirality-state projections can be simplified with the use of the relationships \(f(\hat{\mathbb{M}})\Pi_{\varsigma}=f(\mathbb{M}_{\varsigma})\Pi_{\varsigma}\) and \(f(\hat{\mathbb{M}}^{\dagger})\Pi_{\varsigma}=f(\mathbb{M}_{-\varsigma})\Pi_{\varsigma}\), which hold true for any function \(f(x)\) (\(\varsigma=\pm 1\)). The projection of \(\hat{\mathbb{S}}_{F}(p)\) involved in the creation of a neutrino with chirality \(-\varsigma\) has the form:

$$\hat{\mathbb{S}}_{F}(p)\Pi_{\varsigma}=(\hat{p}+\mathbb{M}_{-\varsigma})\frac{1}{p^{2}-\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}}\Pi_{\varsigma}.$$
(11)

The advantage of this representation is that in the denominator there is a \(n\times n\) matrix multiplied by a unit matrix in the space of bispinors, while the original expression contains a \(4n\times 4n\) matrix with the mixed bispinor and flavor indices. The mass spectrum is determined by poles of the propagator, in this instance by zeros of the characteristic polynomial

$$p(\lambda)=\det||\lambda-\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}||=0.$$
(12)

For three flavors, \(p(\lambda)\) is a polynomial of the third degree that has, respectively, either three real roots, or one real and two complex roots. The mass spectrum does not depend on \(\varsigma\), because \(\det||\lambda-\mathbb{M}_{L}\mathbb{M}_{R}||=\det||\lambda-\mathbb{M}_{R}\mathbb{M}_{L}||\). Given that \((\mathbb{M}_{L}\mathbb{M}_{R})^{*}=\mathbb{M}_{R}\mathbb{M}_{L}\), the complex conjugation of \(p(\lambda)\) gives the characteristic equation \(p(\lambda^{*})=0\). Thus, either all three roots are real, or one root is real, and two roots have imaginary parts of opposite sign. The product \(\mathbb{M}_{L}\mathbb{M}_{R}\), however, is a Hermitian matrix. Hence all eigenvalues of \(\mathbb{M}_{L}\mathbb{M}_{R}\) are real. Let \(\lambda_{i}=m_{i}^{2}\) be one of them. There exists a normalized vector \(|i\varsigma\rangle\) such that

$$\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}|i\varsigma\rangle=\lambda_{i}|i\varsigma\rangle.$$
(13)

Multiplying on the left by \(\langle i\varsigma|\) yields \(\lambda_{i}=\langle i\varsigma|(\mathbb{M}_{-\varsigma})^{\dagger}\mathbb{M}_{-\varsigma}|i\varsigma\rangle>0\). Furthermore, if \(|i\varsigma\rangle\) is a normalized eigenvector of \(\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}\), then

$$|i-\varsigma\rangle=\frac{1}{m_{i}}\mathbb{M}_{-\varsigma}|i\varsigma\rangle$$
(14)

is the normalized eigenvector of \(\mathbb{M}_{-\varsigma}\mathbb{M}_{\varsigma}\) with the same eigenvalue \(\lambda_{i}\). In the basis of flavor states \(|\alpha\rangle\), the coordinates of these eigenvectors appear to be complex conjugate: \(\langle\alpha|i\varsigma\rangle^{*}=\langle\alpha|i-\varsigma\rangle\). One can check the validity of this relation by writing Eq. (13) in the flavor basis. After performing a complex conjugation, an equation follows for the eigenvectors of \(\mathbb{M}_{-\varsigma}\mathbb{M}_{\varsigma}\), which means the eigenvectors of \(\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}\) and \(\mathbb{M}_{-\varsigma}\mathbb{M}_{\varsigma}\) differ by phase factors only. The phase factors of \(|i\varsigma\rangle\), however, can always be chosen so that the relationship

$$\langle\alpha|i\varsigma\rangle^{*}=\langle\alpha|i-\varsigma\rangle$$
(15)

is fulfilled. The mixing of Majorana neutrinos is described thereby in terms of a single mixing matrix. Equation (14) is equivalent to the condition \(\langle i-|\mathbb{M}_{L}|j+\rangle=m_{i}\delta_{ij}\), which shows that for \(n\) flavors the mass matrix \(\mathbb{M}_{L}\) is diagonalized by two unitary matrices:

$$\mathbb{V}^{\dagger}\mathbb{M}_{L}\mathbb{U}=\textrm{diag}(m_{1},\ldots,m_{n}),$$
(16)

with \(V_{\alpha i}=\langle\alpha|i-\rangle\) and \(U_{\alpha i}=\langle\alpha|i+\rangle\). Moreover, by virtue of Eq. (15), \(\mathbb{V}\) is expressed in terms of \(\mathbb{U}\): \(\mathbb{V}=\mathbb{U}^{*}\), so that

$$\mathbb{U}^{T}\mathbb{M}_{L}\mathbb{U}=\mathbb{V}^{T}\mathbb{M}_{R}\mathbb{V}=\textrm{diag}(m_{1},\ldots,m_{n}).$$
(17)

By means of the \(4n\times 4n\) matrix \(\hat{\mathbb{U}}=\mathbb{U}\Pi_{-}+\mathbb{V}\Pi_{+}\), the mass term of the Lagrangian transforms into a diagonal matrix,

$$\hat{\mathbb{M}}\rightarrow\gamma^{0}\hat{\mathbb{U}}^{\dagger}\gamma^{0}\hat{\mathbb{M}}\hat{\mathbb{U}}=\mathbb{V}^{\dagger}\mathbb{M}_{L}\mathbb{U}\Pi_{-}+\mathbb{U}^{\dagger}\mathbb{M}_{R}\mathbb{V}\Pi_{+}$$
$${}=\hat{1}\otimes\textrm{diag}\left(m_{1},...,m_{n}\right),$$
(18)

where \(\hat{1}\) is a unit matrix in the space of bispinors. Since \(\hat{\mathbb{U}}\) is real, the neutrality of bispinors in the flavor basis is a necessary and sufficient condition for the neutrality of bispinors diagonalizing the mass matrix. When the PMNS matrix \(\mathbb{U}\) is employed, Eq. (17) yields a mass matrix \(\mathbb{M}_{R}\) with non-zero phases on the diagonal. These phases are further eliminated by chiral rotations, as indicated below Eq. (6), after which the mass matrix takes on the form (7).

Using trace identities [12], the characteristic polynomial (12) for three neutrino flavors can be written as follows

$$p(\lambda)=(\lambda-\lambda_{0})^{3}$$
$${}-\frac{1}{2}(\lambda-\lambda_{0})\textrm{Tr}[\mathbb{C}^{2}]-\frac{1}{3}\textrm{Tr}[\mathbb{C}^{3}],$$
(19)

where

$$\mathbb{C}=\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}-\frac{1}{3}\textrm{Tr}[\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}]$$
(20)

is a traceless Hermitian matrix, and

$$\lambda_{0}=\frac{1}{3}\textrm{Tr}[\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}]$$
(21)

is a real positive number that represents the average of neutrino masses squared. The coefficients of the characteristic polynomial are traces of the powers of \(\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}\).

The matrix \(\mathbb{C}\) can be decomposed over the Gell–Mann–Okubo matrices: \(\mathbb{C}=\omega_{a}\lambda^{a}\), where \(\omega_{a}\) are real numbers (\(a=1,\ldots,8\)). Suppose we reduce \(\mathbb{C}\) to a diagonal form: \(\mathbb{C}=O\Lambda O^{\dagger}\), where \(O\) is unitary. The diagonal traceless matrix \(\Lambda\) has the form: \(\Lambda=a\lambda^{3}+b\lambda^{8}\). The eigenvalues of \(\mathbb{C}\) are \(a+b/\sqrt{3}\), \(-a+b/\sqrt{3}\), and \(-2b/\sqrt{3}\). For quantities \(a\) and \(b\), solving the cubic equation by the Viet method, we find

$$a=\sqrt{\frac{\textrm{Tr}[\mathbb{C}^{2}]}{2}}\sin\Theta,$$
(22)
$$b=-\sqrt{\frac{\textrm{Tr}[\mathbb{C}^{2}]}{2}}\cos\Theta,$$
(23)

where

$$\Theta=\frac{1}{3}\arccos\left(\frac{\sqrt{6}\textrm{Tr}[\mathbb{C}^{3}]}{(\textrm{Tr}[\mathbb{C}^{2}])^{3/2}}\right).$$
(24)

The masses of diagonal neutrinos are finally expressed in terms of \(\mathbb{M}^{\dagger}\mathbb{M}\):

$$m_{1}^{2}=\lambda_{0}+a+\frac{b}{\sqrt{3}},$$
$$m_{2}^{2}=\lambda_{0}-a+\frac{b}{\sqrt{3}},$$
$$m_{3}^{2}=\lambda_{0}-\frac{2b}{\sqrt{3}}.$$
(25)

The condition of \(m_{i}^{2}\) to be real also follows from the inequality \((\textrm{Tr}[\mathbb{C}^{2}])^{3}\geq 6(\textrm{Tr}[\mathbb{C}^{3}])^{2}\), which holds for all \(3\times 3\) traceless Hermitian matrices. The eigenvalues of \(\mathbb{M}_{R}\mathbb{M}_{L}\) coincide with the eigenvalues of \(\mathbb{M}_{L}\mathbb{M}_{R}\). In general, there are three mass eigenstates for three Majorana neutrinos. The ordering of \(m_{i}^{2}\) depends on the angle (24).

Using Sylvester’s matrix theorem [11], the propagator (10) can be written in the form:

$$\hat{\mathbb{S}}_{F}(p)=\sum_{i\varsigma}\frac{1}{p^{2}-m_{i}^{2}}\left(\hat{p}+\mathbb{M}_{-\varsigma}\right)|i\varsigma\rangle\langle i\varsigma|\Pi_{\varsigma}.$$
(26)

The neutrino states are projected onto the eigenstates of \({\mathbb{M}}_{\varsigma}{\mathbb{M}}_{-\varsigma}\) with the use of the Frobenius covariants:

$$|i\varsigma\rangle\langle i\varsigma|\equiv\mathbb{F}_{i\varsigma}=\prod_{j\neq i}\frac{\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}-\lambda_{j}}{\lambda_{i}-\lambda_{j}}.$$
(27)

The mass eigenstates with \(\varsigma=\pm 1\) form complete bases in flavor space:

$$\sum_{i}|i\varsigma\rangle\langle i\varsigma|=1.$$
(28)

To illustrate the statements made, consider

Example 1. Let the mass matrix be given by the equation

$$\hat{\mathbb{M}}=\lambda^{3}+i\gamma_{5}\lambda^{6},$$
(29)

where \(\lambda^{a}\) are the Gell-Mann matrices. The components of \(\hat{\mathbb{M}}\) are measured in mass units. The standard representation for \(\lambda^{a}\) is assumed to be identical with the flavor basis, so that

$$\mathbb{M}_{L}=\mathbb{M}_{R}^{\ast}=\lambda^{3}-i\lambda^{6}=\left|\left|\begin{matrix}1&0&0\\ 0&-1&-i\\ 0&-i&0\end{matrix}\right|\right|.$$
(30)

The eigenvalues and eigenvectors can readily be found using the MAPLE symbolic computation software package [18]:

$$(m_{Li})^{T}=(m_{Ri}^{*})^{T}=\left|\left|\begin{matrix}1\\ -1/\eta\\ -\eta\end{matrix}\right|\right|,$$
$$\tilde{W}_{\alpha i}=\ \left|\left|\begin{matrix}1&0&0\\ 0&-i/\eta&-i\eta\\ 0&1&1\end{matrix}\right|\right|,$$
(31)

where \(\eta=e^{i\pi/3}\), the index \(\alpha\) numbers the rows, the index \(i\) numbers the columns. After normalizing the eigenvectors of \(\mathbb{M}_{L}\):

$$\tilde{W}_{\alpha i}\rightarrow W_{\alpha i}=\left(\sum_{\beta}\tilde{W}_{\beta i}^{2}\right)^{-1/2}\tilde{W}_{\alpha i},$$
(32)

the matrix \(\mathbb{W}\) becomes a complex orthogonal matrix obeying

$$\sum_{\alpha}{W}_{\alpha i}{W}_{\alpha j}=\delta_{ij}$$

or, in matrix notation, \(\mathbb{W}^{T}\mathbb{W}=1\). The mass matrix is converted to the diagonal form in two steps. On the vector space of dimension \(4n\), we define first

$$\hat{\mathbb{W}}=\mathbb{W}\Pi_{-}+\mathbb{W}^{\ast}\Pi_{+}.$$
(33)

The action of \(\hat{\mathbb{W}}\) on the mass matrix gives

$$\hat{\mathbb{M}}\rightarrow\gamma^{0}\hat{\mathbb{W}}^{T}\gamma^{0}\hat{\mathbb{M}}\hat{\mathbb{W}}={\mathbb{W}}^{T}\mathbb{M}_{L}{\mathbb{W}}\Pi_{-}+{\mathbb{W}}^{{\dagger}}\mathbb{M}_{R}{\mathbb{W}^{*}}\Pi_{+}$$
$${}=\textrm{diag}\left(e^{i\gamma_{5}\vartheta_{1}}|m_{L1}|,...,e^{i\gamma_{5}\vartheta_{n}}|m_{Ln}|\right),$$
(34)

where \(e^{-2i\vartheta_{i}}=m_{Li}/m_{Ri}\). At the second stage, the chiral rotations \(\nu_{i}\rightarrow\exp(-i\gamma_{5}\vartheta_{i}/2)\nu_{i}\) eliminate the phase factors. \(\hat{{\mathbb{W}}}\) is a real matrix, therefore, it converts Majorana fields to Majorana fields. The drawback of \({\hat{\mathbb{W}}}\) is its non-unitarity.

The matrices included in the denominator of the neutrino propagator can be found to be

$$\mathbb{M}_{R}\mathbb{M}_{L}=\left(\mathbb{M}_{L}\mathbb{M}_{R}\right)^{\ast}=\left|\left|\begin{matrix}1&0&0\\ 0&2&i\\ 0&-i&1\end{matrix}\right|\right|.$$
(35)

The corresponding eigenvalues and eigenvectors are equal

$$(m_{i}^{2})^{T}=\left|\left|\begin{matrix}1\\ 2-g\\ 1+g\end{matrix}\right|\right|,$$
$$\tilde{U}_{\alpha i}=\tilde{V}_{\alpha i}^{*}=\left|\left|\begin{matrix}1&0&0\\ 0&-i/g&ig\\ 0&1&1\end{matrix}\right|\right|,$$
(36)

where \(g=(1+\sqrt{5})/2\). The normalization gives

$$\langle\alpha|i+\rangle=\langle\alpha|i-\rangle^{\ast}=\left(\sum_{\beta}|\tilde{U}_{\beta i}|^{2}\right)^{-1/2}\tilde{U}_{\alpha i}.$$

It can be verified that Eq. (14) applied to \(\langle\alpha|i+\rangle\) yields vectors proportional to \(\langle\alpha|i-\rangle\). Calculation of the parameters entering Eqs. (21)–(24) gives \(\lambda_{0}=4/3\), \(a=\frac{2}{\sqrt{3}}\cos(\frac{1}{3}\arccos(\frac{11}{16}))\), \(b=-\frac{2}{\sqrt{3}}\times\sin(\frac{1}{3}\arccos(\frac{11}{16}))\). Using relations (25) leads to values \(m_{i}^{2}\) that coincide with those given above.

Integrating the time-like component of the four-momentum in Eq. (8), we find neutrino propagator in the coordinate representation

$$i\hat{\mathbb{S}}_{F}(x-y)=\theta(x^{0}-y^{0})\sum_{i\varsigma}\int\frac{d\mathbf{p}}{(2\pi)^{3}2E_{i}(\mathbf{p})}$$
$${}\times e^{-iE_{i}(\mathbf{p})(x^{0}-y^{0})+i\mathbf{p}(\mathbf{x}-\mathbf{y})}\left(\hat{p}_{i}+\mathbb{M}_{-\varsigma}\right)\mathbb{F}_{i\varsigma}\Pi_{\varsigma}$$
$${}-\theta(y^{0}-x^{0})\sum_{i\varsigma}\int\frac{d\mathbf{p}}{(2\pi)^{3}2E_{i}(\mathbf{p})}$$
$${}\times e^{+iE_{i}(\mathbf{p})(x^{0}-y^{0})-i\mathbf{p}(\mathbf{x}-\mathbf{y})}\left(\hat{p}_{i}-\mathbb{M}_{-\varsigma}\right)\mathbb{F}_{i\varsigma}\Pi_{\varsigma},$$
(37)

where \(\hat{p}_{i}=E_{i}(\mathbf{p})\gamma_{0}-\boldsymbol{\gamma}\mathbf{p}\) and \(E_{i}(\mathbf{p})=\sqrt{m_{i}^{2}+\mathbf{p}^{2}}\).

Massive neutrinos change chirality during oscillations. The right antineutrinos take part in weak interaction. Vertices with right neutrinos should be taken into account because Majorana neutrinos are neutral. The matrix elements shown in the first row of Table 1 control the emission and absorption of neutrinos. The charged leptons are represented by \({\ell}_{\alpha}^{-}=(e^{-},\mu^{-},\tau^{-})\). The next rows of Table 1 contain three more pairs of matrix elements produced by permuting fermionic fields. The oscillations of the left neutrino created at point \(y\) into the left neutrino absorbed at point \(x\) are described by the matrix elements in the first row of Table 1. The oscillations of the right neutrino generated at point \(x\) into the right neutrino absorbed at point \(y\) occur with the reverse ordering of time \(x^{0}-y^{0}<0\).

Table 1 Four pairs of the lepton currents at points \(x\) and \(y\), leading to neutrino oscillations of various types, are presented in the first two columns; the next column shows the time ordering, the last column shows the processes described by the lepton currents

It is not difficult to identify the component of the neutrino propagator that triggers each of the processes shown in the table because projection operators \(\Pi_{\pm}\) enter the vertices. The component of \(\hat{\mathbb{S}}_{F}(x-y)\Pi_{+}\) proportionate to \(\hat{p}_{i}\Pi_{+}\) is suitable in the case of \(\nu_{\beta L}\leftarrow\nu_{\alpha L}\). When there are oscillations \(\nu_{\beta R}\rightarrow\nu_{\alpha L}\), the term proportional to \(\mathbb{M}_{L}\Pi_{+}\) is relevant, and so on. The vertices listed in the table only partially act independently. For instance, the transformation of the time inversion relates the process \(\nu_{\beta L}\leftarrow\nu_{\alpha L}\) to the process \(\nu_{\beta R}\rightarrow\nu_{\alpha R}\) in the second row: \(\mathcal{T}\hat{A}(\nu_{\beta L}\leftarrow\nu_{\alpha L})\mathcal{T}^{-1}=\hat{A}(\nu_{\beta R}\rightarrow\nu_{\textrm{R}\alpha})\). In relation to processes \(B^{c}\to A^{c}\), all processes \(A\to B\) share this feature.

The component of the amplitudes, \(\Pi_{\upsilon}i\hat{\mathbb{S}}_{F}(x-y)\Pi_{\varsigma}\) (\(\upsilon,\varsigma=\pm 1\)), that describes the neutrino propagation, is represented by a \(4\times 4\) matrix in the space of bispinors of the charged leptons. Neglecting the spread in the neutrino momenta, the neutrino component can be represented as follows

$$\hat{A}(\nu_{\beta L}\leftarrow\nu_{\alpha L})=-\hat{A}(\nu_{\beta R}\rightarrow\nu_{\alpha R})$$
$${}\propto-\sum_{i}e^{-iE_{i}(\mathbf{p})\tau}\langle\beta|i+\rangle\langle i+|\alpha\rangle\frac{\hat{p}_{i}}{E_{i}(\mathbf{p})}\Pi_{+}$$
$${}=-\sum_{i}e^{-iE_{i}(\mathbf{p})\tau}\langle\beta|\mathbb{F}_{i+}|\alpha\rangle\frac{\hat{p}_{i}}{E_{i}(\mathbf{p})}\Pi_{+},$$
(38)

where \(\tau=|x^{0}-y^{0}|\) and \(|\alpha\rangle\) are the neutrino flavor states. The matrix element

$$\langle\alpha|i+\rangle=U_{\alpha i}$$
(39)

can be recognized as the PMNS matrix element. The ratio \({\hat{p}_{i}}/{E_{i}(\mathbf{p})}\) is sensitive to the neutrino masses for \(|\mathbf{p}|\lesssim m_{i}\), the exponential factor \(\exp(-iE_{i}(\mathbf{p})\tau)\) is sensitive to the neutrino masses in a broad range of \(\mathbf{p}\) due to the macroscopic time factor \(\tau\), the Frobenius covariant (27) entering the second line depends on the neutrino masses from the start, while the other components of the amplitude have purely kinematic origin and are the same for all types of neutrinos.

In recent decades, the theory of neutrino oscillations in a vacuum and in a medium has been studied in many details [3, 7, 9, 19–23]. In the limit of large values of neutrino momenta, the kinematic factors of reactions involving neutrinos are factorized, and an expression for the probability of oscillations can be found:

$$P(\nu_{\beta L}\leftarrow\nu_{\alpha L})$$
$${}=\sum_{ij}e^{-i(E_{i}(\mathbf{p})-E_{j}(\mathbf{p}))t}U_{i\beta}^{\ast}U_{i\alpha}U_{j\beta}U_{j\alpha}^{\ast},$$
(40)

where \(t\) is the oscillation time. It is assumed that the momentum \(\mathbf{p}\) of neutrinos is fixed. This equation becomes

$$P(\nu_{\beta L}\leftarrow\nu_{\alpha L})$$
$${}=\sum_{ij}e^{-i(E_{i}(\mathbf{p})-E_{j}(\mathbf{p}))t}\langle\beta|\mathbb{F}_{i+}|\alpha\rangle\langle\alpha|\mathbb{F}_{j+}|\beta\rangle$$
(41)

once formulae (27) and (39) are taken into account. This equation is also a consequence of Eq. (38).

The probability has a similar form in the spatial picture with fixed energy \(E\). The exponential factor changes to

$$\exp(-i(E_{i}(\mathbf{p})-E_{j}(\mathbf{p}))t)$$
$${}\approx\exp\left(-i\frac{m_{i}^{2}-m_{j}^{2}}{2|\mathbf{p}|}t\right)$$
$${}\rightarrow\exp\left(-i\frac{m_{i}^{2}-m_{j}^{2}}{2E}L\right),$$
(42)

where \(L\) is the distance between the source and the detector. A description of the complete kinematic structure of the oscillation amplitudes in the spatial representation can be found in [21].

In this paper, the analysis of oscillations is based on the plane-wave approximation [22]. Despite the successful application of the plane-wave approximation, modeling processes using wave packets of the involved particles encompasses some challenges. Wave packets are employed in quantum field theory to accurately define asymptotic states, while scattering theory includes wave packets from the outset. When the problem is reformulated using the wave packet techniques, the standard formulas for the probability of neutrino oscillations are reproduced; however, neutrino oscillations are expected to be suppressed at distances greater than the coherence and dispersion lengths typical for neutrinos arriving from galactic sources [23].

Chirality-changing transitions of neutrinos are of interest in the context of double-\(\beta\) decay. The process \(\nu_{\beta R}\leftarrow\nu_{\alpha L}\) is described by the currents shown in the second raw of Table 1. The same matrix elements are responsible for the process \(\nu_{\beta R}\rightarrow\nu_{\alpha L}\). The component of the amplitude associated with the neutrino propagation equals \(\Pi_{+}i\hat{\mathbb{S}}_{F}(x-y)\Pi_{+}\), so we obtain

$$\hat{A}(\nu_{\beta R}\leftarrow\nu_{\alpha L})=-\hat{A}(\nu_{\beta R}\rightarrow\nu_{\alpha L})$$
$${}\propto-\sum_{i}e^{-iE_{i}(\mathbf{p})\tau}\ m_{i}\langle i+|\beta\rangle\langle i+|\alpha\rangle\frac{1}{E_{i}(\mathbf{p})}\Pi_{+}$$
$${}=-\sum_{i}e^{-iE_{i}(\mathbf{p})\tau}\langle\beta|\mathbb{M}_{L}\mathbb{F}_{i+}|\alpha\rangle\frac{1}{E_{i}(\mathbf{p})}\Pi_{+}.$$
(43)

Initially, under the sum sign, the expression \(\langle\beta|\mathbb{M}_{L}|i+\rangle\langle i+|\alpha\rangle\) occurs. In order to obtain the first line, we use the identities \(\langle\beta|\mathbb{M}_{L}|i+\rangle\) = \(m_{i}\langle\beta|i-\rangle\) = \(m_{i}\langle\beta|i+\rangle^{*}\) = \(m_{i}\langle i+|\beta\rangle\), which also allow calculating the mass matrix \(\langle\beta|\mathbb{M}_{L}|\alpha\rangle\) in terms of the neutrino masses and the mixing matrix:

$$\langle\beta|\mathbb{M}_{L}|\alpha\rangle=\sum_{i}m_{i}\langle i+|\beta\rangle\langle i+|\alpha\rangle.$$
(44)

This equation is equivalent to Eq. (17). The mass matrix of the right chirality is determined by \(\langle\beta|\mathbb{M}_{R}|\alpha\rangle=\langle\alpha|\mathbb{M}_{L}|\beta\rangle^{\ast}\). We remark out that, as can be seen from Example 1, the eigenvalues of the matrix \(\mathbb{M}_{L}\) do not match the masses \(m_{i}\) obtained by diagonalizing the matrix \(\mathbb{M}_{R}\mathbb{M}_{L}\). Section 4 gives further details of the double-\(\beta\) decay process.

Thus, the amplitude of processes associated with the emission and absorption of Majorana neutrinos is expressed explicitly in terms of the mass matrix using Eqs. (20)–(25) and (27).

2.2 Dirac Neutrino

The results can be generalized to the case of Dirac neutrinos. The Lagrangian for Dirac neutrinos has the form

$$\mathcal{L}_{D}=\bar{\nu}_{\alpha}i\hat{\nabla}\nu_{\alpha}-\bar{\nu}_{\alpha}\left(M_{\alpha\beta}+i\gamma_{5}N_{\alpha\beta}\right)\nu_{\beta},$$
(45)

where \(\nu_{\alpha}\) are the four-component complex spinors of Dirac neutrinos. Hermiticity of the Lagrangian implies that the matrices \(\mathbb{M}\) and \(\mathbb{N}\) are Hermitian, whereas \(\mathbb{M}_{\pm}=\mathbb{M}\pm i\mathbb{N}\) are arbitrary complex matrices satisfying \(\mathbb{M}_{+}^{\dagger}=\mathbb{M}_{-}\). The mass matrix is defined by \(\hat{\mathbb{M}}=\mathbb{M}+i\gamma_{5}\mathbb{N}\). The neutrino propagator has the form of Eq. (10). Because \(\mathbb{M}_{R}\mathbb{M}_{L}\) is Hermitian, Eqs. (8)–(28) apply to Dirac neutrinos, except for (15) and (17).

In the Weyl representation of two-component spinors \(\nu_{\alpha L}\) and \(\nu_{\alpha R}\), the Lagrangian (45) takes the form

$$\mathcal{L}_{M}=\bar{\nu}_{\alpha L}i\hat{\nabla}\bar{\nu}_{\alpha L}-\bar{\nu}_{\alpha L}M_{\alpha\beta}^{{D}}\nu_{\beta R}+\textrm{h.c.},$$
(46)

where \(M^{{D}}_{\alpha\beta}=M_{\alpha\beta}+iN_{\alpha\beta}\) is a complex \(n\times n\) matrix \(\mathbb{M}_{R}\).

The matrix \(\mathbb{M}_{R}\mathbb{M}_{L}\) for Majorana and Dirac neutrinos is controlled by \(n^{2}\) independent real parameters. Majorana neutrinos have \(n(n+1)\) and Dirac neutrinos have \(2n^{2}\) parameters in the matrix \(\mathbb{M}_{L}\) (or \(\mathbb{M}_{R}\)). The diagonalization of \(\mathbb{M}_{R}\mathbb{M}_{L}\) does not automatically imply that in the derived basis \(\mathbb{M}_{L}\) (or \(\mathbb{M}_{R}\)) is diagonal. The only result that follows from arguments used in the previous subsection is Eq. (18), where in the case of Dirac neutrinos \(\mathbb{U}\) and \(\mathbb{V}\) are two distinct unitary matrices. The unitary matrices \(\mathbb{V}\) and \(\mathbb{U}\) in Eq. (16) can be expressed as \(\mathbb{V}=\mathbb{P}_{1}\mathbb{V}^{\prime}\mathbb{P}_{2}\) and \(\mathbb{U}=\mathbb{P}_{3}\mathbb{U}^{\prime}\mathbb{P}_{4}\), respectively, where \(\mathbb{P}_{i}=\textrm{diag}(e^{i\chi_{i1}},\ldots,e^{i\chi_{in}})\). Without reducing generality, \(\chi_{21}=\chi_{41}=0\) can be used, resulting in the unitary matrices \(\mathbb{V}^{\prime}\) and \(\mathbb{U}^{\prime}\) being defined by \(n^{2}-2n+1\) independent parameters. Since the right neutrinos do not participate in interactions, it is possible to proceed to their linear combination by a unitary transformation. Rotating the right neutrinos with \(\mathbb{P}_{3}^{\dagger}\mathbb{V}\) removes \(n^{2}\) parameters. By phase rotations \(\nu_{L}\to\mathbb{P}_{4}^{\dagger}\nu_{L}\), the \(n-1\) additional parameters can be removed. The number of independent parameters of the mass matrix becomes \(n^{2}-n+1\). This result is the number of diagonal masses \(m_{i}\) plus the number of parameters of \(\mathbb{U}^{\prime}\). For \(n=3\) the mass matrix is determined by \(9-3+1=7\) parameters.

To illustrate the similarities and differences between the Dirac and Majorana mass matrices, consider

Example 2. Let the mass matrix be given by the equation

$$\hat{\mathbb{M}}=\lambda^{3}+i\gamma_{5}\lambda^{7}.$$
(47)

The components of \(\hat{\mathbb{M}}\) are measured in mass units. The matrix is not symmetric, but it satisfies the condition of Dirac self-conjugacy \(\gamma^{0}\hat{\mathbb{M}}^{{\dagger}}\gamma^{0}=\hat{\mathbb{M}}\). The left and righ chiral projections of \(\hat{\mathbb{M}}\) read

$$\mathbb{M}_{L}=\lambda^{3}-i\lambda^{7}=\left|\left|\begin{matrix}1&0&0\\ 0&-1&1\\ 0&-1&0\end{matrix}\right|\right|,$$
$$\mathbb{M}_{R}=\lambda^{3}+i\lambda^{7}=\left|\left|\begin{matrix}1&0&0\\ 0&-1&-1\\ 0&1&0\end{matrix}\right|\right|.$$
(48)

Unlike Majorana neutrinos, \(\mathbb{M}_{L}\) and \(\mathbb{M}_{R}\) are not complex conjugate. The eigenvalues and unnormalized eigenvectors can be found to be:

$$(m_{Li})^{T}=(m_{Ri})^{T}=\left|\left|\begin{matrix}1\\ -1/\eta\\ -\eta\end{matrix}\right|\right|,$$
$$\tilde{W}_{\alpha i}^{L}=\left|\left|\begin{matrix}1&0&0\\ 0&1/\eta&\eta\\ 0&1&1\end{matrix}\right|\right|,\quad\tilde{W}_{\alpha i}^{R}=\left|\left|\begin{matrix}1&0&0\\ 0&-1/\eta&-\eta\\ 0&1&1\end{matrix}\right|\right|,$$

where \(\eta=\) \(e^{i\pi/3}\). The matrices entering the denominators of the neutrino propagator are equal to

$$\mathbb{M}_{R}\mathbb{M}_{L}=\left|\left|\begin{matrix}1&0&0\\ 0&2&-1\\ 0&-1&1\end{matrix}\right|\right|,$$
$$\mathbb{M}_{L}\mathbb{M}_{R}=\left|\left|\begin{matrix}1&0&0\\ 0&2&1\\ 0&1&1\end{matrix}\right|\right|.$$
(49)

The corresponding eigenvalues and eigenvectors become

$$(m_{i}^{2})^{T}=\left|\left|\begin{matrix}1\\ 2-g\\ 1+g\end{matrix}\right|\right|,\quad\tilde{U}_{\alpha i}=\left|\left|\begin{matrix}1&0&0\\ 0&1/g&-g\\ 0&1&1\end{matrix}\right|\right|,$$
$$\tilde{V}_{\alpha i}=\left|\left|\begin{matrix}1&0&0\\ 0&-1/g&g\\ 0&1&1\end{matrix}\right|\right|,$$
(50)

where \(g=(1+\sqrt{5})/2\). Calculation of the parameters entering Eqs. (21)–(24) gives \(\lambda_{0}=4/3\), \(a=\frac{2}{\sqrt{3}}\sin(\frac{1}{3}\arccos(\frac{11}{16}))\), \(b=-\frac{2}{\sqrt{3}}\cos(\frac{1}{3}\arccos(\frac{11}{16}))\). Using the relations (25) leads to values \(m_{i}^{2}\) that coincide with those given above. A significant difference from the Majorana mass case is that two chiral projections require two mixing matrices,

$$\langle\alpha|i+\rangle=\left(\sum_{\beta}|\tilde{U}_{\beta i}|^{2}\right)^{-1/2}\tilde{U}_{\alpha i},$$
(51)
$$\langle\alpha|i-\rangle=\left(\sum_{\beta}|\tilde{V}_{\beta i}|^{2}\right)^{-1/2}\tilde{V}_{\alpha i},$$
(52)

to diagonalize the propagator. The relationship (15) for Dirac neutrinos is no longer valid, as can be seen from the example considered.

Table 2 shows the possible types of neutrino oscillations caused by the weak lepton currents at the space–time points \(x\) and \(y\). The combinations are also shown that could lead to the oscillations of Dirac neutrinos in the presence of right currents. One part of the amplitudes is expressed through the other by means of time inversion. The transition amplitudes represented in the last row are obtained by permuting the arguments \(x\leftrightarrow y\) of the amplitudes in the first row. The oscillation amplitudes corresponding to the first row are determined by Eqs. (38) with the right Majorana neutrinos \(\nu_{\alpha R}\) identified as the left Dirac antineutrinos \((\nu_{\alpha L})^{c}\).

Table 2 First two columns show the lepton currents at points \(x\) and \(y\), leading to neutrino oscillations; the next column shows the time ordering, the last column shows the process corresponding to the lepton currents; the second and third rows show the lepton currents that produce oscillations of Majorana neutrinos, but not Dirac neutrinos

The conclusion that is crucial to this discussion is that the neutrino mass matrix explicitly enters the amplitude of the oscillations of Dirac neutrinos in general.

Dirac neutrino is a superposition of two Majorana neutrinos with masses that are equal in absolute values and opposite in signs (see, e.g., [24]). Although the case described here is essentially a modified version of the one previously discussed, it is highlighted because of its significance in application. The mixed case of Majorana and Dirac neutrinos is reduced to Majorana neutrinos and is not discussed here.

3 NEUTRINO OSCILLATIONS IN MEDIUM

The electroweak mean-field potential created by the \(Z^{0}\) boson and the exchange potential of the \(W^{\pm}\) boson generated by charged leptons in the medium have both an effect on neutrinos. The mean field in the rest frame of matter has a time-like element: \((V_{\mu})_{\alpha\beta}=(V_{\alpha\beta},\mathbf{0)}\). The potential in the flavor basis is diagonal: \(V_{\alpha\beta}=V_{(\alpha)}\delta_{\alpha\beta}\), as illustrated in Fig. 1. Hereafter, \(\mathbb{V}\) stands for a matrix containing the components \(V_{\alpha\beta}\).

Fig. 1
figure 1

Diagrams contributing to the mean-field potential of neutrino in matter. \(f\) is a Standard Model fermion, \(\ell^{-}_{\alpha}\) is a charged lepton, and \(\nu_{\alpha}\) are neutrinos with the \(\alpha\) flavor. The \(Z\) boson exchange in the \(t\)-channel (\(a\)) produces a diagonal potential independent of neutrino flavor. Charged leptons and neutrinos interact additionally through the \(s\)-channel exchange of the \(W^{-}\) boson, as shown in the diagram (\(b\)). The individual components of the mean-field potential (\(b\)) are determined by the partial lepton densities. The total potential is diagonal in the neutrino flavor. The initial and final states of the medium particles are identical, since the mean-field potential is determined by the forward scattering amplitude at zero transmitted momentum.

3.1 Majorana Neutrino

The effective Lagrangian of Majorana neutrinos in the medium has the form

$$\mathcal{L}_{M}^{\prime}=\frac{1}{2}\bar{\nu}_{\alpha}\left(i\hat{\nabla}\delta_{\alpha\beta}+\hat{V}_{\alpha\beta}\gamma_{5}\right)\nu_{\beta}$$
$${}-\frac{1}{2}\bar{\nu}_{\alpha}\left(M_{\alpha\beta}+i\gamma_{5}N_{\alpha\beta}\right)\nu_{\beta}.$$
(53)

The mean-field potential \(\hat{V}_{\alpha\beta}=(\gamma_{0}V_{\alpha\beta},\mathbf{0})\) is determ̄ined by the forward scattering amplitude of neutrinos and depends on the composition of the substance. The initial and final states of the scattering particles are identical. Neutrino exchange with the medium particles via the \(Z\) boson is independent of the neutrino flavor. \(W^{\pm}\) bosons contribute to \(\mathbb{V}\) through the \(s\)-channel exchange involving lepton component of the substance. Figure 1 shows both types of diagrams, each generating a flavor-diagonal mean-field potential. Consequently, the mean-field potential is generally represented as \(V_{\alpha\beta}=V_{(\alpha)}\delta_{\alpha\beta}\). For example, there are no \(\tau^{-}\) mesons in interiors of neutron stars, therefore, the component of \(V_{(\alpha)}\) associated with the \(\tau^{-}\) meson is zero. The mean-field potential influences the left neutrino component. The mean field is a symmetric matrix in flavor indices, so its vector component disappears (cf. Eq. (3)) and only its axial component remains in Eq. (53).

In the momentum representation, the Feynman propagator takes the form

$$\hat{\mathbb{S}}_{F}(p)=\frac{1}{\hat{p}-\hat{\mathbb{M}}+\hat{\mathbb{V}}\gamma_{5}}=\left(\hat{p}+\hat{\mathbb{M}}^{\dagger}+\hat{\mathbb{V}}\gamma_{5}\right)$$
$${}\times\frac{1}{p^{2}-\hat{\mathbb{M}}\hat{\mathbb{M}}^{\dagger}-\hat{\mathbb{V}}^{2}+\left(\hat{p}\hat{\mathbb{V}}-\hat{\mathbb{V}}\hat{p}\right)\gamma_{5}}.$$
(54)

Unlike the case discussed in Section 2, the chiral projections do not eliminate \(\gamma\)-matrices completely from the denominator. Given \(\hat{p}\hat{\mathbb{V}}-\hat{\mathbb{V}}\hat{p}=2\mathbb{V}\boldsymbol{\alpha}\mathbf{p}\) in the rest frame of matter, it is helpful to decompose \(\hat{\mathbb{S}}_{F}(p)\) further using projections into states where the scalar product of the matrix \(\boldsymbol{\alpha}\) and the momentum \(\mathbf{p}\) equals \(\pm|\mathbf{p}|\). The corresponding projection operators are defined by \(P_{\pm}=(1\pm\boldsymbol{\alpha}\mathbf{n})/{2}\), where \(\mathbf{n}=\mathbf{p}/|\mathbf{p}|\), \(P_{\sigma}\) commute with \(\Pi_{\varsigma}\). In the standard representation and in the Weyl representation of the \(\gamma\)-matrices, the relationship \(\boldsymbol{\alpha}=\gamma_{5}\otimes\boldsymbol{\sigma}\) holds true, so that \(\boldsymbol{\alpha}\mathbf{n}\) can be interpreted as chirality times helicity (\(\sigma=\varsigma\eta\), \(\eta\) is helicity).

The neutrino propagator is represented by sum of the four projections:

$$\hat{\mathbb{S}}_{F}(p)=\sum_{\varsigma\sigma}\left(\hat{p}+{\mathbb{M}}_{-{\varsigma}}+{\varsigma}\hat{\mathbb{V}}\right)$$
$${}\times\frac{1}{p^{2}-{\mathbb{M}}_{{\varsigma}}{\mathbb{M}}_{-{\varsigma}}-{\mathbb{V}}^{2}+2{\varsigma}\sigma\mathbb{V}|\mathbf{p}|}\Pi_{\varsigma}P_{\sigma}.$$
(55)

The denominator is a Hermitian matrix, so the eigenvalues are real. They are determined from the characteristic equation

$$\det\left|\left|\lambda-{\mathbb{M}}_{{\varsigma}}{\mathbb{M}}_{-{\varsigma}}-{\mathbb{V}}^{2}+2\eta\mathbb{V}|\mathbf{p}|\right|\right|=0.$$
(56)

The value of \(\lambda\) depends on the momentum. The complex conjugation of the characteristic equation permutes the matrices \({\mathbb{M}}_{R}\) and \({\mathbb{M}}_{L}\), so the poles of the expression (55) depend on the helicity \(\eta\) only. The dispersion law takes the form

$$E_{i\eta}(\mathbf{p})=\sqrt{\lambda_{i\eta}+\mathbf{p}^{2}},$$
(57)

where \(\lambda_{i\eta}\) is one of the solutions. Entering the determinant (78) is a \(3\times 3\) matrix. This suggests that we are dealing with a cubic equation. To find the roots, one can use the Vieta method, which has already been employed to analyze neutrino oscillations in vacuum.

The procedure is equivalent to that used in Section 2. We redefine the matrix \(\mathbb{C}\) and the parameter \(\lambda_{0}\):

$$\mathbb{C}=\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+{\mathbb{V}}^{2}-2\eta\mathbb{V}|\mathbf{p}|$$
$${}-\frac{1}{3}\textrm{Tr}\left[\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+{\mathbb{V}}^{2}-2\eta\mathbb{V}|\mathbf{p}|\right],$$
(58)
$$\lambda_{0}=\frac{1}{3}\textrm{Tr}\left[\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+{\mathbb{V}}^{2}-2\eta\mathbb{V}|\mathbf{p}|\right].$$
(59)

The formulas (20)–(25) and (27) allow to find the eigenvalues \(\lambda_{i\eta}\) and the in-medium dispersion law.

The eigenstates determined by the equation

$$\left(\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+{\mathbb{V}}^{2}-2\eta\mathbb{V}|\mathbf{p}|\right)|i\varsigma\eta\rangle=\lambda_{i\eta}|i\varsigma\eta\rangle$$
(60)

depend on the momentum. The values of \(E_{i\eta}(\mathbf{p})\) are real, since

$$E_{i\eta}^{2}(\mathbf{p})$$
$${}=\langle{i\varsigma\eta}|\left((\mathbb{M}_{-\varsigma})^{\dagger}\mathbb{M}_{-\varsigma}+\left(|\mathbf{p}|-\eta\mathbb{V}\right)^{2}\right)|i\varsigma\eta\rangle\geq 0.$$

The eigenvalues of Eq. (60) are real and coincide for \(\varsigma=\pm 1\), while the coordinates of \(|i\pm\eta\rangle\) in the flavor basis are related by the equation

$$\langle\alpha|i\varsigma\eta\rangle^{*}=\langle\alpha|i-\varsigma\eta\rangle.$$
(61)

In the simplest case of the mean-field potential proportional to the unit matrix, which occurs when the diagram in Fig. 1\(b\)) is turned off, the momentum dependence of the dispersion law is quite transparent. Assuming \(V_{(\alpha)}=V\) for \(\alpha=1,2,3\), we find

$$E_{i\eta}(\mathbf{p})=\sqrt{m_{i}^{2}+\left(|\mathbf{p}|-\eta V\right)^{2}},$$
(62)

where \(m_{i}^{2}\) are the eigenvalues of \(\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}\). Since \(\mathbb{V}\) is commutative with \(\mathbb{M}_{\varsigma}\), the dependence on \(\eta\) disappears from \(|i\varsigma\eta\rangle\).

The equation for eigenvalues in an arbitrary potential admits a simple solution for neutrinos with a large momentum. In the left side of Eq. (60), the leading term is \(2\eta\mathbb{V}|\mathbf{p}|\) with \(\eta=\pm 1\). Since the potential is diagonal in flavor, the eigenstates are also diagonal in flavor. The first two terms on the left side of the equations give a correction to \(\lambda_{i\eta}\). In the leading order, the eigenvalues, the dispersion law, and the eigenstates take the form

$$\lambda_{\alpha\eta}=-2\eta V_{\alpha}|\mathbf{p}|+m_{\alpha}^{2}+{V}_{\alpha}^{2}+\ldots,$$
(63)
$$E_{\alpha\eta}=|\mathbf{p}|-\eta V_{\alpha}+\frac{m_{\alpha}^{2}}{2|\mathbf{p}|}+\ldots,$$
(64)
$$|\alpha\varsigma\eta\rangle=|\alpha\rangle+\frac{\eta}{2|\mathbf{p}|}$$
$${}\times\sum_{\beta\neq\alpha}|\beta\rangle\frac{1}{V_{(\beta)}-V_{(\alpha)}}\langle\beta|\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}|\alpha\rangle+\ldots,$$
(65)

where \(m_{\alpha}^{2}=\langle\alpha|\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}|\alpha\rangle\). In the limit of \(|\mathbf{p}|\to\infty\), the oscillations disappear for \(V_{(\beta)}\neq V_{(\alpha)}\). The degeneracy of the potential \(V_{(\beta)}\approx V_{(\alpha)}\) for any pair of \(\alpha\neq\beta\) restores the mixing [19, 20]. The mixing matrix (39) is split in the medium into two components depending on the parameter \(\eta\):

$$U^{\dagger}_{i\alpha}\to\langle i+|\alpha\varsigma\eta\rangle=U^{\dagger}_{i\alpha}+\frac{\eta}{2|\mathbf{p}|}$$
$${}\times\sum_{\beta\neq\alpha}U^{\dagger}_{i\beta}\frac{1}{V_{(\beta)}-V_{(\alpha)}}\langle\beta|\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}|\alpha\rangle+\ldots$$
(66)

The correction described by the second term can be significant for low values of \(|\mathbf{p}|\) or in the presence of degeneracy \(V_{(\beta)}\approx V_{(\alpha)}\).

Integrating the time-like component of four-momentum in the coordinate representation, we obtain

$$i\hat{\mathbb{S}}_{F}(x-y)=\theta(x^{0}-y^{0})$$
$${}\times\sum_{i\varsigma\eta}\int\frac{d\mathbf{p}}{(2\pi)^{3}2E_{i\eta}(\mathbf{p})}e^{-iE_{i\eta}(\mathbf{p})(x^{0}-y^{0})+i\mathbf{p}(\mathbf{x}-\mathbf{y})}$$
$${}\times\left(\hat{p}_{i\eta}+{\mathbb{M}}_{-\varsigma}+\varsigma\hat{\mathbb{V}}\right)|i\varsigma\eta\rangle\langle i\varsigma\eta|\Pi_{\varsigma}P_{\sigma}$$
$${}-\theta(y^{0}-x^{0})\sum_{i\varsigma\eta}\int\frac{d\mathbf{p}}{(2\pi)^{3}2E_{i\eta}(\mathbf{p})}$$
$${}\times e^{+iE_{i\eta}(\mathbf{p})(x^{0}-y^{0})-i\mathbf{p}(\mathbf{x}-\mathbf{y})}$$
$${}\times\left(\hat{p}_{i\eta}-{\mathbb{M}}_{-\varsigma}-\varsigma\hat{\mathbb{V}}\right)|i\varsigma\eta\rangle\langle i\varsigma\eta|\Pi_{\varsigma}P_{\sigma},$$
(67)

where \(\sigma=\varsigma\eta\).

The Frobenius covariants are given by

$$|i\varsigma\eta\rangle\langle i\varsigma\eta|\equiv\mathbb{F}_{i\varsigma\eta}$$
$${}=\prod_{j\neq i}\frac{\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+{\mathbb{V}}^{2}-2\eta\mathbb{V}|\mathbf{p}|-\lambda_{j\eta}}{\lambda_{i\eta}-\lambda_{j\eta}}.$$
(68)

Like in the previous section we consider first neutrino oscillations \(\nu_{\beta L}\leftarrow\nu_{\alpha L}\), described by the current \(j_{\mu}^{\textrm{CC}}=\bar{\nu}_{\alpha}\gamma_{\mu}(1-\gamma_{5})\ell_{\alpha}^{-}\) at the neutrino creation point \(y\) and by the current \(j_{\mu}^{\textrm{CC}}=\bar{\ell}_{\beta}^{-}\gamma_{\mu}(1-\gamma_{5})\nu_{\beta}\) at the neutrino absorption point \(x\). A process is also possible in which the right neutrino \(\nu_{\alpha R}\) is created at the point \(x\) and propagates to the point \(y\) without chirality flip, where it is absorbed (\(x^{0}-y^{0}<0\)). Neglecting the neutrino momentum spread, the neutrino component of the amplitude, \(\Pi_{-}i\hat{\mathbb{S}}_{F}(x-y)\Pi_{+}\), becomes

$$\hat{A}(\nu_{\beta L}\leftarrow\nu_{\alpha L})\propto-\sum_{i\eta}e^{-iE_{i\eta}(\mathbf{p})\tau}\langle\beta|\left(\hat{p}_{i\eta}+\hat{\mathbb{V}}\right)$$
$${}\times\mathbb{F}_{i\varsigma\eta}|\alpha\rangle\frac{1}{E_{i\eta}(\mathbf{p})}\Pi_{\varsigma}P_{\eta},$$
(69)
$$\hat{A}(\nu_{\beta R}\rightarrow\nu_{\alpha R})\propto\sum_{i\eta}e^{-iE_{i\eta}(\mathbf{p})\tau}\langle\beta|\left(\hat{p}_{i\eta}-\hat{\mathbb{V}}\right)$$
$${}\times\mathbb{F}_{i\varsigma\eta}|\alpha\rangle\frac{1}{E_{i\eta}(\mathbf{p})}\Pi_{\varsigma}P_{\eta},$$
(70)

where \(\varsigma=+1\).

The mixing matrices with chirality \(\varsigma=\pm 1\) are connected by Eq. (61), however to characterize the mixing of neutrinos with helicity \(\eta=\pm 1\), two distinct matrices \(\langle\alpha|i++\rangle\) and \(\langle\alpha|i+-\rangle\) are required. In the limit of vanishing matter density, the matrices \(\langle\alpha|i+\pm\rangle\) become the PMNS matrix.

The process \(\nu_{\beta R}\leftarrow\nu_{\alpha L}\) with chirality flip is associated with the lepton current \(j_{\mu}^{\textrm{CC}}=\bar{\nu}_{\alpha}\gamma_{\mu}(1-\gamma_{5})\ell_{\alpha}^{-}\) at the creation point \(y\) and by the lepton current \(j_{\mu}^{\textrm{CC}}=\bar{\nu}_{\beta}\gamma_{\mu}(1-\gamma_{5})\ell_{\beta}^{-}\) at the absorption point \(x\). The process \(\nu_{\beta R}\rightarrow\nu_{\alpha L}\) is described for \(x^{0}-y^{0}>0\) by the same currents. The component of the amplitudes responsible for chirality-flip propagation of neutrinos equals \(\Pi_{+}i\hat{\mathbb{S}}_{F}(x-y)\Pi_{+}\), so we obtain

$$\hat{A}(\nu_{\beta R}\leftarrow\nu_{\alpha L})=-\hat{A}(\nu_{\beta R}\rightarrow\nu_{\alpha L})$$
$${}\propto-\sum_{i\eta}e^{-iE_{i\eta}(\mathbf{p})\tau}\langle\beta|\mathbb{M}_{L}\mathbb{F}_{i\varsigma\eta}|\alpha\rangle$$
$${}\times e^{-i\phi_{\alpha}-i\phi_{\beta}}\frac{1}{E_{i\eta}(\mathbf{p})}\Pi_{\varsigma}P_{\eta},$$
(71)

where \(\varsigma=+1\). The explicit dependence on the neutrino mass matrix is obvious. The arguments given after Eq. (43) prove the equation \(\langle\beta|\mathbb{M}_{L}|i\varsigma\eta\rangle=m_{i}\langle i\varsigma\eta|\beta\rangle\) for \(\mathbb{V}\) commutative with \(\mathbb{M}_{\pm}\).

3.2 Dirac Neutrino

The in-medium Lagrangian for Dirac neutrinos is written as

$$\mathcal{L}_{D}^{\prime}=\bar{\nu}_{\alpha}\left(i\hat{\nabla}\delta_{\alpha\beta}-\hat{V}_{\alpha\beta}\Pi_{-}\right)\nu_{\beta}$$
$${}-\bar{\nu}_{\alpha}\left(M_{\alpha\beta}+i\gamma_{5}N_{\alpha\beta}\right)\nu_{\beta}.$$
(72)

The mean-field term \(\hat{V}_{\alpha\beta}\Pi_{-}\) distinguishes it from the vacuum Lagrangian, \(\mathcal{L}_{D}\). The mass matrices \(\mathbb{M}\) and \(\mathbb{N}\) are Hermitian. The Feynman propagator has the following form in the momentum representation:

$$\hat{\mathbb{S}}_{F}(p)=\frac{1}{\hat{p}-\hat{\mathbb{M}}-\hat{\mathbb{V}}\Pi_{-}}.$$
(73)

Example 3. Given that \(\hat{\mathbb{V}}\) is determined by time-like component of the four-vector, the substitutions \(E\to E-\mathbb{V}\) and \(\mathbb{V}\to\mathbb{V}/2\) reduce the problem to the case of Majorana neutrinos provided \(\mathbb{V}\) is proportional the unit matrix. In such a case, the dispersion law reads

$${E}_{i\eta}^{D\pm}(\mathbf{p})=\frac{1}{2}V\pm E_{i\eta}(\mathbf{p}),$$
(74)

where \(E_{i\eta}(\mathbf{p})\) simplifies to (62) with \(V\) replaced by \(V/2\). In a neutron star, the mean-field potential generated by neutrons is \(V\sim-20\) eV (see, e.g., [25]). At zero temperature and for \(V/2+m_{i}<0\) the potential well is filled first by Dirac neutrinos with a negative energy at \(|\mathbf{p}|\sim|V|/2\), then states with higher energy are filled. At a low density of neutrinos, the design of the Fermi sphere turns out to be quite unusual: the neutrino states form in the momentum space a hollow sphere with a radius close to \(|V|/2\). Figure 2 depicts the dispersion law for the Majorana neutrino for \(V=-10\) eV. The positive-energy Dirac neutrino dispersion law is similar, but it is displaced vertically and has a modified scale of the horizontal axis.

Fig. 2
figure 2

Two branches of dispersion law for Majorana neutrinos in a medium with a mean-field potential of \(V=-10\) eV, typical for nuclear matter at saturation, and a neutrino mass of \(m_{i}=0.1\) eV. The upper and lower lines correspond to neutrinos with left and right helicity \(\eta=\mp 1\), respectively. The neutrino energy takes a minimum value \(E_{i+}(\mathbf{p})=m_{i}\) at \(|\mathbf{p}|=10\) eV. The dispersion law of positive-energy Dirac neutrinos takes the same form for \(V=-20\) eV and is shifted along the vertical axis accordingly.

The BCS mechanism causes a gap to appear in the neutrino excitation spectrum at ultra-low temperatures [26]. The pairing of Dirac neutrinos is caused by the exchange of the Higgs boson, which generates an attraction between fermions. In the case under consideration, the BCS mechanism is applicable to neutrinos adjacent to both the outer surface of the hollow Fermi sphere and the inner. At the temperature of old neutron stars of the order of 100 eV, condensation has not yet formed, so the appearance of superfluidity can be expected only in the very distant future as a result of further cooling of neutron stars.

The bottom curve of the dispersion law reaches the border of the continuum of states with negative energy at a depth of the potential well of \(-4m_{i}\). This leads to the creation of neutrino-antineutrino pairs from vacuum and filling of levels with energy \(E_{i}^{D+}(\mathbf{p})<-m_{i}\).

The technique for determining the dispersion law in a general scenario for in-medium Dirac neutrinos is more involved. We decompose first the neutrino propagator (73) in a series by powers of the potential. The chirality-state projections of each of the sum terms can be simplified by using Eq. (10). The summation of the outcome results in

$$\hat{\mathbb{S}}_{F}(p)=\left(\hat{p}+\mathbb{M}_{L}\right)\frac{1}{p^{2}-\mathbb{M}_{R}\mathbb{M}_{L}-\hat{\mathbb{V}}\hat{p}}\Pi_{+}$$
$${}+\left(\hat{p}+\mathbb{M}_{R}\right)\frac{1}{p^{2}-\mathbb{M}_{L}\mathbb{M}_{R}}\Pi_{-}+\left(\hat{p}+\mathbb{M}_{L}\right)$$
$${}\times\frac{1}{p^{2}-\mathbb{M}_{R}\mathbb{M}_{L}-\hat{\mathbb{V}}\hat{p}}\hat{\mathbb{V}}\mathbb{M}_{R}\frac{1}{p^{2}-\mathbb{M}_{L}\mathbb{M}_{R}}\Pi_{-}.$$
(75)

The second and third lines, which refer to \(\hat{\mathbb{S}}_{F}(p)\Pi_{-}\), have two types of poles. The first corresponds to a neutrino in a medium, while the second to a neutrino in a vacuum. The second type pole is inactive. To demonstrate this, we employ the chirality-state projectors to act on \(\hat{\mathbb{S}}_{F}(p)\Pi_{-}\) from the left. Basic algebra gives

$$\hat{\mathbb{S}}_{F}(p)\Pi_{-}$$
$${}=\mathbb{M}_{L}\frac{1}{p^{2}-\mathbb{M}_{R}\mathbb{M}_{L}-\hat{\mathbb{V}}\hat{p}}\left(\hat{p}-\hat{\mathbb{V}}\right)\mathbb{M}_{L}^{-1}\Pi_{-}$$
$${}+\frac{1}{p^{2}-\mathbb{M}_{R}\mathbb{M}_{L}-\hat{p}\hat{\mathbb{V}}}\mathbb{M}_{R}\Pi_{-}.$$
(76)

The transformations are carried out under the assumptions that the matrix \(\mathbb{V}\) only has a time-like component, all matrices are non-singular, and the matrix \(\mathbb{V}\) generally does not commute with the mass matrix. Since the poles of the free neutrino propagator vanish, the vacuum poles of the expression (75) are fictitious. With the substitutions \(E\rightarrow E-\mathbb{V}\) and \(\mathbb{V}\to\mathbb{V}/2\), the propagator poles for Majorana neutrino and Dirac neutrino formally coincide.

The denominators of propagators turn into \(n\times n\) matrices in the flavor space using projectors \(P_{\sigma}\). Using the commutativity of \(P_{\sigma}\) and \(\Pi_{\varsigma}\), the relationships \(\hat{\mathbb{V}}\hat{p}P_{\sigma}=\mathbb{V}(E-\sigma|\mathbf{p}|)P_{\sigma}\), \(P_{\sigma}\gamma_{0}=\gamma_{0}P_{-\sigma}\) and \((\hat{p}-\hat{\mathbb{V}})P_{\sigma}=P_{-\sigma}(\hat{p}-\hat{\mathbb{V}})\), the denominators of \(\hat{\mathbb{S}}_{F}(p)\Pi_{\varsigma}P_{\sigma}\) can be transformed into \(p^{2}-\mathbb{M}_{R}\mathbb{M}_{L}-\mathbb{V}(E-\eta|\mathbf{p}|)\). The excitation spectrum is determined thereby from the equation:

$$\left(\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+\mathbb{V}(E-\eta|\mathbf{p}|)\right)|i\sigma\varsigma\rangle=\lambda_{i\eta}(E)|i\varsigma\eta\rangle.$$
(77)

It is useful to compare this equation to Eq. (60).

The characteristic equation of Dirac neutrinos differs from that of Majorana neutrinos in terms depending on the potential:

$$\det\left|\left|\lambda-{\mathbb{M}}_{{\varsigma}}{\mathbb{M}}_{-\varsigma}-\mathbb{V}(E-\eta|\mathbf{p}|)\right|\right|=0.$$
(78)

For completeness, we also give the expressions of the matrix \(\mathbb{C}\) and the parameter \(\lambda_{0}\):

$$\mathbb{C}=\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+\mathbb{V}(E-\eta|\mathbf{p}|)$$
$${}-\frac{1}{3}\textrm{Tr}\left[\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+\mathbb{V}(E-\eta|\mathbf{p}|)\right],$$
(79)
$$\lambda_{0}=\frac{1}{3}\textrm{Tr}\left[\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+\mathbb{V}(E-\eta|\mathbf{p}|)\right].$$
(80)

The relations (22)–(24) with \(n=3\) are fully applicable for Dirac neutrinos in matter, in the relations (25) it is necessary to replace \(m_{i}^{2}\to\lambda_{i\eta}(E)\). Unlike the Majorana neutrino, \(\lambda_{i\eta}(E)\) depends on the energy. The characteristic equation, therefore, is a polynomial of degree \(6\) with respect to the energy, and for a matrix \(\mathbb{V}\) of general form it is not solved in quadratures.

The Frobenius covariants read:

$$|i\varsigma\eta\rangle\langle i\varsigma\eta|\equiv\mathbb{F}_{i\varsigma\eta}(E,\mathbf{p})$$
$${}=\prod_{j\neq i}\frac{\mathbb{M}_{\varsigma}\mathbb{M}_{-\varsigma}+\mathbb{V}(E-\eta|\mathbf{p}|)-\lambda_{j\eta}(E)}{\lambda_{i\eta}(E)-\lambda_{j\eta}(E)}.$$
(81)

Let us consider the case of large momenta. The dominant term \(\mathbb{V}(-\eta|\mathbf{p}|)\) in Eq. (77) if not entirely compensated by the term \(\mathbb{V}E\) implies that the zeroth approximation is formed by the flavor states. In the leading order, we obtain

$$\lambda_{\alpha\eta}(E_{\alpha\eta}^{D+})=(1-\eta)V_{\alpha}|\mathbf{p}|+m_{\alpha}^{2}$$
$${}+\frac{1}{2}(1-\eta)V_{\alpha}^{2}+\ldots,$$
(82)
$$E_{\alpha\eta}^{D+}=|\mathbf{p}|+\frac{1}{2}(1-\eta)V_{\alpha}+\frac{m_{\alpha}^{2}}{2|\mathbf{p}|}+\ldots,$$
(83)
$$|\alpha\varsigma\eta\rangle=|\alpha\rangle$$
$${}-\frac{1}{2|\mathbf{p}|}\sum_{\beta\neq\alpha}|\beta\rangle\frac{1}{V_{(\beta)}-V_{(\alpha)}}\langle\beta|\mathbb{M}_{R}\mathbb{M}_{L}|\alpha\rangle+\ldots$$
$$(\eta=-1),$$
(84)

where \(m_{\alpha}^{2}=\langle\alpha|\mathbb{M}_{R}\mathbb{M}_{L}|\alpha\rangle.\) As noted earlier, the poles of the Majorana neutrino and Dirac neutrino propagators formally coincide after the energy shift \(E\rightarrow E-\mathbb{V}\) and redefinition \(\mathbb{V}\to\mathbb{V}/2\); this is reflected in the dispersion laws (64) and (83). According to Eq. (83), a Dirac neutrino with right helicity does not feel the external potential, which physically distinguishes a Dirac neutrino from a Majorana neutrino. Moreover, for \(\eta=+1\), the term \(\mathbb{V}(E-\eta|\mathbf{p}|)\) in Eq. (77) is not dominant, so the validity of the expressions (82)–(84) is restricted to the case \(\eta=-1\).

We omit the intermediate steps and give the final expression for neutrino component of the oscillation amplitudes of Dirac neutrinos

$$\hat{A}(\nu_{\beta L}\leftarrow\nu_{\alpha L})=-\hat{A}((\nu^{c})_{\beta R}\rightarrow(\nu^{c})_{\alpha R})$$
$${}\propto-\sum_{i\eta}e^{-iE_{i\eta}^{D+}(\mathbf{p})\tau}\langle\beta|\mathbb{F}_{i\varsigma\eta}(E_{i\eta}^{D+},\mathbf{p})|\alpha\rangle$$
$${}\times\frac{2\hat{p}_{i\eta}^{D+}}{2{E}_{i\eta}^{D+}(\mathbf{p})-\lambda_{i\eta}^{\prime}({E}_{i\eta}^{D+}(\mathbf{p}))}\Pi_{\varsigma}P_{\eta},$$
(85)

where \(\varsigma=+1\), \(\lambda_{i\eta}^{\prime}(E)={\partial}\lambda_{i\eta}(E)/{\partial E}\), and \(E\) is replaced by \(E_{i\eta}^{D+}\) in \(\mathbb{F}_{i\varsigma\eta}(E,\mathbf{p})\). By virtue of Eq. (77) \(\lambda_{i\eta}^{\prime}(E)=\langle i\varsigma\eta|\mathbb{V}|i\varsigma\eta\rangle\). Expressions in quadratures complicate the analysis of the analytical properties of the propagator as a function of energy. Returning to expressions (75) and (76), we note that the inversion of matrices in the denominators results in the denominators becoming a polynomial of degree \(2n=6\), which, in the absence of degeneracy, has simple zeros on the real axis of energy. The expression shown in (85) is the result of integration with respect to the time-like variable \(p^{0}\), while the integration contour is closed in the lower half-plane.

4 OTHER NEUTRINO-RELATED PROCESSES

It was demonstrated in the preceding sections that experiments involving neutrino oscillations in a medium and in a vacuum can both be explained solely in terms of the neutrino mass matrix without the need for the mixing matrix. This is demonstrated via Sylvester’s theorem, which makes use of Frobenius covariants. Additional neutrino-related processes are covered below. It turns out that the probabilities are also directly expressed using the neutrino mass matrix in all of the scenarios that are considered to be the most typical.

4.1 \(\beta\)-Decay of Tritium

Tritium is used as a \(\beta\) emitter in advanced experiments on the direct measurement of the mass of an electron antineutrino [27–29]. For the last decade, precise measurements of the electron antineutrino mass have been carried out by the KATRIN collaboration. As a result of the analysis of the edge of the experimental spectrum, a restriction on the neutrino mass \({<}0.8\) eV at \(\textrm{CL}=90\%\) was obtained [29].

The edge of the \(\beta\) electron energy spectrum is not sensitive to decay channels associated with the electron shell excitation. However, when the \(\beta\) electron energy exceeds the binding energy of an electron in a tritium atom, it is necessary to take into account new channels available in energy. The loss of energy for ionization leads to a distortion of the shape of the \(\beta\) spectrum. In the \(\beta\) decay of tritium, the Coulomb interaction of the knocked-out electron with the nucleus and the channels associated with the population of higher discrete levels by an atomic electron are also important to construct a realistic \(\beta\) decay spectrum over the entire range of the energy spectrum of \(\beta\) electrons. Theoretical corrections arising on the atomic and nuclear level are also discussed in [30, 31].

Neglecting electron shell structure effects the \(\beta\)-decay probability of tritium gives (see, e.g., [32])

$$d\Gamma_{\beta}=\frac{G_{\textrm{F}}^{2}\cos^{2}\Theta_{\textrm{C}}}{2\pi^{3}}|\mathcal{M}_{\textrm{nucl}}|^{2}F(E,Z=2)$$
$${}\times\sum_{i}\left|\langle i+|\alpha\rangle\right|^{2}p_{\nu_{i}}E_{\nu_{i}}p_{e}E_{e}dE_{e},$$
(86)

where \(G_{\textrm{F}}\) denotes the Fermi constant, \(\Theta_{\textrm{C}}\) is the Cabibbo angle, \(\mathcal{M}_{\textrm{nucl}}\) is the nuclear matrix element, \(F(E,Z=2)\) is the Fermi function, \(p_{\nu_{i}}\) and \(E_{\nu_{i}}\) are the momentum and energy of the neutrino with a mass of \(m_{i}\), \(p_{e}\) and \(E_{e}\) are the momentum and energy of the \(\beta\) electron, \(Q=E_{\nu_{i}}+E_{e}\) is the decay energy, \(|\alpha\rangle\) is the electron-flavor neutrino state, and \(\langle i+|\alpha\rangle\) is the mixing matrix. The kinematic factor is the phase space volume. Using the formula (27), we obtain

$$d\Gamma_{\beta}=\frac{G_{\textrm{F}}^{2}\cos^{2}\Theta_{\textrm{C}}}{2\pi^{3}}|\mathcal{M}_{\textrm{nucl}}|^{2}F(E,Z=2)$$
$${}\times\sum_{i}\langle\alpha|\mathbb{F}_{i+}|\alpha\rangle p_{\nu_{i}}E_{\nu_{i}}p_{e}E_{e}dE_{e}.$$
(87)

The decay probability is expressed in terms of the neutrino mass matrix in the flavor basis.

4.2 Electron Capture

The inverse to the \(\beta\) decay of an atom is the capture of an orbital electron. The following formula for the total probability is used in the limit of zero width of excited electron shells [33]:

$$\lambda_{x}=\frac{G_{\textrm{F}}^{2}\cos^{2}\Theta_{\textrm{C}}}{(2\pi)^{2}}n_{x}\mathcal{B}_{x}\beta_{x}^{2}C_{x}$$
$${}\times\sum_{i}\left|\langle i+|\alpha\rangle\right|^{2}p_{\nu_{i}}E_{\nu_{i}},$$
(88)

where \(n_{x}\) is the relative occupation number of the \(x\) shell, \(\beta_{x}\) is proportional to the wave function of the captured electron inside the parent nucleus. The overlap-exchange correction \(\mathcal{B}_{x}\) is a square of the atomic transition matrix element which takes into account collective electron effects. The phase-space volume is proportional to \(p_{\nu_{i}}E_{\nu_{i}}\). The shape function of the transition, \(C_{x}\), contains information about the nuclear structure and transformation properties of the decay amplitude. Other notations as in the previous subsection.

In terms of the Frobenius covariant (27), the decay probability becomes

$$\lambda_{x}=\frac{G_{\textrm{F}}^{2}\cos^{2}\Theta_{\textrm{C}}}{(2\pi)^{2}}n_{x}\mathcal{B}_{x}\beta_{x}^{2}C_{x}$$
$${}\times\sum_{i}\langle\alpha|\mathbb{F}_{i+}|\alpha\rangle p_{\nu_{i}}E_{\nu_{i}}.$$
(89)

Electron capture is used or is planned to be used to measure the electron neutrino mass [34–38].

4.3 Neutrinoless Double-\(\beta\) Decay

Experimental observation of neutrinoless double-\(\beta\) decays (\(0\nu\beta^{-}\beta^{-}\), \(0\nu\beta^{+}\beta^{+}\)), neutrinoless double electron capture (\(0\nu\)2EC) or mixed decays \(0\nu\beta^{-}\)EC could prove the non-conservation of the total lepton number and, accordingly, Majorana nature of neutrino mass. These processes are also of interest for determining the absolute scale and the hierarchy type of neutrino masses. In view of the value of processes sensitive to the non-conservation of the total lepton number, extensive literature is devoted to the physics of \(0\nu\beta\beta\) decay, \(0\nu\)2EC, and the corresponding models of the nuclear structure [39–45].

The \(0\nu\beta^{-}\beta^{-}\), \(0\nu\beta^{+}\beta^{+}\), \(0\nu\beta^{+}\)EC decays and \(0\nu\)2EC are caused by the chirality flip of a massive Majorana neutrino. The amplitudes of the \(0\nu\beta^{+}\beta^{+}\), \(0\nu\beta^{+}\)EC and \(0\nu\)2EC processes are proportional to the propagator projected onto the left state of the emitted and the left state of the absorbed neutrino. The projection of the propagator can be written as follows

$$\langle\alpha|\Pi_{+}\hat{\mathbb{S}}_{F}(p)\Pi_{+}|\alpha\rangle$$
$${}=\sum_{i}\frac{1}{p^{2}-m_{i}^{2}}\langle\alpha|\mathbb{M}_{L}\mathbb{F}_{i+}|\alpha\rangle\Pi_{+},$$
(90)

where \(|\alpha\rangle\) is the electron-flavor state of the neutrino in the initial and final states.Footnote 1 We use the condition of completeness of mass eigenstates \(|i\varsigma\rangle\) and take into account that the transmitted momentum of neutrino has the order of the Fermi momentum at saturation \(p_{\textrm{F}}=270\) MeV. Since this momentum is much higher than the neutrino masses \(m_{i}\lesssim 1\) eV, one can replace \(m_{i}^{2}\) by the average value of the squares of the neutrino masses, \(\langle m_{i}^{2}\rangle\). By performing the summation, one has

$$\langle\alpha|\Pi_{+}\hat{\mathbb{S}}_{F}(p)\Pi_{+}|\alpha\rangle$$
$${}\approx\frac{1}{p^{2}-\langle m_{i}^{2}\rangle}\langle\alpha|\mathbb{M}_{L}|\alpha\rangle\Pi_{+}.$$
(91)

The decay rate is determined by a single element of the mass matrix:

$$d\Gamma\varpropto\left|\langle\alpha|\mathbb{M}_{L}|\alpha\rangle\right|^{2}=\left|\sum_{i}m_{i}\langle i+|\alpha\rangle^{2}\right|^{2}.$$
(92)

Contributions to the double-\(\beta\) processes of Majorana neutrinos with a mass greater than \(p_{\textrm{F}}\) are also discussed [46]. For heavy neutrinos, it is necessary to perform summation in Eq. (90) with the exact Frobenius covariants (27).

Equation (90) shows that the \(0\nu\beta\beta\), \(0\nu\beta^{-}\)EC and \(0\nu\)2EC probabilities are expressed in terms of the mass matrix alone.

4.4 Quasi-Elastic Neutrino Scattering

As the next example, consider a neutrino scattered on a lepton: \(\nu_{i}+\ell_{i}^{-}\to\nu_{j}+\ell_{f}^{-}\). According to the scattering theory, there is an asymptotic state \(|i+\rangle\otimes|\ell_{i}^{-}\rangle\) of a neutrino and a lepton with masses \(m_{i}\) and \(\mu_{i}\) at infinity \(t\to-\infty\), and an asymptotic state \(|j+\rangle\otimes|\ell_{f}^{-}\rangle\) of a neutrino and a lepton with masses \(m_{j}\) and \(\mu_{f}\) at infinity \(t\to+\infty\). Since the interaction is formulated in electroweak basis, the amplitude \(\mathcal{A}_{\alpha\beta}\) is calculated first for flavor neutrinos, then converted to the basis of the diagonal in mass asymptotic states:

$$\mathcal{A}_{ji}=\langle j+|\alpha\rangle\mathcal{A}_{\alpha\beta}\langle\beta|i+\rangle.$$
(93)

The summation is implied over the repeated indices. The cross section is represented as

$$d\sigma=\frac{1}{4\sqrt{(p_{i}p_{\ell^{-}_{i}})^{2}-m_{i}^{2}\mu_{i}^{2}}}\left|\mathcal{A}_{ji}\right|^{2}$$
$${}\times\frac{1}{(2\pi)^{2}}\Phi_{2}(\sqrt{s},m_{j},\mu_{f}),$$
(94)

where \(p_{\ell^{-}_{i}}\) is the initial four-momentum of lepton, \(s=(p_{i}+p_{\ell^{-}_{i}})^{2}\), and \(\Phi_{2}(\sqrt{s},m_{j},\mu_{f})\) is the relativistic phase space of the final particles. The identity \(|\mathcal{A}_{ji}|^{2}=\textrm{Tr}\left[\mathcal{A}|i+\rangle\langle i+|\mathcal{A}^{{\dagger}}|j+\rangle\langle j+|\right]\) allows the cross-section for the non-polarized particles to be written in terms of the Frobenius covariants (27):

$$d\sigma=\frac{1}{4\sqrt{(p_{i}p_{\ell^{-}_{i}})^{2}-m_{i}^{2}\mu_{i}^{2}}}$$
$${}\times\textrm{Tr}\left[\mathcal{A}\mathbb{F}_{i+}\mathcal{A}^{{\dagger}}\mathbb{F}_{j+}\right]\frac{1}{(2\pi)^{2}}\Phi_{2}(\sqrt{s},m_{j},\mu_{f}).$$
(95)

The trace can be calculated in flavor basis. Explicit use of a basis diagonal in neutrino masses is not mandatory. The corresponding expressions under the trace sign contain only the neutrino mass matrix. These conclusions are valid for both Majorana neutrinos and Dirac neutrinos.

After completing this paper, our attention was drawn to the papers [47–51], which discuss the problem of constructing the mixing matrix from the mass matrix of neutrinos, including the use of analogues of the Frobenius covariants in a less general framework. Equation (28) of the paper [49] is reproduced from Eqs. (60) and (77) with the use of the perturbation theory in the limit of small mean-field potential \(|\mathbf{p}||V_{\alpha}|\ll m_{i}^{2}\) and \(|\mathbf{p}|\gg m_{i}\). In Eqs. (65) and (84) we give similar expressions in the limit of \(|\mathbf{p}||V_{\alpha}|\gg m_{i}^{2}\) and \(|\mathbf{p}|\gg m_{i}\), which is interesting for physics of newly born neutron stars. In the paper [51], an interesting theorem is discussed on the relashionship of the module of the mixing matrix elements with the eigenvalues of the mass matrix and the eigenvalues of the diagonal minors of the mass matrix. The Frobenius covariants, apparently, successfully solve the problem of constructing the mixing matrix, including determining the phases of matrix elements. In the simplest case of neutrinos in vacuum, Eqs. (27) and (39) give, in particular,

$$U_{\alpha i}U_{\beta i}^{*}=\langle\alpha|\mathbb{F}_{i+}|\beta\rangle.$$
(96)

The main message of our paper, on the other hand, is that in neutrino physics, all physical characteristics of the processes are expressed through the neutrino mass matrix without reference to the mixing matrix or assumptions about the hierarchy of neutrino masses. The proposed procedure works for the Majorana mass term, a complex symmetric matrix (7) specified by nine parameters for three neutrino flavors. Of course, an observation of the \(0\nu\beta\beta\)-decay or an assumption about the value of \(M^{M}_{ee}\) is necessary.

5 CONCLUSION

In light of the lack of an explanation for the observed violation of \(CP\) symmetry, the abundance of baryons in the universe, the presence of dark matter, etc., the limitations of SM are obvious. The characteristics of neutrinos, such as their nature (Majorana or Dirac), and the mass matrix with a large number of parameters, can point the right path when looking for SM generalizations. This paper demonstrates that the amplitudes of neutrino-related processes can be constructed by methods, which do not require a mixing matrix and instead allow explicit representations to be built purely in terms of the neutrino mass matrix. Sections 2 and 3 analyze in detail the involved processes of neutrino oscillations in a vacuum and in a medium. Section 4 demonstrates that the probabilities of the most typical processes described by diagrams with external neutrino lines are also expressed in terms of only the neutrino mass matrix without reference to the mixing matrix. We discuss quasi-elastic neutrino scattering, and single- and double-\(\beta\) decays. Although the scattering theory is formulated in terms of asymptotic states, the transformation of scattering amplitudes into the basis of massive diagonal states is not mandatory.

Compared to methods using a mixing matrix, the neutrino mass matrix can be sensitive to alternative experimental data sets. Indirect determination of the mass matrix through the mixing matrix involves loss of information, which leads typically to an increase in statistical errors. The components of the neutrino mass matrix can apparently be determined with reduced errors if the measurement data are directly compared with the mass matrix. Using such a fitting, it is potentially possible to gain access to the structure of the mass matrix — Majorana or Dirac.