Inclusive Nucleon Decay Searches as a Frontier of Baryon Number Violation

Proton decay, and the decay of nucleons in general, constitutes one of the most sensitive probes of high-scale physics beyond the Standard Model. Most of the existing nucleon decay searches have focused primarily on two-body decay channels, motivated by Grand Unified Theories and supersymmetry. However, many higher-dimensional operators violating baryon number by one unit, ∆B = 1, induce multi-body nucleon decay channels, which have been only weakly constrained thus far. While direct searches for all such possible channels are desirable, they are highly impractical. In light of this, we argue that inclusive nucleon decay searches, N → X+anything (where X is a light Standard Model particle with an unknown energy distribution), are particularly valuable, as are model-independent and invisible nucleon decay searches such as n → invisible. We comment on complementarity and opportunities for such searches in the current as well as upcoming largescale experiments Super-Kamiokande, Hyper-Kamiokande, JUNO, and DUNE. Similar arguments apply to ∆B > 1 processes, which kinematically allow for even more involved final states and are essentially unexplored experimentally.

I. Introduction 1 II. Nucleon decay operators 2 A. Operator dimension d = 6 2 B. Operator dimension d > 6 3 III. Exclusive nucleon decay searches 5 IV. Inclusive nucleon-decay searches for ∆B = 1 processes 7 A. Model-independent and invisible mode searches 8 1. Isotope fission products 8 2. Nuclear de-excitation emission 8 B. N → X+anything 10 1. N → (e ± , µ ± , π ± , K ± , ρ ± , K * ,± ) + anything 10 2. N → (π 0 , K 0 , η, ρ, ω, K * ,0 ) + anything 10 3. N → γ + anything Baryon number B and lepton number L are seemingly accidentally conserved in the Standard Model (SM), which makes searches for their violation extremely important. So far we have not observed any B or L violating processes despite decades of experimental investigation, yet there are many reasons to expect that these symmetries could be broken. The linear combination B + L is in principle already violated by 3 + 3 units within the SM itself through non-perturbative instanton effects, albeit highly suppressed [1]. From a more fundamental topdown perspective, global symmetries such as U (1) B and U (1) L are expected to be violated at some level by quantum gravity effects [2,3], which is however difficult to quantify. Furthermore, baryon-number violation is one of the key prerequisites for successful baryogenesis [4], which would address the observed baryon-antibaryon asymmetry of our Universe. Explicit B-violation and associated proton decay is a defining prediction of Grand Unified Theories (GUTs) [5,6] that unify the three forces of the SM into a single gauge group, offering an explanation for the observed charge quantization as well as gauge coupling unification. GUTs typically lead to effective dimension-six (d = 6) operators with ∆B = ∆L = 1 that induce two-body nucleon decays such as p → e + π 0 and n → e + π − , mediated by heavy gauge bosons with family-universal couplings [7][8][9].
It has to be stressed, however, that the significance of nucleon decays stretches far beyond GUTs. ∆B = 0 processes generically appear in numerous theoretical extensions of the SM, such as supersymmetry (SUSY) [8], typically mediated at the renormalizable level by leptoquarks and diquarks [10][11][12]. A somewhat model-independent approach to study ∆B = 0 operators is through the SM effective field theory (SMEFT), neglecting modeldependent interference effects [13][14][15][16][17]. This allows one to arXiv:1910.07647v1 [hep-ph] 16 Oct 2019 identify higher-dimensional operators with a flavor structure that can make rather unconventional nucleon decay channels dominant. Further, higher-dimensional d > 6 operators typically induce nucleon decays with multibody final states and even in simple UV-complete models one can encounter more complicated nucleon decay channels such as n → e + e − ν [18] or p → π + π + e − νν [12,15]. Hence, there is a vast landscape of possible motivated nucleon decay modes of varying complexity.
Experimentally, an extensive nucleon decay search program has been carried out over multiple decades, covering more than 60 decay channels [19]. The most sensitive searches, coming from the Super-Kamiokande (SK) experiment [20,21] (see Ref. [22] for a review), have already pushed the nucleon lifetime limits for certain channels above 10 34 yr [23], twenty-four orders of magnitude beyond the age of our Universe. The frontier of baryon number violation searches will be spearheaded by next-generation large-scale underground neutrino experiments, namely the Jiangmen Underground Neutrino Observatory (JUNO) [24], the Deep Underground Neutrino Experiment (DUNE) [25], and Hyper-Kamiokande (HK) [26]. It is of paramount importance to take full advantage of these considerable efforts and to identify potential new signals in order to ensure that interesting channels are not overlooked due to theoretical biases.
In this work we revisit nucleon decay channels arising from higher-dimensional operators and discuss some of the possible resulting final states. While systematically searching through all of the kinematically allowed nucleon decay channels with increased final state complexity would constitute the strongest probes, this approach quickly becomes highly impractical beyond the simplest of the modes. In view of this, we highlight the importance of inclusive nucleon decay searches. Although these searches are not as sensitive as exclusive ones looking at a particular channel, they allow to cover very broad parameter space in a model-independent manner and are practically far more feasible. This approach is particularly fruitful to revisit in view of the upcoming large-scale experiments.
This paper is organized as follows: in Sec. II we discuss higher-dimensional operators that lead to nucleon decay and argue in particular that many of them lead to multibody final states that are not covered in current searches. In Sec. III we provide a brief overview of current and upcoming detectors as well as existing exclusive nucleon-decay searches. We discuss and propose possible inclusive searches as well as model-independent signatures in Sec. IV. Sec. V is devoted to a short discussion of ∆B > 1 processes such as dinucleon decay, which would also profit from inclusive searches. Finally, we conclude in Sec. VI.

II. NUCLEON DECAY OPERATORS
Several well-motivated theoretical models such as GUTs or R-parity-violating SUSY lead to nucleon decay, typically with specific two-body channels such as p → e + π 0 or p →νK + being dominant [8,9]. In order to discuss nucleon decay in its generality without being restricted to a certain model, we instead consider various possible higher-dimensional operators that can mediate these processes, an approach 1 that goes back to Refs. [13][14][15][16][17]. We aim to determine which operators lead to twobody and which lead to multi-body nucleon decays.
A. Operator dimension d = 6 In the SMEFT, operators that exhibit ∆B = 1 start to appear at operator mass dimension d = 6. Keeping the flavor structure, these ∆B = ∆L = 1 operators can be written as where α, β, γ denote the color, i, j, k, l the SU(2) L , and a, b, c, d the family indices [13][14][15][16][17]. u, d, and are the right-handed up-quark, down-quark, and lepton fields, while Q and L are the left-handed quark and lepton doublets, respectively. The y j couplings have massdimension −2 and the first-generation entries are constrained to be < (O(10 15 ) GeV) −2 due to the induced two-body nucleon decays. Specifically, all of the above operators generate the well-constrained decay p → e + π 0 with a rate of order A variety of other two-body nucleon decay channels are induced as well, including muon and kaon modes once we consider second-generation flavor indices. Three-body decay modes with similar rates are induced as well [28] but ultimately lead to weaker constraints.
Operators in L d=6 involving either charm, top, bottom or tau are seemingly unconstrained by nucleon decay since these particles are heavier than the proton; it is however possible to go through heavy off-shell particles and still induce nucleon decay, as emphasized in Ref. [29] for operators involving a tau and more generally in Ref. [30] (see also Ref. [31] for a UV-complete example with scalar leptoquarks). As an extreme example, Ref. [30] has considered the coupling y 4 3333 and shown that at two-loop level the simple decay n →ν τ π 0 is induced with an estimated rate Γ(n →ν τ π 0 ) 1 10 33 yr An even stronger limit y 4 3333 (10 11 GeV) −2 has been estimated from p →ν τ K + in Ref. [32]. Despite the suppression by loop factors, Cabibbo-Kobayashi-Maskawa mixing angles, and the Fermi constant G F , these limits are far stronger than any constraints from ∆B = 1 top or tau decays on the same operator, making nucleon decays clearly the best search channels. Notice that any operator involving tau leptons brings a missing-energy tau-neutrino in the final state, reducing the detection efficiency somewhat; turning this around, it is then crucial to search for (two-body) nucleon decays involving one neutrino, as these are the best channels of ∆B = 1 operators that involve tau leptons. As it is conceivable that the UV completion that generates the ∆B = 1 operators also singles out tau leptons (or any other lepton flavor for that matter) one must not rely exclusively on searches involving electrons and muons [33].
The main conclusion of the above exercise is that any ∆B = 1 operator will lead to nucleon decay, no matter the flavor structure. In the d = 6 case of Eq. (1) it is furthermore always possible to close SM loops in order to generate a two-body nucleon decay. These might not always be the dominant decay modes, but in light of the clean decay channels they are clearly the preferred way to search for the operators of Eq. (1). This picture changes once we consider ∆B = 1 operators of mass dimension d > 6 as we will discuss below.

B. Operator dimension d > 6
Assuming a ∆B = 1 operator of mass dimension d ≥ 6 with coefficient Λ 4−d we can estimate the amplitude for an k-body nucleon decay as M ∼ m d−k−1 p Λ 4−d and the decay width [34] neglecting the final state masses and possible symmetry factors (see Ref. [35] for loop-suppressed nucleon decays). The decay is suppressed for large d and large k, effectively lowering the probed scales Λ. Assuming conservatively that Λ should lie above TeV in order to evade LHC constraints on the underlying colored mediator particles with quark couplings, we find that d > 14 (for k = 2) or d > 6 (for k = 15) in order to push the lifetime above 10 30 yr (which corresponds to a reasonable lower bound on the total nucleon lifetime, see discussion in Sec. IV). Clearly, nucleon decays can probe very high-dimensional operators and very high multiplicity, making it possible and necessary to go beyond two-body decays mediated by d = 6 operators. As already pointed out by Weinberg [15], ∆B = 1 operators with d > 6 can carry different total lepton number ∆L, which can be used to make them dominant over the d = 6 operators. Interesting connections between ∆B, ∆L, and the mass dimension of the operator d were proven in Refs. [36,37]: d is odd ↔ |∆(B − L)| = 2, 6, 10, 14, . . . , as well as the weak inequality for the minimum dimension d min of an operator with ∆B and ∆L. For the cases of practical interest we give d min in Fig. 1, adapted from Ref. [37]. Below we will discuss ∆B = 1 operators with d > 6 in order to establish the importance of multi-body nucleon decay searches. At d = 7 one finds ∆B = −∆L = 1 operators that induce e.g. n → e − K + or p → e − π + K + [15]. In analogy to the d = 6 operators from above one can show that all of these d = 7 operators induce two-body nucleon decays at some loop level, irrespective of their flavor structure. While these might not necessarily be the dominant decay modes, they are clearly far easier to constrain experimentally. The currently best constrained channels are p → νK + and n → − π + , while many other two-body final states with |∆(B − L)| = 2 have unfortunately not been updated for twenty years (listed below in Tab. I).
At d = 8 one finds ∆B = ∆L = 1 operators of the form uuQLφφ, ddQLφφ, and dQQ φφ in addition to simply dressing the d = 6 operators of Eq. (1) with |φ| 2 , where φ is the SM scalar doublet. 2 Once again two-body nucleon decays are the best search channels.
Starting at d = 9 the phenomenology becomes vastly more interesting. Suppressing all indices and using a very compact notation, the ∆B = 1 operators can be written in the form FIG. 1. Processes with baryon and lepton number violation by ∆B and ∆L units, respectively. We only show one example process, others are implied (e.g. nn → ππ also give n-n oscillation, pp → π + π + , and many more). "Instanton" refers to processes such as 3n → 3ν that break the same quantum numbers as non-perturbative electroweak instantons. 0ν2β (0ν4β) refers to neutrinoless double [38,39] (quadruple [40]) beta decay. Final states with neutrinos make an experimental determination of ∆L impossible, but are shown here for the sake of brevity. Also shown is the minimal mass dimension d of the underlying effective operator following Ref. [37]. In addition to total lepton number L, all operators and processes also carry lepton flavor, which turns the above 2D plot into a 4D lattice [33,41]. 1 -O 9 6 have recently been discussed in Ref. [33]; they fulfill ∆B = −∆L = 1, but since they involve three lepton fields they can carry non-trivial lepton flavor number, i.e. |∆L α | > 1 for one α ∈ {e, µ, τ }, which makes it possible to make them dominant over d = 7 operators (which carry at most |∆L α | = 1) and then lead to three-body [44] and four-body nucleon decays [33]. O 9 7 -O 9 21 will lead to ∆B = −∆L = 1 two-body nucleon decays reminiscent of the d = 7 operators, although O 9 13 -O 9 21 will dominantly induce multi-meson final states, e.g. n → − K + π 0 , due to the large number of quark fields. Finally, O 9 22 and O 9 23 violate ∆B = 1 3 ∆L = 1 [15] and thus lead to nucleon decays with at least three antileptons in the final state, typically accompanied by one or more mesons. We stress that contrary to many statements in the literature these operators do indeed generate nucleon decays despite the fact that they necessarily contain charm or top quarks, in complete analogy to the discussion in Ref. [30]. At loop level the operators O 9 22,23 with arbitrary flavor structure will induce the decays p → +νν and n → 3ν, which are experimentally constrained already.
At d = 10 there are ∆B = ∆L = 1 as well as ∆B = − 1 3 ∆L = 1 operators. The former are generically expected to be suppressed compared to d = 6 operators, except for operators that involve three lepton fields with non-trivial flavor, in particular those operators that give rise to the very clean channels p → e + µ − µ − and p → µ + e − e − , discussed at length in Ref. [33]. The ∆B = − 1 3 ∆L = 1 operators are of the form dddLLLφ and lead to multi-body nucleon decays such as n → − ννK + , p → − ννπ + π + [12,15,45], and n → 3ν, only the latter being explicitly constrained so far.
For examples of ∆B = 1 operators with d > 10 we refer to Refs. [37,45,46]. If such operators are to dominate over lower-dimensional ones they should carry a different lepton (flavor) number, e.g. |∆L| > 3. It is clear that these will lead to multi-lepton final states, typically accompanied by mesons. We note that the semi-inclusive invisible neutron decay, n → neutrinos, can carry away an arbitrary amount of lepton number and flavor via neutrinos and is thus a particularly powerful decay to probe, not least because the experimental signature is the same no matter the number or flavor of the final state neutrinos. We will come back to this channel in Sec. IV.
The calculation and comparison of all possible nucleon final states induced by a given d 6 operator is of course impractical and will not be attempted here. From the examples discussed above we can however already gleam our simple main point: current nucleon-decay searches cover but a fraction of relevant modes and should by all means be extended in order not to miss new physics. From our discussion it is clear that two-body nucleon decays are powerful probes of d = 6-8 operators, whereas d ≥ 9 operators, in particular those with |∆L| > 1 or |∆L α | > 1, lead to multi-body final states that mostly lie outside of current searches. Since nucleon decay is an extremely sensitive probe of new physics, even higherdimensional operators with d 6 can give testable signals and should be investigated. To constrain these operators one then has to either go through all kinematically allowed multi-body final states (in general without knowing the angular and spectral distributions) or focus on inclusive searches, to be discussed in Sec. IV. Before, we must first discuss the status of exclusive searches.

III. EXCLUSIVE NUCLEON DECAY SEARCHES
Several major directions have been pursued to experimentally study nucleon decays. The main points of focus are the ability of an experiment to achieve high signal-detection efficiency as well as scalability to large fiducial masses/volumes, which is particularly critical for such rare-event searches. Some of the earliest searches were performed with scintillators, pioneered by Reines et. al. [47][48][49][50][51]. These exclusive searches, focusing on some specific nucleon decay channels, assumed that a final state of the decay contains an electromagnetically interacting particle. Current large-volume scintillator neutrino experiments based on carbon 12 C, such as KamLAND [52] and Borexino [53], have fiducial masses of 1 kton and can be also utilized for nucleon-decay searches (e.g. Refs. [54,55]). The main distinguishing feature of these experiments is a very low sub-MeV energy threshold.
Other early searches looked for fragments of fission induced by nucleon decays within radioactive ore. As we further discuss below, a particular advantage of such searches is that they are completely model-independent. However, the amount of source material is limited and hence they are significantly less sensitive than dedicated state-of-the-art exclusive searches. For exclusive searches, fine-grained (typically iron-based) calorimeters, such as in the Kolar Gold Field (KGF) [56], NUSEX [57], Soudan [58], and Fréjus [59] experiments, have been also utilized. In these experiments one can achieve high signal efficiency due to precise reconstruction and tracking. However, these configurations do not scale well and hence suffer from a limited fiducial mass ( 1 kton).
WC experiments can identify charged particles with high efficiency and are readily scalable, enabling fiducial volumes far greater than any other technique. The current most sensitive proton-decay limits, already exceeding lifetimes of 10 34 yr, come from the state-of-the-art WC experiment SK [23], which boasts a fiducial volume of 22.5 kton and has been collecting data for two decades.
WC detection relies on the presence of electromagnetically-interacting particles in the final state of nucleon decays. The observable signature from a particle traversing the detector is a Cherenkov radiation ring, classified either as "showering" for e ± and γ or "non-showering" for µ ± and π ± . WC searches are particularly sensitive when the final state is fully visible and one can efficiently reconstruct the original parent nucleon, such as in the case of the leading GUT modes p → e + π 0 and p → µ + π 0 , where the pion is identified via π 0 → γγ. To produce Cherenkov radiation, charged particles of mass m traveling in a medium of refraction index n must exceed the Cherenkov momentum threshold of p th = m/ √ n 2 − 1. In water, n 1.33 and thus p th = 1.14m, which translates to a required minimum momentum of 0.58 MeV, 121 MeV, 159 MeV, and 563 MeV for e ± , µ ± , π ± , and K ± , respectively. Hence, kaons from single nucleon decays, which have energy below m p /2 470 MeV where m p is the proton mass, are always invisible in WC detectors and one must rely on reconstructing products from their subsequent decays. Thus the typically leading SUSY GUT decay channel p → νK + suffers from a low detection efficiency in WC detectors. Scintillation detectors have an advantage over WC in that there is no Cherenkov threshold and a higher light yield. However, because emission of scintillation light is nearly isotropic the directional information is lost. Further, while scintillators are scalable, they are still behind in volume compared to leading WC experiments.
Several next-generation large-scale experiments will allow to push nucleon decay searches even further. The upcoming successor of SK, the HK WC experiment, is expected to have a 187 kton fiducial volume in its initial configuration [26]. This will allow HK to probe nucleon lifetimes up to 10 35 yr. Since a large fiducial volume benefits all the decay modes, improved limits across the board are expected. On the scintillator front, the liquid scintillator experiment JUNO will have 20 kton fiducial volume [24] and will start taking data in a few years. The DUNE [25] experiment with ∼ 40 kton fiducial volume based on liquid-argon timeprojection-chamber (LArTPC) technology will allow to track both light deposits as well as charge. Despite their smaller size compared to HK, the JUNO and DUNE detectors are particularly promising for decay modes involving kaons, such as p → νK + , which are fully visible and can be reconstructed with high efficiency. Additionally, for the multi-body decay channels discussed in this article the efficiency in WC experiments could be decreased, as it becomes difficult to reconstruct multiple overlap-ping Cherenkov rings, making alternative technologies as those in JUNO and DUNE crucial. For example, in ∆B = 2 n-n oscillations the resulting n is subsequently captured on n or p and produces a slew of decaying pions [64,65]. In DUNE, the multi-track nature of such decays could be utilized with high efficiency [66].
We briefly note that planned next-generation darkmatter direct-detection experiments such as Argo [67] (based on liquid Argon) or Darwin [68] (Xenon) could be able to achieve ultra-low sub-keV energy thresholds combined with a O(10-100) ton fiducial mass and could in principle also be utilized for nucleon decay searches. However, their fiducial volume is still multiple orders of magnitude smaller than dedicated large-scale neutrino experiments and would have a hard time competing with their exclusive searches, although their low threshold could be beneficial for model-independent searches. Despite significant experimental efforts proton decay or any other ∆B = 0 process has not been observed thus far. Almost all kinematically-allowed two-body nucleon decay channels have been searched for in various experiments, although several of the limits are outdated. We collect the strongest limits in Tab. I and also specify the violation in (B − L) symmetry that such decays will induce. We note that a neutrino or an antineutrino of any flavor in the final state will escape the detector before interacting and is hence invisible, which can result in a variation of |∆(B − L)| by 2 units. We do not include ∆B = 1 multi-nucleon decays such as pn → e + n or pp → e + ∆ + here, some of which were searched for in Ref. [59]. We also omit showing processes that violate Lorentz or charge symmetry.
In two-body decay modes, kinematics and phase-space considerations uniquely determine the resulting energymomentum distribution of the resulting final-state particles, making these searches highly model independent. For multi-body decays, additional model-dependence from dynamics comes into play, as discussed for example in Ref. [77]. Multi-body searches can furthermore suffer from lowered detection efficiencies and enhanced systematic uncertainties, e.g. from multiple hadronic/nuclear interactions. Only a subset of all kinematically allowed three-body nucleon decay channels have been searched for so far, collected in Tabs. II and III along with the |∆(B − L)| structure for the modes. To remain agnostic regarding theoretical models, a uniform phase-space energy and momentum distribution for final-state particles is typically assumed in such searches.
We are not aware of any searches involving more than three particles in the final state so we will not provide a table listing all kinematically allowed modes. Instead, we will discuss inclusive searches, which could in principle be sensitive to arbitrary multi-body final states. TABLE I. Kinematically allowed two-body nucleon decays with 90% C.L. upper limit on the partial decay width Γ. Here, ν can be a neutrino or antineutrino of any flavor, which does not change the observation signature. The column |∆(B −L)| indicates the violation of B − L in the decay, which depends on whether ν is a neutrino or an antineutrino. An asterisk denotes a limit that has been translated using the properties of WC detectors (e − ∼ e + ∼ γ, µ − ∼ µ + , π − ∼ π + ) but was not given explicitly in the references. Unconstrained channels are still subject to limits from inclusive searches, discussed in Sec. IV. See text for more details.   As we have shown in Sec. II, without focusing on specific models many ∆B = 1 operators induce nucleon decays with a multi-body final state. More so, these channels could be dominant and often the simpler two-body decay modes are completely forbidden, e.g. by a lepton number or flavor symmetry. It is therefore necessary to broaden the nucleon-decay searches, which have primarily been focused on two-body channels, in order to also cover multi-body channels. While this can be done by performing exclusive searches looking at specific multibody channels, an exhaustive search program quickly becomes impractical as the more complex modes are considered. Hence, fully model-independent as well as inclusive N → X + anything channel searches that are less dependent on the number and details of final-state particles are particularly important. We outline such searches below.

A. Model-independent and invisible mode searches
Model-independent searches allow us to probe all possible nucleon decay channels simultaneously and are not constrained to any specific final-state particle. The most stringent model-independent limits come from nuclear de-excitation emission searches. Such limits apply to all seemingly unconstrained channels described in Tabs. I, II, and III.
Model-independent searches for neutron decays also constitute the main method to study invisible decays that leave no other trace in the detector, such as n → invisible, where the invisible final state could consist of an arbitrary combination of neutrinos and antineutrinos. Invisible modes decaying to SM particles, such as n → 3ν, could become significant in models based on extra dimensions (e.g. [81,82]) or partially-unified Pati-Salamtype constructions [83]. Proposals for invisible decays to more exotic states, such as dark fermions [84] or unparticles [85], have been also suggested. Applying these searches to processes like p → 3ν can test electric-charge conservation.

Isotope fission products
A model-independent search for nucleon decays can be performed by analyzing fission isotope products due to nucleon decay within an element [86]. If the fission fragment is stable (unstable), a geochemical (radiogenic/radiochemical) search method can be used.
In the radiogenic/radiochemical approach one looks for radioactive remnants due to nucleon decay within some favorable and abundant isotope, such as 58 Ni → 57 Co [87] or 39 K → ( 38 Ar, 38 K) → 37 Ar [88]. These experiments must be placed deep underground to suppress background isotope production due to cosmic ray interactions. Further, as radioactive decay happens during experimental observation large quantities of material are required. Typically, these experiments are not easily scalable, with limits on nucleon decay lifetimes of τ p ∼ 10 26-27 yr already requiring ton-scale fiducial mass.
Due to their ultra-low keV-level thresholds and well understood backgrounds, large DM experiments could also constitute efficient detectors for radiogenic nucleon decay searches. Using ∼ 7 kg of fiducial mass and ∼ 6 kg · yr of exposure, a limit of τ p ∼ 10 24 yr was obtained with the DAMA liquid Xenon scintillator [89]. A future large-scale Xenon experiment with fiducial mass of O(10-100) tons, such as DARWIN [68], might thus achieve an improved sensitivity reach by several orders of magnitude.
A variation on the radiogenic approach is to search for appearance of fragments from spontaneous fission of a large nucleus into several sizable components. Here, using 7 kg · hr exposure of 232 Th a limit of τ p ∼ 10 21 yr was set in Ref. [90]. However, using this method effectively requires good background understanding and large quantities of very heavy elements.
In the geochemical approach one searches for stable nuclear residue associated with nucleon decays that have accumulated over billions of years within naturally abundant ore. This approach has been employed for double beta-decay studies, especially when the resulting isotope is a noble gas. When applied to nucleon decay by considering processes such as 23 Ne → 22 Ne, 39 K → 38 Ar, and 133 Cs → 132 Xe, a nucleon decay lifetime sensitivity of ∼ 10 23 yr can be achieved [91], assuming spectroscopic resolution of 1 part in 10 8 (10 ppb). An improved lifetime limit of ∼ 10 25 yr was achieved in Ref. [92] by considering beta decays in 130 Te → 129 Sb → 129 Xe, measured double beta-decay lifetime and known abundances as well as atmospheric contamination of various xenon isotopes in telluride ore. Due to the very long lifetime of the source material, these searches are not limited by fiducial mass as radiochemical searches are. These searches require precise impurity identification. While purification levels of material, such as 130 Te crystals produced for double beta-decay studies, in current experiments has greatly improved in recent years and approaches 1 part in 10 14 [93], obtaining accurate understanding of accumulated background contributions is highly non-trivial.
While we envision that some improvement on such searches is feasible, it will be difficult for these searches to compete with nuclear de-excitation emission searches.

Nuclear de-excitation emission
Aside from analyzing fission fragments, another general and model-independent signature of nucleon decay is nuclear de-excitation emission. Every nucleon decay results in a hole within the host nucleus, e.g. in oxygen 16 O (as relevant for WC-based detectors) or carbon 12 C (as relevant for scintillator-based detectors). Such decays from an inner nuclear shell will typically leave the residual nucleus in an excited state, which subsequently de-excites by emission of secondary (n, p, α, γ) particles and also possibly further via β-decay of the residual radioactive nuclei. Emission steps can be identified with nuclear-shell models and energetic considerations. The de-excitation particles accompanying nucleon decay could be searched for as a model-independent signature [94][95][96][97]. Such signals, however, could be plagued by typical radioactive background within experiments as well as neutrino interactions that result in similarly excited nuclear states [98].
For WC experiments, de-excitation γ's are especially important for these searches. De-excitation γ emission with energies of a few MeV is expected with a large branching ratio (e.g. 44% for 6.18 MeV γ-ray for 15 O [95]). However, this energy range is plagued by backgrounds associated with solar neutrinos and cosmogenics. Placing the experiment very deep underground (e.g. 6000 meter water equivalent for SNO+) and using a directional cut to suppress solar neutrino contributions, the SNO+ scintillator experiment was able to achieve a limit on invisible proton decay of with only 0.9 kton of water and 0.58 kton·yr exposure obtained during its water-filled WC phase [99].
Large WC experiments that have significant exposure can take advantage of a relatively clean signature from the higher energy γ-rays, which however have suppressed branching ratios. In particular, search for γrays in the 19-50 MeV energy range has been used to set a lower partial-lifetime limit associated with n → 3ν of τ /Br ∼ 2 × 10 31 yr by the Kamiokande experiment using ∼ 8.5 kton·yr exposure [100]. With the associated cumulative branching ratio for such energetic γ's of (1.4 ± 0.7) × 10 −4 [96], the resulting lifetime limit is τ ∼ 5 × 10 26 yr [100]. With over 370 kton·yr of exposure available for SK, a limit of order 10 28 yr can be achieved using this method.
A coincidence multi-event signature associated with several spatially and temporally correlated secondary particles from de-excitation can be employed to differentiate signal from background (see discussion in Ref. [96]). Such low energy MeV-level signals associated with de-excitation can be particularly advantageous for scintillator-based experiments, as considered by Borexino [54] and KamLAND [80]. In particular, a limit on invisible nucleon decay of τ inv. n > 0.58 × 10 30 yr (10) was obtained with the 3.2 kton KamLAND scintillator detector and 0.84 kton·yr of data [80]. In WC-based experiments, a de-excitation proton or an α would be below the Cherenkov threshold and thus invisible. A kicked-out neutron will quickly thermalize and diffuse, subsequently being captured on a hydrogen atom with an accompanying emission of a low-energy 2.2 MeV γ-ray. A low-efficiency neutron tagging is currently being used in SK to gather this signal [101], which is below the energy threshold of SK. The planned SK upgrade involving gadolinium doping [102] will allow to detect neutrons with high efficiency due to ∼ 8 MeV γ-ray emission accompanying neutron capture on gadolinium. Since de-excitations will generally result in neutrons 3 , they could be beneficial for WC searches. However, without additional coincidence signatures the signal-background discrimination is non-trivial as neutrons will also generically appear from neutrino interactions such as inverse β-decay ν e + p → n + e + . Previously, the importance of detecting neutrons associated with proton decay have been emphasized for heavy water (D 2 O) searches (e.g. in SNO) [103]. In particular, the proton decay search d → n+? requires less stringent implicit hypotheses regarding the stability of the resulting daughter nuclei and is hence even less model-dependent than some of the other searches discussed in this section.
Coincidence tagging of de-excitation emission together with some visible final state can provide additional discrimination of nucleon decay signal versus backgrounds, as first suggested by Totsuka for SK [94]. However, this is often difficult to utilize in practice. For example, in WC searches the nearly-prompt emission of a low-energy deexcitation γ could be hard to resolve when super-imposed with an energetic Cherenkov ring e.g. from a final state e ± or µ ± . A time separation between the energy depositions, for example in p → νK + due to the 12 ns kaon lifetime, is then beneficial and can be used to suppress background -as is done in SK analyses [74,[104][105][106].
Future detectors such as JUNO, DUNE, and HK are expected to improve on these important channels. While we are not aware of dedicated sensitivity studies, JUNO liquid-scintillator detector, which is analogous to Kam-LAND but with a 20 times larger size (20 kt), could be expected to test invisible-nucleon decays beyond 10 31 yr lifetime. Further, with low detection thresholds and sizable 40 kt mass, DUNE could also be very promising. To our knowledge no dedicated studies for de-excitation emission associated with nucleon decay from argon, as would be relevant for DUNE, have been performed thus far.
Despite being derived as limits on invisible nucleon decay, the discussed limits can be interpreted as the best current model-independent limits on the total nucleon lifetime and are likely applicable to all otherwise unconstrained modes in Tabs. I, II, and III. We note that there is still some implicit dependence on the actual nucleon decay mode. For example, decay modes with highly energetic pions can potentially destroy the daughter nucleus rather than just leaving a hole, but we envision that such events would be captured with other searches.
Finally, we comment that the invisible searches presented here could also be interpreted as tests of forbidden nuclear transitions, as was done in Ref. [100].

B. N → X+anything
Above we have discussed nucleon-decay searches that are nearly model-independent and only make use of the properties of the daughter nucleus that is created after nucleon decays within a nucleus. The current most sensitive searches in this category are from de-excitation emission, with resulting lifetime limits approaching 10 30 yr. This is more than four orders of magnitude below the best exclusive search limits (Tab. I) and one or two orders of magnitude below the predominant majority of other existing exclusive searches (e.g. Tabs. II and III).
In an effort to keep searches as model-independent as possible while utilizing the excellent particle identification of exclusive searches, we will now discuss inclusive nucleondecay searches, which correspond to N → X + anything with a light SM particle X of unknown energy. Without knowledge of the underlying energy distribution these are in general background-dominated searches.
1. N → (e ± , µ ± , π ± , K ± , ρ ± , K * ,± ) + anything In general, from size and exposure considerations, WC experiments are expected to dominate inclusive searches with a charged particle. However, the primary charged particle of interest could be invisible in WC experiments due to its Cherenkov threshold. Hence, additional detection efficiency can be expected in experiments with low-energy thresholds and where they are clearly visible, such as scintillators, especially for particles with a higher mass and a larger Cherenkov threshold.
From electric-charge conservation it is clear that any proton decay will eventually result in at least one positron, albeit potentially space-time-delayed if it originated from the decay of a heavier charged final-state particle (e.g. muon or pion). Hence, p → e ± + anything is an important general inclusive search. Since the positron will almost always have its momentum above the Cherenkov threshold, unless the number of final state particles is very high, a particularly promising sensitivity for this general inclusive search is achievable from a re-analysis of a 1-ring e-like data sample in SK and subsequently HK. The current best lifetime limit on the N → e + + anything search is 0.6 × 10 30 yr, from a 1979 analysis of Ref. [107]. We expect that SK and later HK could significantly improve on this result, due to their large detector volume and exposure. An estimate of such limit can be obtained from the recent SK search that placed a limit of 170 × 10 30 yr [78] on the p → e + νν channel, which bears close similarity to the inclusive p → e + + anything mode. Lowering the positron-momentum threshold employed in that analysis and varying the energy distribution, or even adapting a pure counting analysis, is expected to yield an inclusive N → e ± + anything mode limit of similar order with available SK data, i.e. O(100) × 10 30 yr.
In the case of inclusive N → µ ± + anything channel with a muon, if the muon has enough momentum it will lead to a µ-like ring in WC detectors. The current limit on this search is τ (p → µ + + anything) > 12 × 10 30 yr, converted from a 1981 Homestake experiment analysis of Ref. [60] that used a 0.3 kton WC detector. We again expect SK and future HK to improve on this limit. Similar to the inclusive search with an electron, the recent SK search of p → µ + νν that resulted in a lifetime limit of 220 × 10 30 yr [78] serves as a good indication of the potential limit on N → µ ± + anything. Inclusive searches for charged mesons are feasible as well, even though the mesons will eventually decay into charged leptons and thus in principle will be covered by the searches outlined above. Nevertheless, designated charged-meson searches can provide additional useful information. In WC detectors the charged pion π ± will look similar to a µ ± , albeit with a somewhat higher Cherenkov threshold, so the resulting limit can be expected to be of similar order. An estimate of p → π ± +anything can be obtained from the SK search for the two-body decay p → ν + π + that yielded a limit of 390 × 10 30 yr [73]. Unlike this two-body decay, multibody decays will result in a smeared pion momentumcausing an additional detection-efficiency penalty when the pion is below Cherenkov threshold.
Charged mesons more massive than the pion (i.e. K ± , ρ ± , K * ,± ) will always be below Cherenkov threshold when originating from decays of a single nucleon. Thus, in WC experiments such mesons can only be identified by their decay products. Due to their low momentum and short lifetime, they effectively decay at rest. There is thus not much dependence of the associated inclusive limits on momentum of the heavy meson and if the decay is two-body or multi-body. Hence, limits from exclusive searches with such mesons and other states being invisible, such as τ (p → ν +K + ) > 5900×10 30 yr [74] and τ (p → ν + K * ,+ ) > 130 × 10 30 yr [75] from SK as well as τ (p → ν + ρ + ) > 162 × 10 30 yr from IMB [62], can be readily re-interpreted as inclusive limits on p → (K ± , ρ ± , K * ,± )+anything processes.
We expect that SK can improve by more than an order over the existing IMB limits and even further improvement can be achieved with future HK. Further, JUNO and DUNE are particularly promising for complementary studies of these modes since they allow for an identification of the low-energy heavy charged mesons that are below Cherenkov threshold, with a particularly good efficiency achievable for kaons.
2. N → (π 0 , K 0 , η, ρ, ω, K * ,0 ) + anything Inclusive searches involving neutral mesons M 0 rely on the electromagnetically-interacting daughter particles of M 0 for detection. For example, in the rest frame π 0 → γγ decay results in photon lines of energy m π /2, clearly visible in WC detectors. For slow π 0 these are back-toback, otherwise more collimated. Already in the two-body nucleon decay the pion is not very energetic and carries energy below m N /2. Hence, the inclusive limit n → π 0 +anything is not particularly sensitive to pion momentum and one can estimate it to be of similar order as the limit from existing SK search for n → π 0 ν of 1100 × 10 30 yr [73].
As before, we expect that a re-analysis of the above modes with the SK data-set can improve IMB's limits by more than an order of magnitude and even better limits can be expected from HK.

N → γ + anything
Since we are envisioning multi-body nucleon decays, the emission of a photon will be typically suppressed compared to the same channel without a photon. It is still useful to perform such an analysis for final states that are otherwise difficult to detect, such as n → γ + neutrinos. As emphasized in Ref. [108], neutron decay into neutrinos corresponds to the sudden disappearance of the neutron's magnetic moment and should thus lead to electromagnetic radiation. From a particle-physics perspective this effect can be readily obtained from attaching a photon to the initial quarks in any diagram leading to n → neutrinos, in analogy to radiative τ decays [109]. Compared to low-energy de-excitation photons, these photons can reach energies of up to m N /2 and thus provide a far cleaner signature. However, with the probability that a photon emitted from this process has energy > 100 MeV being only ∼ 5 × 10 −5 [108], it is unlikely that such a search will be more sensitive than the n → invisible limit of Eq. (10). Further, the photon spectrum here depends on the details of the final state and number of emitted neutrinos, implying that obtained limits from such a search are less general than Eq. (10).
Nevertheless, this search would provide some complimentary to the low-energy γ searches from de-excitation and is also close to the inclusive search on e ± , as both result in an e-like ring in WC experiments.

N → ν + anything
As neutrinos escape the detector before interacting, this search is effectively invisible if nucleon decays are considered to occur within the experimental fiducial volume. However, it has been suggested that by considering such decays on macroscopic Earth-sized scales the additional cumulative neutrino flux from decays can be significant and observable. By attributing the experimentally observed ν µ /ν e fluxes to nucleon decays in Earth a limit of ∼ 10 26 yr was estimated [107]. However, a further reanalysis by the Fréjus experiment using 2 kton × yr of data yielded a lower limit of ∼ 10 25 [59] (see their Sec. 7 for flux estimation). Since SK has already collected over 370 kton × yr of data (e.g. [69]), we anticipate that an improved sensitivity of up to ∼ 10 28 yr could be achieved by re-analysis of this channel in SK using the full data set and even further improved upon in future HK. We note that exact search details will depend on the assumptions about the flavor of emitted neutrinos.
V. PROCESSES WITH ∆B > 1 So far we have focused on processes violating baryon number by one unit, ∆B = 1, which are described by the lowest-dimension ∆B = 0 SMEFT operators. However, nucleon decays also constitute sensitive probes of dimension d 6 operators beyond ∆B = 1. Such "multinucleon decays" with ∆B > 1 can be treated in complete analogy to the processes discussed in previous sections. The main difference is that ∆B > 1 processes have more available energy, which allows for the on-shell production of tau leptons [110] as well as heavier mesons such as D and φ, thus increasing the number of possible final states. Appearance of visible heavier charged mesons such as K ± , which are always below Cherenkov threshold when originating from single nucleon decays, also becomes possible. Only a few channels out of all possible kinematically allowed two-body ∆B = 2 di-nucleon decays have been experimentally studied thus far (17 out of 118), as we summarize in Tabs. IV, V, and VI. Even more limited are searches for ∆B > 2 processes, as we briefly comment on in Sec. V C.
A. ∆B = 2, ∆L = 0 As illustrated in Fig. 1, ∆B = 2 processes start at d = 9 and do not involve leptons. The corresponding operators, such as uddudd, have been discussed at length in the literature [113]. Experimentally, they have been probed through di-nucleon decays like nn → ππ as well as neutron-antineutron (n-n) oscillations [114]. Taking into account the possible quark-flavor structure of these operators it is clear that they can induce a variety of nn → M M processes, with M and M being mesons or photons, and similarly np, pp → mesons decays. While more complicated multi-meson final states are also expected, calculations of their branching ratios have not been extensively performed. Experimentally, only a small subset of such possible two-body final states has been searched for, namely nn → 2γ, N N → ππ and N N → KK, with the latter being motivated by    SUSY [115,116]. ∆B = 2 operators with ∆L = 0 involving leptons arise at d = 13 and should be suppressed compared to the d = 9 operators unless they carry lepton flavor [33,41], leading for example to nn → µ + e − or nn → µ + τ − . The tau modes have not been tested thus far, but since these channels have a fully visible final state they should allow for better sensitivity than the already performed np → τ + ν search in SK [27], which gave a limit of 29 × 10 30 yr.
The invisible di-nucleon decay nn → neutrinos can be searched in a way analogous to the ∆B = 1 case by looking for nuclear de-excitation emission [117]. The current best limits come from KamLAND [80], with a limit τ nn > 1.4 × 10 30 yr. Slightly lower limits around a few ×10 28 yr have been obtained by SNO+ [99] for np → invisible and pp → invisible modes. We envision that JUNO will be able to test these modes with lifetime sensitivity over an order of magnitude stronger than KamLAND, approaching 10 31 yr. Dedicated studies of sensitivity to these modes in DUNE, JUNO, and HK would provide valuable complementary information.

B. ∆B = 2, ∆L = ±2
At dimension d = 12, we find ∆B = ∆L = 2 operators, an example being uuuuddee. For first-generation quark-flavor indices, this simply induces the highly visible pp → + + , which has been discussed in the context of GUTs [118]. While the electron and muon modes have been searched for in SK [69], the one kinematically available tau mode pp → e + τ + has not. Similar to the case of nn → µ + τ − we expect that SK can probe its lifetime above 10 31 yr and better than the np → τ + ν search [27]. HK will allow to considerably improve on studies of these channels. We note that taking the quark indices of the d = 12 operators into account can lead to additional emission of kaons as well as other mesons in the final state, which we do not list here.
Still at dimension d = 12, we find ∆B = −∆L = 2 operators such as ddddddēē. Here, the dominant decay mode is nn → − − π + π + , potentially with kaons instead of pions. This competes with the loop-induced nn → νν, for which experimental limits exist. Other than this invisible di-neutron decay, no two-body decay fulfills ∆B = −∆L = 2, highlighting the importance of multi-body final states and inclusive searches.
Baryon number violation by more than two units has not been well explored, with sparse theoretical [119] and experimental [120][121][122][123] studies. As mentioned before, one particularly well-motivated case that falls into this category is ∆B = ∆L = 3, as mediated in the SM by non-perturbative electroweak instantons. At low temperatures and energies as in nuclear decays these processes are highly suppressed, while rate calculations in highenergy collider setups are notoriously difficult [124,125]. While we will not discuss the various ∆B > 2 final states and operators (starting at d = 15) further, we stress that inclusive searches can easily cover these processes as well, as we describe below. D. Inclusive searches for ∆B > 1 As discussed above, there is still significant untested parameter space for ∆B ≥ 2 processes, including about a hundred ∆B = 2 two-body final states with fixed kinematics (see Tabs. IV, V, and VI). While exclusive searches for these channels would clearly provide the best sensitivity, inclusive searches again offer an interesting alternative. A particularly relevant aspect of the proposed inclusive ∆B = 1 searches from Sec. IV B to stress is that they will automatically provide limits for ∆B > 1 channels. For example, an inclusive search for p → e + + anything with positron energy between ∼ MeV and m p /2 would also constrain np → e + + anything or nnp → e + + anything. In comparison, dedicated ∆B > 1 inclusive searches would allow to widen the search window to positron energies above m p /2 or look for on-shell tau leptons.

VI. CONCLUSIONS
Baryon number violation is strongly motivated from many distinct theoretical considerations, including GUTs, supersymmetry, and baryogenesis. Thus, probing proton decay and other baryon-number-violating processes is of fundamental importance in order to learn about physics beyond the SM. A slew of upcoming large-scale experiments, in particular JUNO, DUNE, and Hyper-Kamiokande, will be able to explore these processes with unparalleled sensitivity. Hence, it is crucial to make the most of these efforts (as well as existing detectors, most notably SK) and probe as much parameter space for baryon-number violation as possible.
As we highlighted in this work, nucleon decay constitutes a premier probe of baryon-number violation because ∆B = 1 operators will result in nucleon decay regardless of their quark or lepton flavor structure. While most nucleon decay searches have focused on two-body decay channels (e.g. p → e + π 0 ), as motivated by minimal GUT models and simplified EFT considerations, a vast landscape of possible baryon-number-violating processes exists beyond them.
Nucleon decay provides a unique opportunity to test baryon-number-violating operators of very high dimensionality and even beyond ∆B = 1, achieving sensitivity not available with other techniques. These operators can dominate over the conventionally considered lower-order dimension-six terms due to their lepton number and flavor quantum numbers. Going beyond the lower order, a vast amount of additional possible operators can start to contribute. Generically, such operators will result in complicated nucleon decays with a multi-body final state.
Most of the conducted nucleon decay searches have been exclusive, focusing on decay processes with a particular final state. We have identified a slew of single nucleon decay processes yet to be tested with exclusive studies. Further, exclusive searches for relatively clean channels like p → + − + or n → − K + π 0 could provide information indicative of d = 9 and 10 operators.
Exclusive studies rapidly become unfeasible to test nucleon decay beyond the simplest channels. As we argued in this work, a particularly important venue to further explore baryon-number violation in the future is through studies of inclusive nucleon decays N → X+anything (where X is a light SM particle with unknown energy distribution), as well as model-independent and invisible decay searches such as n → neutrinos (which can carry away an arbitrary lepton number). Such searches allow us to constrain multiple complicated nucleon decay channels simultaneously and broadly cover dimension d ≥ 9 baryon-number-violating operators. More-so, inclusive ∆B = 1 searches will allow us to gain some insight into ∆B > 1 processes as well.
For the searches outlined above, we foresee that the WC detectors SK and HK will in general provide the furthest reach in sensitivity due to their unparalleled size. However, we stress that searches involving heavier final state charged particles, e.g. kaons, can significantly benefit from studies in JUNO and DUNE. Further, the lowenergy sensitivity of JUNO and DUNE could be particularly beneficial for nuclear de-excitation emission searches (e.g. invisible decay channels). We are not aware of dedicated sensitivity studies related to these aspects. Finally, let us stress that even though we envision nuclear decays as the best probe for ∆B = 0, other baryon-numberviolating searches as in meson or tau decays could provide important complementary information.