Dipole portal and neutrinophilic scalars at DUNE revisited: the importance of the high-energy neutrino tail

We estimate the sensitivity of the DUNE experiment to new physics particles interacting with neutrinos, considering the dipole portal to heavy neutral leptons and a neutrinophilic scalar with lepton-number $2$ as examples. We demonstrate that neutrinos from the high-energy tail of the DUNE flux, with energies $E_{\nu}\gtrsim 5-10\text{ GeV}$, may significantly improve the sensitivity to these models, allowing to search for particles as heavy as $\simeq 10\text{ GeV}$. We also study the impact of the so-called tau-optimized neutrino beam configuration, which slightly improves sensitivity to the new physics models considered here. For both models, we consider new production channels (such as deep-inelastic scattering) and provide a detailed comparison of different signatures in the detector.

In a broad class of Standard Model (SM) extensions, an interaction of SM neutrinos with new particles in the few GeV mass range is introduced. Such models may be probed at neutrino experiments, where a highintensity neutrino beam is produced, see e.g. [1]. Examples of currently running and future neutrino beam experiments are T2K [2,3], MiniBooNE [4], Micro-BooNE [5] and DUNE [6]. Other promising facilities to look for such new particles are LHC-based experiments such as SND@LHC [7], FASERν [8], or beam-dump experiments like SHADOWS [9] and SHiP [10].
The goal of this paper is to demonstrate the potential of the DUNE experiment to search for interactions of new physics particles with neutrinos. We do this by revising the sensitivity of DUNE to two example models -a neutrinophilic scalar [11,12] and the dipole portal to heavy neutral leptons (HNLs) [13,14]. We improve on previous studies [11,12,15,16] in several aspects. First, we include the high-energy neutrino tail in the calculation of the production rate, showing its importance to extend the reach in the mass of new physics particles, and we study the effect of the so-called tau-optimized beam configuration. In particular, we demonstrate that, depending on the model, DUNE may have sensitivity to new physics particles with masses up to O(10 GeV) -a few times larger than what it was obtained in the previous literature. Furthermore, for both models, we consider new production channels for the new particles. For the dipole portal, we discuss additional detection signatures; the ratio of the different signal types is a specific prediction of the dipole portal, allowing us to identify this model in case of detection.   [17,18], where Stransverse, L det , L to det , and n nucl stand for the fiducial cross-sectional area of the detector, the fiducial length of the detector, the distance to the proton collision point, and the nucleon number density, respectively.
The paper is organized as follows. In Sec. II, we discuss the flux of neutrinos at DUNE and stress the importance of its high-energy tail for the production of heavy particles. In Sec. III, we consider the dipole portal to HNLs, reestimate the DUNE sensitivity by considering new production channels such as deep inelastic scattering and new search signatures. In Sec. IV, we discuss the neutrinophilic scalar portal, revising the sensitivity in a similar way. Finally, we conclude in Sec. V. Technical details, expressions for matrix elements and cross sections, as well as further discussion can be found in several appendices.

II. HIGH-ENERGY NEUTRINOS AT DUNE
A. The DUNE experiment DUNE [6] is the next-generation long-baseline neutrino oscillation experiment, primarily aiming at precisely measuring neutrino oscillation parameters, such as the unknown CP violation phase in the lepton sector. On the other hand, it is also a powerful tool to search for a variety of new physics beyond the SM. By using a 120 GeV proton beam onto a graphite target with a beam power of 1.2 MW, this experiment can provide 1.1 ×10 21 PoT/year, generating a large flux of light mesons such as π ± , K ±,0 . A fraction of these mesons then would decay inside a 194 meter long decay pipe, producing neutrinos. As a result, DUNE provides the world's most intense neutrino beam with a wide range of neutrino energies peaking at about 2.5 GeV. The interactions of neutrinos are supposed to be studied in two detectors -the near detector (ND) and the far detector (FD). The former is designed to include a 67.2 ton Liquid Argon Time Projection Chamber (LArTPC), a magnetized gaseous argon time projection chamber, and a large, magnetized monitor, which will be located 574 m downstream of the neutrino target following the decay pipe. The FD will be equipped with four 10-kton LArTPC modules at a distance of 1285 km from the target; one can refer to Fig 1.4 in Ref. [17] for the layout of the modules. The useful parameters of ND and FD are summarized in Table I.
In order to maximize the flux of neutrinos/antineutrinos in the direction of the DUNE detectors, the charged particles need to be collimated with respect to the beam axis by a magnetic horn ν e , ND, CP-opt. ν μ , ND, CP-opt.  Figure 1: Neutrino fluxes as defined in Eq. (1) of νe,µ,τ at the DUNE ND and the flux of ντ at the DUNE FD (due to νµ → ντ oscillations), assuming neutrino operating mode of the focusing horns. ντ fluxes are magnified by a factor 10 5 for better visibility. The results for two beam configurations are shown: the CP-optimized (solid), and the tau-optimized (dashed). Note that the ντ flux at the ND does not depend on the horn configuration, since it originates from decays of promptly decaying Ds mesons and τ leptons.
system. In dependence on the operating mode (the one for neutrino or antineutrino), correspondingly positively or negatively charged particles would be collimated. The flux in antineutrino mode is very similar to the one in neutrino mode (see e.g., Fig. 5.4 in [6]), and we expect similar sensitivities to the new physics models discussed below for neutrino and antineutrino beam modes. To be specific, in this work, we will only consider the neutrino mode.
In addition, there are two horn configurations considered by the DUNE collaboration: the CP-optimized configuration, which maximizes the flux at E ντ < 5 GeV important to study CP violations in neutrino oscillations, and the tau-optimized configuration, for which the flux of τ neutrinos at the FD with E ν < 5 GeV is somewhat lower, but the higher-energy flux between 5 and 10 GeV gets significantly increased, which would result in an order of magnitude higher number of neutrino deepinelastic scattering (DIS) events [6].

B. Neutrino flux at DUNE
We define the neutrino flux at a given detector as where N PoT is the number of protons on target, S transverse is the transverse area of the detector (Table I), dN ν /dE ν is the differential distribution of neutrinos traveling in the direction of the detector. The fluxes at ND and FD for various neutrino flavors, assuming the neutrino operating mode, are shown in Fig. 1. Below, we briefly discuss their main characteristics.
The electron and muon neutrinos at DUNE are produced mainly by decays of light long-lived mesons such as π ± , K ± , K 0 L , and muons. The low-energy part of the spectra of ν µ (E νµ 6 GeV) and ν e (E νe 10 GeV) originates from decays π + → ν µ + µ + , µ + → ν e +ν µ + e + (2) correspondingly. The high-energy tail comes from decays of kaons K + , K 0 L . The relative suppression ν e /ν µ comes from the fact that muons are long-lived, τ µ /τ π/K ∝ 10 2 , and only a small fraction of them, 10 −2 , decays inside the decay pipe before being scattered/absorbed in the material after the decay pipe. To obtain the fluxes of these neutrinos, we use the publicly available results of the detailed GEANT4 [19][20][21] based simulation (G4LBNF) of the LBNF beamline developed by the DUNE collaboration [6,22]. Technical details are given in Appendix A.
When considering the DUNE sensitivity to new physics interacting with ν e,µ , we will only consider the ND. The reason is that the fluxes of ν e/µ at the FD are much smaller than at the ND due to the much smaller angular coverage of the FD. To argue this point, let us compare the products of the flux times the detector volume V det , Φ × V det , at the ND and FD. Using that V FD /V ND ≈ 800 if assuming four FD module (Table I), and neglecting the oscillations for the moment, one has where (L to ND /L to FD ) 2 ≈ 2 · 10 −7 is the solid angle suppression (see Table I). Thus, the FD is less relevant for searching for hypothetical BSM particles interacting with these flavors. Let us now discuss τ neutrinos. At the DUNE ND, their main production channels are the decays and their charge conjugated channels, where X denotes lepton or hadron final states. Since D s and τ decay promptly, their distribution (and hence the flux of τ neutrinos and antineutrinos at ND) is not affected by the horn configuration. However, the flux of these D s /τoriginated ν τ s orν τ s at the FD would be too small compared to the flux of ν τ originated from the oscillations ν µ → ν τ , and therefore the former can be safely ignored.
Only the latter (ν τ from oscillations) would be present at the FD when we calculate the sensitivities of neutrino mode to new physics. The DUNE simulations did not include D s mesons and τ leptons. To generate the τ neutrino flux, we have used the spectrum of D s mesons from [23], then simulated the decay chain (4) (see the Appendix A), and selected the τ neutrinos that point to the ND.
The relative flux suppression ν τ /ν µ at the ND parametrically behaves as (5) Here, P pp→Ds 4 · 10 −6 , P pp→π + ≈ 6.3 are multiplicities (the number of decayed mesons per PoT), taken from [24], Br(D s → τ ) ≈ 0.055 [25], and p Eν is a factor depending on the neutrino energy, varying from O(1) at E ν 3 GeV to 10 3 for E ν 50 GeV. The reason of this behavior is that τ leptons (and hence τ neutrinos) have much larger mean energy than pions and kaons (and hence ν e/µ ) decaying inside the decay pipe.
Apart from the production in decays of mesons and τ leptons, ν τ may also be produced via oscillations ν µ → ν τ . This channel is not relevant for the ND, since the typical neutrino oscillation length is much larger than the distance from the target to ND. However, it becomes the main contribution to the ν τ flux at the FD (see Sec. II). To obtain the oscillated flux at the FD, we have first extracted the flux of ν µ at the FD, and then convoluted it with the oscillation probability assuming ∆m 2 = 2.523 × 10 −3 eV 2 and L ≈ 1300 km for all neutrinos. The resulting flux is where P νµ→ντ is O(1) at the FD. Since the oscillation probabilities ν µ → ν τ andν µ →ν τ are the same up to small CP-violating effects, independently on the oscillating neutrino energy, the property of the relative suppression of the antineutrino to neutrino fluxes in the neutrino mode translates to the fluxes of τ neutrinos and antineutrinos at the FD. Therefore, the ν τ fluxes at the ND and FD have a different origin. As a result, the conclusion that the FD is not relevant for the DUNE sensitivities to new physics becomes invalid in the case of interactions with τ flavor. To illustrate this, let us again compare the products of the flux times the detector volume. For the CP-optimized horn configuration, using Eqs. (5), (7), one has Considering the FD and E ν 5 GeV, both P osc and p Eν can be of O(1), leading to the ratio in Eq. (8) being very large. Hence, the FD may provide better sensitivities to new physics coupling to ν τ [15]. However, since P osc (E ν ) ∼ E −2 ν at energies E ν 10 GeV, together with large p Eν they cause the suppression of this ratio at large energies. As a result, with the increase of the neutrino Detector Horn conf.
The total numbers of neutrinos traveling in the direction of the DUNE detectors (with the parameters given in Table I) per PoT are given in Table II. C. High-energy neutrinos and production of new physics particles In Refs. [12,15], which have studied the sensitivity of DUNE to new physics particles produced in scatterings of electron and muon neutrinos, artificial cuts on the neutrino spectrum have been imposed: E ν < 6 GeV and E ν < 10 GeV correspondingly. This is reasonable if one studies the production of particles Y with mass m Y E ν,max . Indeed, first, the production of such particles does not require large energies. Second, only a tiny fraction of ν e , ν µ have energies E ν > 5 − 10 GeV, see Fig. 1.
However, as we show below, these cuts can significantly underestimate the maximal mass of new physics particles that may be searched for at DUNE. Namely, the sensitivity of DUNE estimated in [12,15] rapidly drops at masses m Y = 2.5 − 3 GeV, which is directly related to the cuts. To understand this, let us look closer at the relevant Y production processes (here without specifying the model details). Correspondingly, they are 1 1 Apart from the scatterings off nucleons, [15] considered scatterings off electrons and nuclei, but the former requires much larger neutrino energies to produce a particle with the given mass, while the latter is suppressed due to nuclear form factors for large masses m r −1 nuclear . ν,min , Eq. (11), required to produce a particle Y in the scattering processes in Eq. (9), on the Y mass. The shortand long-dashed black lines denote correspondingly the maximal energies of νe,µ and ντ obtained in the simulation (see Fig. 1).
where n, p are nucleons, and in the second process a threshold for the invisible transverse momentum of p T,Y > p T,min = 0.5 GeV is required to suppress backgrounds (we will discuss this in more detail in Sec. IV).
The minimal energy E ν,min of the neutrino required to produce the particle Y in these processes is, correspondingly, The behavior of E (1), (2) ν,min as a function of m Y is shown in Fig. 2. We see from the figure that at DUNE it actually may be possible to produce much heavier particles -up to m Y 8 GeV (the production from ν e/µ ), or to m Y 10 GeV (from ν τ ) if neutrinos with energies up to 100 GeV are taken into account, offering a potential trade-off to the reduced flux at high energies.
In the next two sections, we will study how the DUNE sensitivity to the mentioned models extends in detail. We will consider both ND/FD CP-optimized and tauoptimized horn configurations.

III. NEUTRINO DIPOLE PORTAL
The effective Lagrangian of the neutrino dipole portal below the electro-weak (EW) scale is [13] L dipole = d αN σ µν P L ν α F µν + h.c., (12) where ν α is the SM neutrino of flavor α = e, µ, τ , σ µν = i 2 [γ µ , γ ν ], P L = (1 − γ 5 )/2, F µν = ∂ µ A ν − ∂ ν A µ is the electromagnetic field strength tensor, and N is a Heavy Neutral Lepton (HNL). Note that this Lagrangian is not gauge-invariant and therefore valid only at energies below the EW scale, above which we need to consider UV-completions of the operator in Eq. (12). Here we remain agnostic about the UV origin of this new interaction and study its phenomenological implications at energies below the EW scale.

A. Phenomenology at DUNE
The DUNE sensitivity to HNLs with a dipole portal interaction has been studied previously in Refs. [15,16]. The HNL production mechanism studied there was mainly quasi-elastic (QE) neutrino up-scattering [13], where T = e, n/p, Ar or atomic nuclei in the crust along the trajectory of the neutrino beam.
Here, we will include also the deep-inelastic (DIS) contribution, ν α + p/n → N + X, where X is an arbitrary hadronic state, see Fig. 3 (a). Furthermore, in addition to the neutrino up-scattering, we consider the production of HNLs by decays of short-lived and long-lived mesons [13], π ±,0 , K ±,0 , η, ρ 0 , see Fig. 3 (b), (c) for two example diagrams. A detailed discussion of the relevant cross sections and production rates is given in Appendix B 1. The HNL production may occur either inside the detector (neutrino up-scattering) or outside it (both neutrino up-scattering and meson decays). In the latter case, HNLs need to reach the detector in order to decay inside it and hence be detected. As discussed in detail in Appendix B 1 d, in most cases the production from meson decays plays only a sub-leading role and the main production channel is either inside or outside up-scattering. The only exception is HNL production via d τ at the ND.
The main HNL decay channels are with l = e, µ. The dominant channel is the mono-photon channel, with the decay width being for Dirac neutrinos. Above the di-electron and di-muon mass threshold, the leptonic channel becomes available. This channel has been considered in ref. [35] in the context of testing the dipole portal at FASER and in ref. [40] in the context of the T2K near detector. These processes are sub-dominant. However, they have the advantage that it is possible to reconstruct the decay vertex since we have two charged particles. In the limit m N 2m l , the decay width behaves as The branching ratios of the leptonic decay modes are shown in Fig. 4. At large masses m N m l , the suppression of the leptonic decay width with a factor α EM /3π 10 −3 compared to the photon channel (15) gets partially compensated by the logarithm.
The two decay modes lead to different experimental signatures and imply different requirements for background rejection. Furthermore, if both decay channels can be observed, their ratio is a specific prediction of the model, serving as a smoking-gun signature.
The combination of the different production and decay processes leads to different signatures of the dipole portal at DUNE: 1. Monophoton -an event consisting of a single isolated photon appearing inside the detector. This type of events occur when the HNL is produced outside of the detector (via e.g. decays of mesons or by neutrino up-scatterings), then enters the detector, and decays through N → ν + γ.
2. Double-bang -an event inside the detector consisting of two vertices: the one with recoil matter particles (electrons, nucleons, nuclei) and the one with a displaced monophoton or a pair of charged leptons [16,32]. This type of signature appears if a HNL is produced inside the detector via neutrino up-scattering, and then travels a distance larger than the DUNE spatial resolution, which is of order

Prompt photon/di-lepton -a single event with
a recoil electron/nucleon/nucleus, plus a leptonantilepton pair or a photon. It occurs if an HNL is produced inside the detector and decays within ∆l DUNE from its production point, such that its production and decay points cannot be resolved at DUNE. It is the main signature for heavy HNLs m N 1 GeV.
Technical details about estimates of the number of events for various signatures and comparison of the different production channels are given in Appendix B. Depending on the signature, possible backgrounds include: scatterings ν + T → ν + π 0 + X [15], with subsequent decays π 0 → 2γ, where 2 photons cannot be resolved and the final state X either was not (for the monophoton signature) or was detected (for the prompt HNL decays and double bang signature); processes ν + T → n + X where X has been detected while the neutron, being undetected, produced a displaced monophoton track at lengths comparable to its absorption length (for the double bang signature [16]). In the case of prompt events, a possible way to discriminate the background is to look at the energy of the recoil particles. Namely, considering the neutrino scatterings, this energy is typically much larger since the scattering medi-ators are heavy W/Z which do not restrict the transferred momentum, while for the dipole portal it is a photon γ, which prefers q 2 → 0. A detailed study of background discrimination goes beyond the scope of this paper. Instead, we will show the reach of DUNE in the form of iso-contours corresponding to N events = 2 events per 1 year of DUNE operation.

B. Discussion of results
Our sensitivity results are summarized in Fig. 5, which shows the iso-event contours for various signatures at DUNE ND and FD assuming the CP-optimized horn configuration. (We comment on the tau-optimized configuration later in this section.) From the figure, we conclude that depending on the signature and the flavor of the neutrino coupled to the HNL, DUNE may probe masses m N up to 10 GeV. In addition, the lower bound of the sensitivity may even touch the supernova bounds.
Let us discuss the results in more details. At the ND, the mono-photon signature (caused by decays of mesons and the outside up-scatterings) dominates at masses m N < 300 MeV − 1.7 GeV, depending on the neutrino flavor coupled to N . The reason is that (i) the outside medium is more dense and (ii) more length is available for neutrinos to scatter compared to the ND (see a discussion in Appendix B 1 d). For heavier HNLs, however, the number of monophoton events rapidly decreases, and the signature becomes sub-dominant. The reason for this is that the HNL lifetime shortens with increasing HNL mass, l decay, and HNLs that are produced outside the detector decay before reaching the ND. To increase l decay,N , one needs to decrease d 2 . However, the decrease leads to the suppression of the number of produced HNLs. As a result, the conditions of having large enough decay length and large enough production rate become inconsistent at m N 1 GeV.
The di-lepton signatures, with both prompt and displaced pairs of leptons, are sub-dominant due to the small branching ratio Br(N → l + l − ν). However, being combined with the events with prompt photons, they may help in discriminating signal events from background, and would provide additional information about HNL dγ LEP 95% CL dZ,γ LEP 95% CL  12), from top to bottom, respectively. Iso-contours correspond to 2 events of decaying HNLs per year at the DUNE ND for de and dµ or ND and FD for dτ . Various signatures are considered: monophoton (green solid), "double-bang" (green dashed), di-leptons (red solid and red dashed), and all these signatures together (blue). The excluded parameter space is shown in gray, taken from [15,46]. The DUNE sensitivity obtained in [15] corresponding to the combination of all signatures except for di-leptons is shown by the light gray curves for comparison. The bounds from LEP depend on the UV completion of the model, which we illustrate by showing the constraints assuming zero couplings to Z bosons ("dγ") and assuming dZ = dγ tan(θW ) ("dZ,γ"). properties.
The sensitivity of the DUNE FD to d = d e,µ is strongly suppressed compared to the ND. This is because of the much smaller neutrino flux of ν e,µ at the FD compared to the ND, according to Eq. (3). Therefore, we do not show the sensitivity of the FD to these couplings.
For d = d τ , there is a complementarity between the ND and FD, as illustrated in Fig. 5 (bottom right). As we have already discussed in Sec. II, ν τ s at the ND are much less numerous than at the FD, but much more energetic on average. As a result, the FD allows probing much smaller couplings than the ND for masses m N 1 GeV, where the energy of neutrinos is not important. However, for m N 1 GeV, where the HNL production requires large neutrino energies, the sensitivity of FD quickly drops. In contrast, this is not the case for ND, which allows probing much larger masses in the domain that is currently not excluded by past experiments.
An important remark has to be made. The region of new parameter space of heavy HNLs that may be probed by DUNE is model-dependent, especially for the case of the coupling to ν τ . Namely, bounds from LEP depend on the UV completion of the model (12). In particular, a class of the completions introduces a coupling d Z to Z bosons, which is related to d in Eq. (12) via the algebra of electroweak symmetry generators [13]. In Fig. 5, we indicate the model dependence of the bounds by showing the LEP constraints assuming d Z = 0 and d Z = 0. If d Z = 0 and assuming coupling to ν τ , DUNE may only probe HNLs with masses m N 4 GeV, while if d Z = 0, the sensitivity extends to m N 9 GeV.
For comparison we show in Fig. 5 also the sensitivity from Ref. [15], which combines all signatures except for the di-lepton one. We can appreciate the improved sen-sitivity due to additional production channels and the inclusion of the high-energy tail of the neutrino flux. 3 We note also that when adopting the same assumptions about exposure and experimental parameters, our results for the double-bang signature are in good agreement with Ref. [16].
Above, we have considered the HNLs as the Dirac particles. The analysis for the Majorana case is completely similar. The only difference in the sensitivity is purely numerical: given the same coupling d, the Majorana HNLs would have twice larger decay width Γ N . The effect on the sensitivity depends on the considered parameter space and the signature. Namely, in the regime where the HNL decay length is much larger than the scale of DUNE, the number of events is N events ∝ d 2 Γ N , hence being twice larger for Majorana HNLs. Since Γ N ∝ d 2 , in the plane d − m N , the lower bound of the sensitivity would be 2 1/4 lower for the Majorana particles. In the opposite regime of short decay lengths, the number of events behaves as N events ∝ d 2 × exp[−l min Γ N /cγ N ], where l min is the minimal displacement for the given signature. For the prompt signatures, l min = 0, and the decay width does not affect the sensitivity at all (the case of prompt signatures). This is not the case for the case of the outside production and the double bang. The upper bound of the sensitivity to these signatures would be a factor of √ 2 lower for Majorana particles. Finally we note that apart from these numerical differences, the presence of a Majorana mass for the HNL would lead to theoretical inconsistencies, as in this case, the dipole operator with d in the relevant range would induce a large contribution to the light neutrino mass matrix via a loop diagram, in conflict with bounds on neutrino masses, see e.g., [13]. Therefore, a consistent theory of large d in case of Majorana NHLs would require more elaborate model building.
Let us now consider the impact of the different horn configurations, in particular the tau-optimized flux. The fluxes of all the neutrinos except for ν τ at the ND are very sensitive to the horn configuration. 4 In Fig. 1, we compare the fluxes for the CP-optimized and tau-optimized configurations. We see that the amount of neutrinos with energies E ν > 5 GeV for the tau-optimized flux may be larger by a factor up to 20. This translates to a comparable increase of the cross section of the production 3 Apart from the analysis improvements discussed above, we have also corrected a few numerical issues in the calculations of Ref. [15], the most relevant being (i) the treatment of the offaxis neutrino flux and (ii) correcting the number of nucleon number density in the detector (roughly a factor 3). Item (i) explains different results concerning outside up-scattering production both for ND and FD sensitivities (see a brief discussion in Appendix B 1 d). 4 τ neutrinos at the ND are produced mainly by decays of Ds mesons and τ leptons, which occur almost instantly after the production of these particles inside the target. Therefore, the ντ flux is insensitive to the horn configuration.
of heavy HNLs with m N 2 GeV, which requires highenergy neutrinos, see Fig. 6 (top left). The comparison of the sensitivities of DUNE assuming the CP-optimized and tau-optimized horn configurations is shown in the remaining panels of Fig. 6. In dependence on the HNL mass, the improvement in the sensitivity may reach a factor as large as 3.

IV. NEUTRINOPHILIC SCALAR PORTAL
The neutrinophilic scalar portal is characterized by the effective Lagrangian [11,12] where ν c means charged conjugated neutrino field, and φ is a massive complex scalar carrying lepton number −2.
The interaction in Eq. (17) can emerge from a dimension-6 operator of the form (LH)(LH)φ/Λ 2 , with L and H denoting the lepton and SM Higgs doublet fields, respectively. Again we remain agnostic about the UV completion and focus on the low-energy implications of the effective interaction (17). The scalar introduced in Eq. (17) can also act as a portal to dark matter [12]. SM precision measurements of meson decays or invisible Higgs decays provide constraints on the parameters m φ and g α [11]; the excluded range is shown by the gray shaded region in Fig. 8.
At neutrino experiments, the neutrinophilic scalar portal in Eq. (17) may be searched for by an excess of events with missing transverse momentum p T [12,49], e.g.
where X is a hadronic final state, and φ leaves the detector invisibly or decays into neutrinos and hence leads to missing p T (see Fig. 7). Indeed, since neutrinos are collimated with respect to the beam axis, they carry negligible p T . Therefore, the products of the reaction (18) have total zero p T . Finite reconstruction efficiencies smear the p T distribution even if all the reaction products are detected. However, the resulting distribution is softer than the distribution of events with invisible particles such as φ. Additionally, due to the lepton number violation in the operator (17), an antimuon is produced together with hadrons. If the muon charge can be identified, this would allow reducing background even further. In particular, the study performed in [12] has shown that the distribution of p T of muons and nucleons produced by the process ν + n → p + + µ − , assuming the muon charge identification, drops to zero at p T 0.5 GeV.
In [12], the sensitivity of DUNE to the neutrinophilic scalar portal has been estimated by considering the quasielastic (QE) scattering requiring that the missing transverse momentum is p T,miss > 0.5 GeV. In this work, we extend the analysis  of [12] in the following ways. First, we include the highenergy tail of the neutrino flux, by considering also neutrinos with E ν > 10 GeV. This will provide sensitivity to somewhat larger scalar masses (remind Fig. 2). Second, we include form-factors to the QE proton or neutron vertex [50] in order to take into account that especially for heavy φ, m φ r −1 p 1 GeV, the nucleons cannot be considered as point-like, leading to a suppression of the QE cross section. Third, we consider in addition the production process via DIS, which dominates when the incoming neutrinos are energetic enough. In Appendix C we provide details on the relevant cross sections and further discussion of the respective impact of each of these improvements in the analysis. The 90% CL sensitivity to the coupling g µ , considering the scalar interactions with ν µ , is shown in Fig. 8. We require p T,φ > 0.5 GeV and assume the absence of background, which, according to Ref. [12], corresponds to the case when the muon charge may be identified (see also the discussion above). Therefore, the sensitivity is calculated by requiring N events > 2.3 for the full DUNE exposure of 5 years in the neutrino mode (5.5×10 21 PoT). If assuming no muon charge ID, the 90% sensitivity would drop by a factor of 2 because of the non-zero background [12].
In this figure, we also show the ND sensitivity from Fig. 2 in [12], where we choose the curve obtained under the assumption of the presence of the muon charge identification. As for the cases of the interactions with ν e and ν τ , DUNE has no sensitivity. There are two reasons for this. First, in case of interactions with ν e,τ , electrons and Excluded [1 90 1. 01 25 9] Quasi-el, CP-opt. Quasi-el+DIS, CP-opt. Quasi-el, τ-opt. Quasi-el+DIS, τ-opt. For comparison, the gray line shows the sensitivity obtained in Ref. [12] for the CP-optimized horn configuration.
τ leptons would be produced in the process (19) instead of muons. It is more difficult to look for missing p T in this case, since it is more complicated to reconstruct the kinematics of electrons and τ (which in addition decays into τ neutrinos). Second, even if assuming a perfect kinematics reconstruction of e and τ , the fluxes of ν e,τ are suppressed by at least two orders of magnitude (see Fig. 1). This translates into at least one order of magnitude suppression of the sensitivity. On the other hand, the currently excluded domain for the e, τ cases remains very similar to the case of the interaction with ν µ or even becomes tighter (see [12]). Let us briefly summarize the results. Compared to the estimate of the sensitivity from [12], which has been performed for the CP-optimized horn configuration, at small masses m φ 1.5 − 2 GeV, our estimate shows slightly worse sensitivity. This is the effect of the inclusion of the nucleon form factors. The inclusion of DIS can only partially compensate for this effect, see blue solid and blue dashed curves. However, due to the inclusion of highenergy neutrinos, we find that DUNE has a sensitivity to somewhat larger scalar masses, up to m φ 4.5 GeV. Compared to the CP-optimized horn configuration, the tau-optimized configuration improves the sensitivity for the whole mass range (green curves). In addition, it extends the probed mass region by 1 GeV.
The point-like nucleon approximation may also affect the analysis of background performed in [12], according to which the SM background vanishes at p T > 0.5 GeV. In particular, form factors suppress large transverse momentum of nucleons (and hence muons), which means that in reality the SM background is even softer and drops to zero at smaller p T . Therefore, the cut p T,miss > 0.5 GeV for the signature (18) may be relaxed. This question is a subject of a separate study; however, we have checked that the relaxation of the cut to p T,miss > 0.2 GeV may improve the sensitivity by up to a factor of 2.

V. CONCLUSIONS
In this paper, we have estimated the sensitivity of the DUNE experiment to new physics particles that interact with neutrinos, considering the neutrino dipole and neutrinophilic scalar portals as examples. Compared to previous studies, we have included new production channels, considered various detection signatures, and taken into account the high-energy tail of the DUNE neutrino flux. Although the fraction of neutrinos with higher energies is strongly suppressed, neutrinos with E ν 10 GeV are important for the production of heavy particles with masses m 1 GeV, see Sec. II C. We have demonstrated that the inclusion of high-energy neutrinos may significantly improve the sensitivity, extending the probed mass range. In particular, DUNE may probe the dipole portal up to HNL masses m 9 GeV (Fig. 5), and the neutrinophilic scalar portal up to masses m 5.5 GeV (Fig. 8), with some dependence on the horn configuration.
For the dipole portal, we have studied in detail various production mechanisms and signatures in the detector. Depending on the parameter region, the mono-photon signal, double-bang signature, or a single event with a shower + prompt photon or lepton-pair may be observable. The relative size of these signals will be a smoking gun signature to identify the dipole portal model. Furthermore, different signatures may require dedicated analysis cuts and background mitigation strategies.

Appendix A: Flux of muon and electron neutrinos at DUNE
In order to get the flux of ν e and ν µ neutrinos and the distribution of decaying mother mesons, we have used the publicly available results of the detailed GEANT4 [19][20][21] based simulation (G4LBNF) of the LBNF beamline developed by the DUNE collaboration [6,22]. It is possible to extract the flux of the electron and muon neutrinos from the simulations at the DUNE ND and FD directly from the simulation files, see Fig. 9. In the simulation data, there are more than 12 million of decaying pions and 1 million of decaying kaons. These numbers are enough to produce the high-energy spectrum for ν µ very well. However, the number of simulated decays into electron neutrinos is a few orders of magnitude smaller, and this is especially the case for the high-energy tail. There are two reasons for this. First, most of the muons ( 99%) from π, instead of decaying, get absorbed outside the decay pipe because of their long lifetime. The decaying muons typically have small energies E µ 10 GeV, which means that they mostly produce low-energy neutrinos. Second, in the DUNE simulation, only a fraction of kaons, namely, Br(K → ν e )/Br(K → ν µ ), has been used to simulate decays into electron neutrinos. This results in large wiggles in the obtained distribution for energies E ν 30 GeV.
Instead of using this direct flux, we have re-simulated the decays of π, K, µ into ν e , allowing all the kaons can decay into neutrinos, and then selecting the neutrinos that point to the ND. To simulate the decay, we first extracted the momenta and decay coordinates of the mesons and muons from the DUNE simulation dataset, then generated the phase space of their decay products at the rest frame of the decaying particles, and then boosted to the lab frame.
We have used the following matrix elements: for the muons, pions, and kaons [51] correspondingly, where f ± are K → π transition form-factors [52]. As a crosscheck, we have verified that the phase space distribution of particles in the process µ + → e + + ν e +ν µ at the rest frame of decaying µ + generated by us coincides with that generated by MadGraph5 [53].
The comparison of the obtained fluxes with the directly extracted fluxes is shown in Fig. 9. They agree with each other very well at all energies in consideration, while the regenerated flux for ν e has fewer wiggles at high energies. Small differences at energies E ν 10 GeV may be explained by slightly different sizes and geometry of the ND considered in the DUNE simulations and in our estimates.  The decay widths are computed as [25] Γ for 2-body decays, where p Y is the momentum of the Y particle in the rest frame of decaying X, and where m ab = (p a + p b ) 2 . In the limit of m N,l m meson , the decay widths are The resulting branching ratios under the assumption |d| = 1 GeV −1 are shown in Fig. 10, where we use the SM decay widths of the mesons from [25]. The largest branching ratio is for the decay K − → N + γ + l − . It is larger than the ratio for the decay π − → N + γ + l − due to the scaling Γ m − →N +γ+l − ∝ m 5 m − . Indeed, the SM decay widths of the charged mesons scale instead according to Γ m − ,SM ∝ m 3 m − . Therefore, we have Br m − →N +γ+l − ∝ m 2 m − . To estimate the flux of HNLs from the mesons in the direction of the DUNE ND, we have followed the same procedure as described in Appendix A.

b. Neutrino up-scattering -inside
Neutrino up-scattering is a process ν + T → N + X, where T is a target particle, and X denotes an arbitrary final state. In this study, we consider the following processes: T = e, X = e; T = n/p, X = n/p; T = Ar, X = Ar; T = n/p, X = hadrons . (B15) The first three sets correspond to quasi-elastic (QE) scattering off electrons, nucleons, and Ar nuclei. The last set corresponds to deep inelastic neutrino up-scattering (DIS). The QE differential cross sections dσ T N,prod (E ν , Q 2 )/dQ 2 , with Q 2 being the modulus-squared of the momentum transferred to the target, have been previously computed, see e.g. [15]. For the estimates of the number of events, it is useful to compute the cross sections averaged over neutrino energies at the DUNE ND, where f ν (E ν ) is the neutrino distribution function: The behavior of these cross sections as a function of the HNL mass and for different target particles is shown in Fig. 11, (top panel). At small HNL masses, the production from scattering off nuclei dominates, which is because of the enhancement of the cross section (per nucleon) by the factor Z 2 /A. However, the nuclear EM form factor exponentially suppresses large momentum transfers (and hence the production of heavy HNLs), and at masses around m N 1 GeV, the scattering off nucleons starts to dominate. At mass m N 8 GeV (12 GeV) for the production from ν e/µ (ν τ at ND), the production via scattering off nucleons reaches the kinematic threshold and instantly turns off. The only remaining production channel is the scattering off nuclei, although it is exponentially suppressed.
To calculate the DIS cross section, we have first implemented the model of the dipole portal in MadGraph5 [53] with the help of FeynRules [55,56]. Then, we generated the process with further interfacing to pythia8 for subsequent showering/hadronization. As a cross-check of the model implementation, we have reproduced the result of the analytic formula from [15] for the production of the HNL in the ν up-scattering off electrons. As a cross-check of showering, we have reproduced the DIS neutrino cross section ν + n → p + X from [57] within the systematic uncertainty for the neutrino energy range 1 GeV < E ν < 40 GeV. The comparison of the DIS cross section with the QE up-scattering cross sections for various HNL masses and neutrino energies, as well as the cross sections averaged over neutrino energies, are shown in Fig. 11. Before averaging over neutrino energies, the DIS cross section dominates for large HNL masses and, simultaneously, if neutrino energies are large enough to allow momenta transfers Q 2 1 GeV 2 to the nucleons (bottom left panel). However, the neutrino high-energy tail is strongly suppressed, and after averaging over neutrino energies, this suppression compensates for the improvement, see Fig. 11 (bottom right).

c. Neutrino up-scattering -outside
In addition to the production by up-scatterings inside the detector, HNLs may be produced outside, when neutrinos scatter off the crust separating the decay pipe and the ND, see Fig. 1.2 in Ref. [6].
For simplicity, we will drop the contribution of DIS, and estimate the outside production via QE processes. In this way, our estimate is conservative.
A simple estimate of the flux of HNLs produced by the outside up-scatterings may be obtained under the assumption that HNLs travel in the same direction as neutrinos. However, for neutrino energies E ν 5 GeV that dominate the flux, this is not the case: the typical angle between the neutrino and the produced HNL may significantly exceed the angular size of the detector seen from the end of the decay pipe, see Fig. 12. The broadening of the HNL flux means that in general an HNL pointing to the ND may be produced by a neutrino that does not travel to the ND, and vice versa -a neutrino that travels to the ND may produce an HNL not reaching the detector.
To estimate the flux of HNLs at the ND, we have followed a procedure similar to the one described in Appendix A. Namely, for all neutrinos produced by meson decays (not only the neutrinos in the direction of the ND), we generated a random scattering point located between the end of the decay pipe and the beginning of the ND. Then, we generated the produced HNL angle/energy according to the angle/energy distribution given by the differential up-scattering cross section, and selected only those HNLs that point to the detector. For simplicity, we assume that the crust is made of 28 Si. Having this dataset, we have computed the total number of HNLs traveling to ND, as well as the distribution function f (upsc) The calculation of the flux of HNLs produced by outside up-scatterings at the FD is technically more complicated, as it requires the flux of neutrinos flying by very small angles θ ∼ 10 −6 − 10 −4 rad. To calculate such flux, we have computed the distribution of muon neutrinos in angles and energies semi-analytically, using the approach of Ref. [58] from Appendix B. We have verified the obtained distribution by comparing the fluxes of neutrinos at ND and FD with the results of the DUNE simulations and finding an excellent agreement [22]. Having the neutrino fluxes, we have followed the same procedure as for ND.

d. Comparison of different production channels
Before calculating the sensitivity, it would be useful to understand which production channel dominates the sensitivity. In this sub-section, we make the comparison between the fluxes of HNLs from up-scatterings and decays of mesons using simple estimates.
a. Outside and inside up-scattering. Let us first compare the number of HNLs produced by up-scattering inside and outside the ND. The comparison is non-trivial, since not all the HNLs produced outside may be able to reach the detector volume: it depends on their decay length l N,decay = cτ N p N /m N , where τ N = Γ decay,N is the HNL lifetime (see Sec. III A). Therefore, instead of the total number of the produced HNLs we consider the (effective) number of produced HNLs: Here, N ν,tot is the total number of neutrinos produced at DUNE. geom is the fraction of neutrinos producing HNLs that travel to the ND; for the inside scatterings, it is just the fraction of neutrinos flying to the ND. σ N,prod is the mean scattering cross section (per nucleon) of these neutrinos; n nucl is the nucleon number density; L scatt,eff is the length available for scatterings, accounting that HNLs have to survive long enough and decay inside the detector. For the outside events, we assume the crust number density n nucl ≈ 1.4 · 10 30 m −3 [25], and where L crust = 304 m or 1300 km is the distance between the end of the decay pipe and the beginning of the ND/FD [6]. It accounts for the fact that only neutrinos scattering in locations closer than l decay,N effectively contribute to the flux of HNLs that can reach the detector. We assume p N = 10 GeV as the characteristic momentum of short-lived HNLs.
The difference in σ The behavior of the ratio (B21) for the case of ND is shown in Fig. 12. We see that the geometric acceptance provides a O(1) correction. For the behavior of this ratio at the FD, the situation is qualitatively similar. The ratios (B21) contradict the results of the paper [15], according to which the outside up-scattering does not contribute to the sensitivity of the ND at all. We note that the outside fluxes have been calculated differently: using the pre-computed off-axis neutrino flux in the GLoBES format and the analytic formula in [15], whereas in this work we re-simulate the neutrino flux and the HNL production event-by-event. It is worth mentioning that the GLoBES files include the high-energy neutrino tail, and hence the difference cannot be attributed to the latter. We have not found an explicit reason for the discrepancy between the results, but there is a simple argument showing that the number of events from outside up-scattering at the ND is at least comparable with the inside flux: let us consider a volume of the crust material having the volume equal to the volume of the ND, and placed just in front of it. Neutrinos up-scatterings inside this volume would produce an HNL flux being a factor of (L det /L scatt,eff ) × (n crust nucl /n ND nucl ) 3 larger than the inside flux. Given the location of the considered crust volume, all of these HNLs would reach the ND and provide a decay signal, in the relevant limit where the decay length of the HNLs is large compared to L det .
b. Production from mesons and outside up-scatterings. The HNL production from mesons is relevant only at ND. The reason is that all mesons decay before the end of the decay pipe, and to be able to reach the FD, the produced HNLs have to travel the distance of O(1300 km). Since the HNL decay length cτ N γ N ∝ d −2 , this is possible only if d is very small, which a priori suppresses the flux of the HNLs.
Let us now compare the amount of HNLs produced by up-scattering, and those produced by mesons and flying in the direction of the ND. Again, we should compare an effective number of particles, defined by Here, N m,tot is the total number of mesons m produced at DUNE, geom is the fraction of produced HNLs that travel to the detector and exp[−L min,eff /l decay,N ] accounts for the exponential suppression of the flux of short-lived HNLs produced by decays of mesons. In particular, for d = d e,µ , HNLs are produced by long-lived mesons, and we fix L min,eff ≈ 221 m, which is the end of the decay pipe [6]. For d = d τ , where the production comes from promptly decaying mesons, we fix L min,eff = L to det = 574 m. The ratio of the produced HNLs resulting from Eq. (B22) is shown in Fig. 12. From the figure, we see that except for the production via d τ , for which the ν τ flux is strongly suppressed, the production from neutrino upscattering dominates. The reason is the relative suppression of the probability to produce HNLs by decays of mesons (L crust = 304 m): in dependence on the meson which may decay into an HNL, see Figures 11, 4.

Number of events at DUNE
The number of HNL events for the case of non-zero coupling d α consists of two parts: which corresponds to the contribution from neutrino up-scattering and decays of mesons. The number of events from inside up-scatterings per one year is Here, N ν,detector = S transverse Φ ν (E ν )dE ν is the total number of neutrinos flying in the direction of the detector, with S transverse being the detector's transverse area (see Table I; note that for the FD we include four modules, which in total have twice larger transverse area than one module). d σ (T ) N prod Eν /dE N is the differential production cross section (per nucleon) averaged over neutrino energies. n det nucl is the nucleon number density (which corresponds to the liquid agron for both ND and FD). P decay (E N ) · L is the HNL decay probability averaged over length of the detector L available for neutrino scatterings: where the factor L comes from the production probability in the neutrino up-scatterings, P prod = σn nucl L, and . The parameter L min is either 0 or ∆l DUNE = 1 cm, depending on the considered signature (displaced or prompt decays). In the limit l N,decay L min , this expression simplifies to Finally, Br(N → channel) is the branching ratio of the HNL decay into the given particle state, which is νγ, νee, νµµ depending on the signature. The approximation made in (B25) is that both neutrinos and HNLs fly along the beam axis. This is justified since the DUNE ND and FD cover a small solid angle and are located on axis. For the number of events from the outside up-scatterings, we have

The shape of the sensitivity curves
The combined sensitivities to inside up-scattering events shown in Fig. 5 behave non-trivially with the HNL mass. Namely at small HNL masses m N 2 GeV, the lower bound d lower decreases, reaches its minimum at some mass m N,peak , and then starts increasing. The reason is the following. Schematically, the number of events is N events = N N,prod × P decay . Since P decay ≤ 1, the sensitivity cannot cover couplings smaller than defined by the condition N N,prod = 2. The number of produced HNLs is N N,prod = f (m N ) · d 2 , where f (m N ) ∝ σ N,prod is the production cross section. As a result, we may define the minimal coupling d N prod =2 = 2/f (m N ). At HNL masses m N 1 GeV, d N prod =2 remains practically constant. The reason is that the cross section weakly depends on m N in this range (Fig. 11).
In the domain of small HNL masses and at d N prod =2 , the HNL decay length l N,decay = cτ N γ N ∝ m −4 N d −2 is parametrically very large, l N,decay L to det . Therefore, the decay probability is P decay ≈ L det /l N,decay 1 (remind Eq. (B27)). Plugging this expression in the expression for the number of events, we get which is a growing function of the HNL mass. This scaling explains why the lower bound decreases. However, with the increase of the HNL mass, (i) the HNL decay length at fixed d grows, (ii) f (m N ) drops, which leads to the increase of d N prod =2 . As a result, P decay (d N prod =2 ) becomes O(1), and mass dependence of the lower bound becomes to be determined only by the behavior of f (m N ), which quickly drops at these masses: d lower ≈ d N prod =2 . The position of the minimum of the lower bound, which we denoted as m N,peak , may be estimated as N events (m N , d N prod =2 ) 2.   Figure 13: Left panel: the dependence of the cross section of the process (C1) and the DIS production cross section on the neutrino energy for different choices of φ masses. To illustrate the effect of the form factors entering the effective nucleon vertex (C3), we also show the cross section in the limit (C9). Right panel: Impact of QE form factors, the high-energy neutrino tail, and the inclusion of DIS production on the the sensitivity of the DUNE ND to the neutrinophilic scalar portal (17). For comparison we also show the sensitivity obtained in [12].

b. Deep inelastic scattering production
The DIS production process is where X denotes any multi-hadron state. To calculate the DIS cross section, we have implemented the Lagrangian from Eq. (17) in MadGraph5 using FeynRules. To check the implementation, we have also included the weak interaction of point-like nucleons (C9) and reproduced the quasi-elastic scattering cross section (19). The behavior of the DIS cross section with neutrino energy is shown in Fig. 13. Similarly to the case of the dipole portal, the DIS cross section is smaller than the quasi-elastic cross section at small neutrino energies but becomes larger at high energies.

Number of events and discussion of sensitivity
The number of events for the neutrinophilic scalar portal at the DUNE ND has the form where n nucl is the number density of nucleons, L det is the DUNE ND detector length, and σ i are the production cross sections for the processes (C1) (C10), where we require p T,φ > 0.5 GeV. To obtain the sensitivity, we require N events > 2.3, corresponding to 90% C.L. in background free regime. The comparison with [12] under the same assumptions is shown in Fig. 13 (blue solid versus grey curves in the right panel); we find good agreement. To illustrate the effects of the inclusion of the form factors and of neutrino energy cut, in the figure, we include sensitivity curves under different assumptions about them. We see that the form factors affect not only the domain of large φ masses but even the domain of small masses m φ 1 GeV. This is caused by the transverse momentum cut p T,φ > 0.5 GeV, which sets a lower bound on the momentum transferred to the nucleons. The suppression caused by form factors increases with the growth of the transferred momentum, and therefore the sensitivity worsens. The comparison of the short-dashed and long-dashed red curves in Fig. 13 (right panel) shows the impact of including the DIS production channel.