Interaction of supernova neutrinos with stochastic gravitational waves

We examine the propagation and flavor oscillations of neutrinos under the influence of gravitational waves (GWs) with an arbitrary polarization. We rederive the effective Hamiltonian for the system of three neutrino flavors using the perturbative approach. Then, using this result, we consider the evolution of neutrino flavors in stochastic GWs with a general energy density spectrum. The equation for the density matrix is obtained and solved analytically in the case of three neutrino flavors. As an application, we study the evolution of the flavor content of a neutrino beam emitted in a core-collapsing supernova. We obtain the analytical expressions for the contributions of GWs to the neutrino fluxes and for the damping decrement, which describes the attenuation of the fluxes to their asymptotic values. We find that the contribution to the evolution of neutrino fluxes from GWs, emitted by merging supermassive black holes, dominates over that from black holes with stellar masses. The implication of the obtained results for the measurement of astrophysical neutrinos with neutrino telescopes is discussed.


I. INTRODUCTION
Neutrinos provide a unique possibility to explore physics beyond the standard model with the help of nonaccelerator methods. Such achievements of neutrino physics became possible after the observation of oscillations of atmospheric and solar neutrinos [1]. These experimental facts are the straightforward indications of the nonzero masses and mixing between different neutrino flavors. External fields, e.g., the neutrino interaction with background matter [2] and electromagnetic fields [3], are known to modify the process of neutrino oscillations.
The gravitational interaction, in spite of its weakness, was found in Refs. [4][5][6] to contribute to the neutrino oscillations dynamics. In the majority of studies, neutrino oscillations in curved spacetime were examined when particles move in static gravitational backgrounds, e.g., in the vicinity of a black hole (BH). It is interesting to analyze the propagation and oscillations of astrophysical neutrinos in time dependent metrics, e.g., induced by a gravitational wave (GW).
The studies of the fermions evolution in GWs were carried out in Refs. [7,8]. Neutrino spin oscillations, i.e. when we deal with transitions between active left polarized and sterile right polarized particles, in background matter under the influence of GW were discussed in Ref. [9]. Neutrino flavor oscillations in GWs, as well as in gravitational fields caused by metric perturbations in the early universe, were considered in Ref. [10].
In this paper, we continue the research in Refs. [11,12], where the influence of stochastic GWs on neutrino flavor oscillations was considered. The main problem of Refs. [11,12] was the consideration of astrophysical neutrinos emitted in decays of charged pions. Although such neutrinos form the major cosmic neutrinos background, their sources are distributed rather uniformly in the universe. Thus, the integral flux in a terrestrial detector should be averaged over the propagation distance of such neutrinos. This fact makes it difficult to separate the contribution of stochastic GWs on the measured flavor composition. To avoid this difficulty, we decide to examine the effect of stochastic GWs on supernova (SN) neutrinos. If an explosion of a core-collapsing SN happens in our galaxy, firstly, it emits a sizable neutrino flux to be measured even by existing neutrino telescopes [13]. Secondly, SN is almost a point-like neutrino source. Hence one should not average over the neutrino propagation distance. In this situation, we expect that the effect of stochastic GWs is not smeared.
The present work is motivated by the recent direct detection of GWs by the LIGO-Virgo collaborations [14]. There are active multimessenger searches of GWs and high energy neutrinos by existing detectors [15,16] and suggestions to implement them in future ones [17]. There are also attempts to observe stochastic GWs [18,19], with various methods for their detection being developed [20].
SN neutrinos were reliably detected in 1987 after the SN explosion in the Large Magellanic Cloud (see, e.g., Ref. [21]). Since then, the experimental techniques in construction of neutrino telescopes made great achievements. Now, if a nearby SN in our galaxy explodes, a huge number of events will be recorded [22]. As mentioned above, a simultaneous detection of GWs and SN neutrinos may be possible. Besides a direct neutrino signal from a certain SN, all collapsing stars in the universe emit neutrinos which form the diffuse SN neutrino background. There are prospects to measure it by existing and future neutrino telescopes (see, e.g., Ref. [23]).
This work is organized in the following way. In Sec. II, we formulate the problem of the propagation of flavor neutrinos in a plane GW with an arbitrary polarization. Then, we derive the equation for the density matrix for flavor neutrinos if we deal with stochastic GWs. This equation is exactly solved for the arbitrary energy spectrum of GWs. Then, in Sec. III, we apply our results for the description of the interaction of SN neutrinos with stochastic GWs. We find the corrections to neutrino fluxes and the damping decrement in an explicit form. Finally, we summarize our results in Sec. IV. The effective Hamiltonian for flavor oscillations under the influence of GW is rederived in Appendix A.

II. EVOLUTION OF FLAVOR NEUTRINOS IN THE GWS BACKGROUND
The system of three active flavor neutrinos ν λ , λ = e, µ, τ , with the nonzero mixing, as well as under the influence of a plane GW with an arbitrary polarization, obeys the following Schrödinger equation: The mixing matrix U can be present in the form, where c ab = cos θ ab , s ab = sin θ ab , θ ab are the corresponding vacuum mixing angles, and δ CP is the CP violating phase. The values of these parameters can be found in Ref. [24]. The Hamiltonian H 1 in Eq. (2.1), which describes the neutrino interaction with GW, has the form is the Hamiltonian in the mass basis, A c = 1 2 sin 2 ϑ cos 2ϕ cos[ωt(1 − cos ϑ)], A s = 1 2 sin 2 ϑ sin 2ϕ sin[ωt(1 − cos ϑ)], h +,× are the amplitudes corresponding to 'plus' and 'cross' polarizations of GW, ω is the frequency of GW, ϑ and ϕ are the spherical angles fixing the neutrino momentum with respect to the wave vector of GW, which is supposed to propagate along the z-axis. To derive Eq. (2.3) we assume that [11] ωL|v a − v b | ≪ 1, a, b = 1, 2, 3, where L is the distance of the neutrino beam propagation and v a is the velocity of a mass eigenstate. Analogously to H (vac) m , we subtract the common diagonal term in H m for a circularly polarized GW with h + = h × was obtained in Ref. [11] based on the exact solution of the Hamilton-Jacobi equation for a test particle in a plane GW. In the present work, we provide a more straightforward perturbative derivation of the same result which is given in Appendix A; cf. Eq. (A7). Of course, the expression for H (g) m coincides with that in Ref. [11] in the limit h + = h × . Now we consider the situation when a neutrino interacts with stochastic GWs. In this case, the angles ϑ and ϕ, as well as the amplitudes h +,× , are random functions of time. To study the neutrino motion in such a background, it is more convenient to deal with the density matrix ρ, which obeys the equation, iρ = [H 0 + H 1 , ρ]. Following Ref. [25], we introduce the density matrix in the interaction picture, ρ int = exp(iH 0 t)ρ exp(−iH 0 t). It satisfies the equation, where H int = exp(iH 0 t)H 1 exp(−iH 0 t). Using the Baker-Campbell-Hausdorff formula and the fact that that both H We assume that stochastic GWs form a Gaussian random process. In this situation, all odd correlators of angle factors A c,s and the amplitudes h +,× are vanishing. Moreover, we take that h + and h × are independent. After averaging, the formal solution of Eq. (2.4) can be present in the form of a series, where we show only two nonzero terms in order not to encumber the text.
We can see that the series in Eq. (2.5) decays into two independent ones corresponding to different polarizations of GW. Each of these series contains only either h + (t)h + (0) or h × (t)h × (0) . In the following, we account for all terms in the expansion in Eq. (2.5). Further analysis of each of the series corresponding to different GW polarizations is identical to that in Ref. [12]. Therefore we omit the details. Now, let us consider the averaging of the angle factors. We should mention that both the amplitudes h +,× (t) and the angles ϑ(t) and ϕ(t) are random functions of time. Indeed, we consider random distribution of GWs sources. It means that, when a neutrino interacts with a certain GW, the angle between a neutrino momentum and the wave vector of GW is randomly distributed from zero to π. However, unlike the correlators h +,× (t 1 )h +,× (t 2 ) , which are taken to be arbitrary, we suppose that both ϑ(t 1 )ϑ(t 2 ) and ϕ(t 1 )ϕ(t 2 ) are proportional to δ(t 1 − t 2 ). This supposition is reasonable since it is based on the assumption of the uniform distribution of the sources of GWs in the universe. The form of the correlators h +,× (t 1 )h +,× (t 2 ) depends on physical processes underlying the GWs production. Thus, it is inexpedient to take that the amplitudes are δ-correlated.
We can study, e.g., the correlator A c (t 1 )A c (t 2 ) , which has the form, where α 1,2 = ωt 1,2 (1 − cos ϑ 1,2 ) and the angles ϑ 1,2 and ϕ 1,2 correspond to t 1,2 . Since the random variables ϑ and ϕ are taken to be δ-correlated, we get that cos α 1 cos α 2 = 1/2. Now, we should average Eq. (2.6) over directions of incoming GWs, Analogously we show that A s (t 1 )A s (t 2 ) = 3 128 . The obtained correlators of the angular factors in Eq. (2.7) should be used in Eq. (2.5). As we mentioned above, Eq. (2.5) splits into two independent series. Accounting for all terms in the expansions and applying the results of Ref. [12], we get that ρ int obeys the equation, where In case of circularly polarized GWs with h + = h × , we reproduce the results of Ref. [12] in Eqs. (2.8) and (2.9).
To proceed with the analysis of Eq. (2.8), we define the new matrix ρ ′ = U † ρ int U , which is the density matrix in the interaction picture for the neutrino mass eigenstates. After some matrix algebra, we get that ρ ′ satisfies the equation,ρ where the initial condition ρ ′ (0) has the form, ρ ′ (0) = U † ρ int (0)U = U † ρ(0)U . It can be also expressed in the components as ρ ′ ab (0) = σ U * σa U σb P σ (0), where we assume that ρ λσ (0) = δ λσ P σ (0). Here the emission probabilities P σ (0) are proportional to the neutrino fluxes at a source: Accounting for Eq. (2.11), we get the expression for the density matrix for flavor neutrinos, The probability to detect a certain flavor, after the neutrino beam passes the distance x ≈ t, reads P λ (x) = ρ λλ (t ≈ x). Using the components of U in Eq. (2.2), it can be represented in the form, are the phases of neutrino vacuum oscillations. Equation (2.12), which takes into account the neutrino interaction with stochastic GWs, should be compared with the analogous probabilities for neutrino vacuum oscillations P (vac) λ , which are derived, e.g., in Ref. [26]. The difference ∆P λ = P , which reveals the effect of GWs on neutrino oscillations, has the form, ab | are the neutrino oscillations lengths in vacuum. If we study the interaction between stochastic GWs and neutrinos emitted by randomly distributed sources, we have to average Eq. (2.13) over the propagation distance. It gives ∆P λ = 0. Therefore, the effect of stochastic GWs on oscillations of such neutrinos is washed out. The fluxes at a source will coincide with these accounting for only vacuum oscillations. Thus, the claim in Ref. [12], that stochastic GWs result in small changes of observed fluxes of neutrinos from randomly distributed sources, is incorrect. The deviation of fluxes, obtained in Ref. [12], is likely to stem from an inexactitude of numerical simulations. We avoid this inexactitude in the present work since we rely on the analytical solution of Eq. (2.8).

III. APPLICATION TO SN NEUTRINOS
In this section, we apply the obtained results for the description of the propagation of SN neutrinos in the background of stochastic GWs.
A huge amount of energy is carried away by neutrinos from a core-collapsing SN. The major neutrino luminosity was reported, e.g., in Ref. [27] to take place during a ν e -burst, which happens because of the direct Urca process e − + p → n + ν e in the neutronazing matter of a protoneutron star (PNS). This burst occurs at ∼ (3 − 4) ms after the core bounce and lasts 0.1 s [27]; see also references therein. The neutrino luminosity can reach ∼ 10 53 erg · s −1 during the burst, with almost all of emitted neutrinos being of the electron type [27].
We start this section with the study of the interaction between stochastic GWs and SN neutrinos emitted in a ν e -burst. In this situation, the fluxes at a source are F νe : F νµ : F ντ S = (1 : 0 : 0). At later moments of time, other neutrino flavors are emitted. The fluxes of different flavors of SN neutrinos become almost equal by t ≈ (0.05 − 0.1) s after the core bounce. Thus, the initial neutrino fluxes are not in the ratio (1 : 0 : 0). The evolution of SN neutrinos with such initial condition in the presence of stochastic GWs is also discussed in this section.
A collapsing star, owing to its relatively small size, can be considered as an almost point-like neutrino source. Indeed, the size of a neutrinosphere, i.e. an effective sphere where neutrinos are trapped inside, is L source 100 km at the moment of a ν e -burst. The energy of SN neutrinos is E ∼ 10 MeV [22]. The oscillations lengths are L 21 ≈ 330 km and L 31 ≈ L 32 ≈ 10 km for this energy and ∆m 2 ab from Ref. [24]. Neutrinos are emitted from any point of a neutrinosphere more or less isotropically. A terrestrial detector can register all SN neutrinos emitted towards it. Thus, we have to integrate the densities of the fluxes over the area of a neutrinosphere, S source ∼ L 2 source , and divide the result by S source . We can call this procedure as the averaging over the emission points. Thus, if we carry out this averaging for ∆P λ in Eq. (2.13), the (31)-and (32)-contributions are smeared. The only nonvanishing contribution is from the solar oscillations channel (21).
As we mentioned above, other neutrino flavors are emitted after a ν e -burst, changing the ratio (1 : 0 : 0) of the initial fluxes. However, the absolute value of the SN neutrino luminosity becomes smaller. The neutrinosphere shrinks at these greater times. Nevertheless, its size remains greater than 10 km. Thus, only the (21)-oscillations channel gives a nonzero contribution to Eq. (2.13) even for t > t burst .
We start by considering neutrinos emitted in a ν e -burst. Accounting for the ratio of the initial fluxes, we get that ∆P λ for such SN neutrinos takes the form, In Eq. (3.1), we do not set the sine and cosine factors to zero despite x ≫ L 21 . The propagation distance is huge, but it is fixed. The correlators of the amplitudes of GWs h +,× (t)h +,× (0) are related to the spectral density S(f ) of GW by [28] ij where f is the frequency measured in Hz. In Eq. (3.2), we use Eq. (A5). Then, we define the function Ω(f ) = f ρc dρGW df , where ρ GW is the energy density of a GW and ρ c = 0.53×10 −5 Gev·cm −3 is the closure energy density of the universe. Using Eq. (3.2), we get that Ω(f ) = πf 3 8ρcG S(f ), where G = 6.9 × 10 −39 GeV −2 is the Newton's constant. The functioñ g(t) has the form,g Now, choosing the source of stochastic GWs, which is fully characterized by Ω(f ), we can evaluate ∆P λ . We suppose that stochastic GWs are emitted by randomly distributed merging supermassive BHs. In the case, we can approximate Ω(f ) by [29] where Ω 0 ∼ 10 −9 , f min ∼ 10 −10 Hz, and f max ∼ 10 −1 Hz. The main contribution to Γ in Eq. (3.1) results from f min .
Hence we can put f max → ∞ since f max ≫ f min . We suppose that 0 < x < L, where L ∼ 10 kpc is the maximal propagation length, which is taken to be comparable with the Galaxy size ∼ 32 kpc. Using Eqs. where τ = x/L and ω min = πLf min are the dimensionless parameters. The function Γ(τ ) is shown in Fig. 1(a) for Ω(ω) corresponding to Eq. (3.4). In Fig. 1, we depict only the normal mass ordering case since the inverted ordering is almost excluded experimentally [30]. We can see in Fig. 1(a) that Γ tends to a constant value at τ → 1. If we study neutrino fluxes at the Earth, we put x = L, or τ = 1. Then, we suppose that the distance between a source, SN, and a detector, the Earth, is great. It corresponds to the limit ω min ≫ 1. Accounting for Eq.  (3.6) If Ω 0 = 10 −9 and f min = 10 −10 Hz, Γ ⊕ = 8 × 10 −2 in full agreement with Fig. 1(a). Although Eq. neutrinos. Indeed, in that case [29], f min ∼ 10 −5 Hz ≫ 10 −10 Hz. Hence Γ ⊕ → 0 and ∆P λ → 0 despite Ω is greater for such sources.
If (F νe : F µ : F ντ ) S = (1 : 0 : 0) for SN neutrinos emitted in a ν e -burst, the probabilities at the Earth for vacuum oscillations are [26] where we accounted for the fact that the contributions of the ( Table I. Since the fluxes F (vac) ν λ are rapidly oscillating on the distance L = 10 kpc, in Table I, we present only the mean values and the amplitudes of oscillations: (mean value ± amplitude). We also give (∆F ν λ ) ⊕ at x L, i.e. the asymptotic values, which correspond to Figs. 1(b)-1(d). The mean value of (∆F ν λ ) ⊕ equals to zero, as explained above. Thus, we show only the amplitudes of oscillations of (∆F ν λ ) ⊕ . One can see in Table I that the relative contribution of stochastic GWs to the measured neutrino fluxes is at the level of (5 − 7) %. Now we turn to the discussion of neutrinos which are emitted after a ν e -burst, i.e. at t > (3 − 4) ms after the core bounce. As we mentioned above, the ratio of the emission fluxes is not equal to (1 : 0 : 0). At these times, the fluxes of different flavors eventually become almost equal. It happens at t ≈ 0.05 s after the core bounce (see, e.g., the numerical simulation, carried out in Ref. [31]). At t > 0.1 s, the absolute values of the fluxes start to decrease [31]. Let us study the influence of stochastic GWs on SN neutrinos emitted at t burst < t < 0.1 s. We can use the general Eq. (2.13) taking that the emission probabilities are time dependent: P σ (0) → P (0) σ (∆t), where ∆t = t − t burst and t burst = (3 − 4) ms is the ν e -burst time. We can approximate P (0) σ (∆t) by the following dependence: where K = 10 s −1 is the fitting factor. The emission probabilities in Eq. (3.8) satisfy σ P (0) σ = 1 at any ∆t. The initial fluxes (F ν λ ) S ∝ P (0) ν λ (∆t) are shown in Fig 2(a). The value of (F νe ) S at ∆t = 0 corresponds to the luminosity ∼ 10 53 erg · s −1 in a ν e -burst. One can see that (F ν λ ) S of different flavors become almost equal at ∆t ≈ K −1 .
The detection probabilities in Eq. (2.13) depend on the propagation length x. We consider the maximal values of ∆P λ when neutrinos arrive to a terrestrial detector. They are where Γ ⊕ = 8 × 10 −2 . Note that ∆P (max) λ in Eq. (3.9) depends on ∆t.
The deviations of the maximal fluxes (∆F , owing to the interaction with stochastic GWs, for different post bounce emission times ∆t are shown in Fig 2(b). One can see that (∆F (max) ν λ ) ⊕ at ∆t = 0 coincide with the values given in the third column in Table I. If ∆t → K −1 , (F ν λ ) ⊕ → 0. Indeed, P

IV. DISCUSSION
In the present work, we have studied the propagation and flavor oscillations of astrophysical neutrinos interacting with stochastic GWs. We have rederived more straightforwardly the effective Hamiltonian for such a system in Appendix A. The analytical expression for the density matrix of flavor neutrinos has been obtained in Sec. II. Then, in Sec. III, we have applied the obtained results for the description of SN neutrinos.
The present research has several advances compared to Ref. [12], where the interaction of astrophysical neutrinos with stochastic GWs was also studied. Firstly, we have accounted for two independent polarizations of GWs contrary to the case of a circularly polarized GW in Ref. [12]. Secondly, now we have found the exact solution of the equation for the density matrix in the general case of three neutrino flavors. This fact allowed us to correct some statements about the asymptotic behavior of neutrino fluxes, which were made basing on numerical simulations in Ref. [12].
In Ref. [12], we studied the interaction between stochastic GWs and astrophysical neutrinos emitted in decays of charged pions. Sources of such neutrinos are distributed more or less uniformly in the universe. Thus one had to average over the neutrino propagation distance to get the fluxes at the Earth. This fact made it difficult to extract the contribution of stochastic GWs to neutrino fluxes. In the present work, we have considered neutrinos emitted in a SN explosion. It enables us not to perform the averaging over the propagation distance except for (31)-and (32)-oscillations channels. It is valid since the oscillations length L 21 is greater than the radius of the neutrinosphere, i.e. the size of the neutrino source. The distance L between a possible SN explosion is great, but it is fixed. Therefore we should not average over L/L 21 in the probabilities for neutrino flavors.
We have also obtained the analytical expression for the damping decrement Γ; cf. Eq. (3.6). This result allowed us to evaluate the contribution of other sources of stochastic GWs, e.g., merging BHs with stellar masses, to the evolution of neutrino fluxes. The straightforward derivation of the effective Hamiltonian for flavor neutrinos oscillations, presented in Appendix A, can be applied to different metrics perturbations besides GWs considered here. Using this result, we can study, e.g., the evolution of neutrinos in perturbations in the early universe [10,32].
In the present work, we have also studied the interaction between stochastic GWs and SN neutrinos emitted in the time interval t burst < t 0.1 s. When t ≈ 0.1 s, the initial fluxes are almost equal, F νe : F νµ : F ντ S = (1 : 1 : 1). We have found in Sec. III that, in this situation, the contribution of stochastic GWs to the evolution of SN neutrinos is washed out; cf. Fig. 2(b). It means that the major effect is for neutrinos emitted at a ν e -burst, which corresponds to F νe : F νµ : F ντ S = (1 : 0 : 0). The neutrino luminosity in a SN explosion remains significant up to 10 s after the core bounce, which is the time scale for the neutrino driven PNS cooling [33]. However, the influence of stochastic GWs on SN neutrinos emitted at t > 0.1 s is vanishing since such neutrinos are emitted with almost equal probabilities. Thus, in Sec. III, it is inexpedient to extend ∆t beyond the 0.1 s interval.
In Table I, we have summarized the contributions of stochastic GWs to the fluxes of flavor neutrinos at the Earth. They are at the level of a few percent. The current neutrino telescopes are able to detect up to several thousand neutrinos from SN in our Galaxy [13]. Future detectors, like the Hyper-Kamiokande, can detect about 7.5 × 10 4 such events [34]. Thus, the interaction with stochastic GWs can results in the change of the SN neutrinos fluxes by ∼ ±350 events, in case of the Super-Kamiokande, and by ∼ ±3750 events, for the Hyper-Kamiokande.

ACKNOWLEDGMENTS
I am thankful to A. V. Yudin and J. W. F. Valle for the communications, as well as to V. A. Berezin for the useful discussion. The work is supported by the government assignment of IZMIRAN.

Appendix A: Derivation of the effective Hamiltonian
The action S a (x, t) of a neutrino mass eigenstate with the mass m a , moving in the curved spacetime with the metric g µν , obeys the Hamilton-Jacobi equation, a (x) = p µ a x µ , where p µ a = (E a , p), E a = m 2 a + p 2 , and p is the constant neutrino momentum.
In the linear approximation in h µν , the first order term S where we use v a ≡ v (0) a = p/E a as the particle velocity. It means that a neutrino is supposed to move along a straight line. In general situation, one has that v a = v (0) a + δv a (t) in curved spacetime. However, since δv a ∝ h µν , this term can be neglected in Eq. (A3) because S (1) a ∝ h µν already. Finally, we get that dS (1) a dt = − 1 2Ea h µν p µ a p ν a . The contribution to the effective Hamiltonian can be obtained as (H m ) aa = dSa dt . One can check the validity of this expression in the vaccum case. Thus, we obtain that the neutrino intraction with a gravitational field, induced by a metric perturbation, contributes to the effective Hamiltonian in the mass basis as Note that we have to take h µν on the particle trajectory, h µν (x, t) = h µν (x(t), t).