Employing nucleon decay as a fingerprint of SUSY GUT models
using SusyTCProton

Stefan Antusch⋆111Email<EMAIL_ADDRESS>Christian Hohl⋆222Email:
<EMAIL_ADDRESS>and Vasja Susič⋆333Email<EMAIL_ADDRESS>
⋆ Department of Physics, University of Basel,

Klingelbergstr. 82, CH-4056 Basel, Switzerland

While the observation of nucleon decay would be a smoking gun of Grand Unified
Theories (GUTs) in general, the ratios between the decay rates of the various
channels carry rich information about the specific GUT model realization. To
investigate this fingerprint of GUT models in the context of supersymmetric
(SUSY) GUTs, we present the software tool SusyTCProton, which is an extension
of the module SusyTC to be used with the REAP package. It allows to calculate
nucleon decay rates from the relevant dimension five GUT operators specified
at the GUT scale, including the full loop-dressing at the SUSY scale. As an
application, we investigate the fingerprints of two example GUT toy models
with different flavor structures, performing an MCMC analysis to include the
experimental uncertainties for the charged fermion masses and CKM mixing
parameters. While both toy models provide equally good fits to the low energy
data, we show how they could be distinguished via their predictions of ratios
for nucleon decay rates. Together with SusyTCProton we also make the
additional module ProtonDecay public. It can be used independently from REAP
and allows to calculate nucleon decay rates from given $D=5$ and $D=6$
operator coefficients (accepting the required SUSY input for the $D=5$ case in
SLHA format). The $D=6$ functionality can also be used to calculate nucleon
decay in non-SUSY GUTs.

###### Contents

  1. 1 Introduction
  2. 2 Choice of two example models for investigating nucleon decay
    1. 2.1 General setup in flavor SUSY GUTs
    2. 2.2 Specific models: their Yukawa and quasi-Yukawa sectors
  3. 3 Considerations and procedure for a numerical analysis
    1. 3.1 Implementation of the Yukawa sector
    2. 3.2 Model parameters
    3. 3.3 Observables
    4. 3.4 Calculational procedure
  4. 4 Results for nucleon decay from SusyTCProton
    1. 4.1 Best-fit points and nucleon decay
    2. 4.2 MCMC results for nucleon decay
  5. 5 Conclusions
  6. A Details on doublet-triplet splitting
    1. A.1 The modified missing partner mechanism
      1. A.1.1 The missing partner mechanism — the usual setup
      2. A.1.2 The perturbativity problem and modifying the usual setup
    2. A.2 Double missing partner mechanism
  7. B Documentation for the proton decay calculation
    1. B.1 Numerical procedure
      1. B.1.1 Definition of the dimension 5 operators
      2. B.1.2 Running of the dimension 5 operators
      3. B.1.3 Dressing of the dimension 5 operators
      4. B.1.4 Running of the dimension 6 operators
      5. B.1.5 Partial decay widths
    2. B.2 Matching with SLHA conventions
  8. C User’s guide to Mathematica packages SusyTCProton and ProtonDecay
    1. C.1 SusyTCProton.m
    2. C.2 ProtonDecay.m

##  1 Introduction

Grand Unified Theories (GUTs) [1, 2, 3] continue to present an interesting
framework for physics Beyond the Standard Model (BSM). The attractive features
of this class of theories are the unification of all gauge interactions at
some high energy scale, (partial) unification of matter into representations
of a bigger symmetry, and an explanation for the quantization of hypercharge
in the Standard Model (SM).

The presence and type of low energy phenomena (e.g. light exotic matter)
typically depends on the GUT model and as such these features are not
universal predictions of gauge unification. An almost universal prediction of
GUT models, however, is proton decay — and more generally, nucleon decay. The
gauge bosons transforming as $(\mathbf{3},\mathbf{2},-5/6)$ under the Standard
Model gauge group
$G_{321}\equiv\mathrm{SU}(3)\times\mathrm{SU}(2)\times\mathrm{U}(1)$ present
in the adjoint representation of $\mathrm{SU}(5)$ violate baryon and lepton
number ($B$ and $L$) symmetries, and mediate the process of nucleon decay,
generating $B$ and $L$ violating four-fermion operators, i.e. $D=6$ proton
decay. Other particles, such as scalar leptoquarks, can of course also mediate
proton decay, but the main point is that as long as the unified group contains
$\mathrm{SU}(5)$ as a subgroup, proton (nucleon) decay will be present.

One advantage of nucleon decay as a probe of new physics is that our
experiments are able to measure decay rates corresponding to very large masses
of the mediators, and current experimental bounds demand that the GUT scale
$M_{\text{GUT}}$ must be higher than around $10^{15}\,\mathrm{GeV}$. Since
unification within field theory should happen below the Planck scale, this
implies a window between roughly
$10^{15}\,\mathrm{GeV}$–$10^{19}\,\mathrm{GeV}$ for unification. The
renormalization group (RG) running of gauge couplings in the SM can be made
consistent with this window in the presence of appropriate intermediate
states.

In supersymmetric (SUSY) theories, assuming the effective theory above the
SUSY scale $M_{\text{SUSY}}$ to be the Minimal Supersymmetric Standard Model
(MSSM), it is possible to construct $B$ and $L$ violating superpotential
operators of four fields of the type $qqql$, where $q$ stands for a quark and
$l$ for a lepton field. Such operators are generated, for example, when
integrating out chiral supermultiplets of heavy color triplets
$(\mathbf{3},\mathbf{1},-1/3)$ under $G_{321}$. The resulting Lagrangian
operators are dimension $5$ operators, consisting of $2$ fermions and $2$
scalar superpartners of SM fermions. The $D=5$ diagrams are then dressed at
the SUSY scale to have SM external states, thus generating four-fermion
nucleon decay operators. Contributions to nucleon decay generated in SUSY in
such a way are referred to as $D=5$ proton decay [4, 5], and they typically
dominate over $D=6$ contributions, also usually requiring a higher GUT scale
$M_{\text{GUT}}\gtrsim 10^{16}\,\mathrm{GeV}$ in SUSY GUTs compared to non-
SUSY GUTs (due to the requirement of sufficiently large triplet masses).

A rough estimate of the total proton decay rate can be obtained using the
formula

$\displaystyle\Gamma$ $\displaystyle\approx|\mathcal{A}|^{2}\,m_{p}^{5},$
(1.1)

where $m_{p}\approx 0.94\,\mathrm{GeV}$ is the proton (nucleon) mass and the
$D=5$ amplitude $A_{D=5}$ for the process can be estimated, see e.g. [6], as

$\displaystyle\mathcal{A}_{D=5}$
$\displaystyle\approx\frac{\alpha}{4\pi}\,y_{u}y_{d}\;\frac{m_{\tilde{w}}}{m_{T}m_{\tilde{f}}^{2}},$
(1.2)

with $\alpha=g^{2}/(4\pi)$, $g$ the unified gauge coupling, $y_{u}$ and
$y_{d}$ the Yukawa couplings of the up- and down-quarks, $m_{\tilde{w}}$ the
gaugino mass, $m_{\tilde{f}}$ the sfermion mass, and $m_{T}$ the mass of the
heavy color triplet mediating the decay. We assumed that the $D=6$
contribution, e.g. from gauge bosons, is negligible, i.e.
$|\mathcal{A}_{D=5}|\gg|\mathcal{A}_{D=6}|$, since its amplitude is
proportional to the square of the inverse heavy mass:

$\displaystyle\mathcal{A}_{D=6}$
$\displaystyle\approx\frac{4\pi\alpha}{m_{X}^{2}},$ (1.3)

with $m_{X}$ the mass of the gauge bosons mediating the $B$ and $L$ violating
interaction. A thorough analysis of nucleon decay, however, should determine
the decay rates into each channel separately, and both $D=5$ and $D=6$
contributions should be taken into account (we assume $R$-parity conservation,
so no $D=4$ proton decay).

It is well known, see e.g. [7], that in SUSY theories the dominant decay
channel is usually the kaon-neutrino channel $p\to K^{+}\bar{\nu}$. The proton
lifetime bound for this channel currently comes from Super-Kamiokande [8]:
$6.6\cdot 10^{33}$ years at $90\,\%$ C.L. Future experiments are set to
improve on this bound further: the expected $90\,\%$ C.L. bounds after 10
years of operation are $3.3\cdot 10^{34}$ years for Dune [9], $3.0\cdot
10^{34}$ years for Hyper-Kamiokande [10] and $2.0\cdot 10^{34}$ years for JUNO
[11].

To investigate the viability of a specific SUSY GUT model and determine its
predictions, at least as far as nucleon decay is concerned, a more precise
estimate of proton/nucleon decay rates is necessary than that of Eq. (1.2).
While the procedure and expressions for the decay rates of the various nucleon
decay channels are known, cf. [12, 7], they are cumbersome to implement from
scratch. A software package capable of the computation of nucleon decay rates
in the various decay channels would thus be most welcome. Ideally, this tool
would be able to handle all the complications in the computation that can
arise, in particular both the sparticle spectrum and the flavor violating
effects in the SUSY sector can play a prominent role. Furthermore, the
effective $D=5$ operators arising from integrating out the heavy triplets may
also carry a flavor imprint from the Yukawa sector, since both the couplings
to the heavy triplets and MSSM Higgs doublets arise from the same GUT
operators. To be able to handle all of this automatically, the tool would thus
ideally be built on top of an existing SUSY spectrum generator of sufficient
sophistication.

In this paper we roll out exactly such a nucleon decay calculation tool for
SUSY GUTs consisting of two in principle independent Mathematica modules
SusyTCProton.m and ProtonDecay.m that we make publicly available:444A User’s
Guide for the installation and use of these packages is provided in Appendix
C.

  1. 1.

The SusyTCProton package is a nucleon decay extension of the existing SusyTC
package [13] built within the framework of REAP [14]. The standard SusyTC
performs the RG evolution of the softly broken MSSM from the GUT scale to the
$Z$-boson scale, automatically matching the MSSM and SM at an appropriate
scale and generating all the SUSY sector information; the nucleon decay
extension of this package performs, in addition, the $1$-loop running of $D=5$
operators from the GUT to the SUSY scale.

  2. 2.

The ProtonDecay package takes as input the $D=5$ decay operators and the
sparticle masses and mixings, all at the SUSY scale, and dresses the relevant
diagrams with SUSY sparticles to determine the effective four-fermion
operators of nucleon decay, and after RG running to the nucleon scale computes
the decay rates in $8$ proton decay channels and $5$ neutron decay channels
(not counting different neutrino flavors, see Table 2 for the list). The input
for this package is compatible both with a SusyTCProton output of internal
values, as well as with an SLHA file [15] amended with the $D=5$ operator
information. This allows ProtonDecay to be used also with any spectrum
generator that can output an SLHA file, such as SPheno [16, 17], SOFTSUSY
[18], SuSeFLAV [19] or SuSpect [20], provided that the user inputs $D=5$
operator values from their own source. In addition, the ProtonDecay package
can handle $D=6$ contributions to nucleon decay, either as standalone terms
(the non-SUSY GUT case) or in tandem with $D=5$ contributions; the user needs
to input them at a scale below which only the SM RG running is applied.

As a demonstration, we employ the above software packages in this paper to
analyze the nucleon decay in two example toy models, cf. [21, 22], both of
which are of the flavor SUSY GUT variety. These models fit the SM fermion
masses and mixings through particular choices of Yukawa textures and GUT-
operators, which generate the Yukawa sector. The type of GUT operators used
implies a choice of the Clebsch-Gordan (CG) coefficients between different
sectors associated by the GUT symmetry. These same features impact also the
$D=5$ proton decay predictions, and it is precisely these flavor effects on
nucleon decay that we focus on here. The magnitude of the decay rates in the
various channels is subject to various potential ambiguities, such as the mass
of the heavy triplet. These are overcome, at least in the models we consider
here, by investigating ratios of decay rates, allowing to identify a nucleon
decay “fingerprint” for each model. If nucleon decay events are seen
experimentally in at least one decay channel, such a fingerprint technique, as
we demonstrate in this paper, can be a powerful tool to either discriminate
within a set of models with equally good Yukawa-sector fits and near-identical
SUSY spectra, or to possibly rule out particular models of this type.

The paper is organized as follows. We define the example flavor GUT models in
Section 2, discuss how we implemented our numerical analysis in Section 3,
present our results regarding nuclear decay in these models in Section 4 and
give our conclusions in Section 5. The paper contains also 3 appendices.
Appendix A presents a discussion on our doublet-triplet splitting setup and a
few of its concrete realizations. Appendix B provides technical details on the
nucleon decay computation and a note on conventions. Finally, Appendix C
contains the documentation for our software packages SusyTCProton and
ProtonDecay.

##  2 Choice of two example models for investigating nucleon decay

###  2.1 General setup in flavor SUSY GUTs

The models under consideration are flavor $\mathrm{SU}(5)$ SUSY GUTs of the
type discussed in [21, 22]. In this paper, we consider $2$ toy models of this
type and analyze their nucleon decay predictions.

The 3 generations of SM fermions are embedded in the standard way into
$\mathrm{SU}(5)$ representations, where the index $i$ runs from $1$ to $3$:

$\displaystyle\mathbf{10}_{Fi}$ $\displaystyle=Q_{i}\oplus u^{c}_{i}\oplus
e^{c}_{i},$ (2.1) $\displaystyle\mathbf{\bar{5}}_{Fi}$
$\displaystyle=L_{i}\oplus d^{c}_{i}.$ (2.2)

The superpotential of the models can be separated into the Yukawa and breaking
part

$\displaystyle W$ $\displaystyle=W_{\text{Yuk}}+W_{\text{breaking}},$ (2.3)

where the breaking part consists of those terms, by definition, which do not
contain the “fermionic” chiral representations $\mathbf{10}_{Fi}$ and
$\mathbf{\bar{5}}_{Fi}$. For each model we assume a particular texture of
Yukawa matrices, to be specified later, and single operator dominance in each
Yukawa entry, i.e. only a single GUT operator is present or dominates every
Yukawa entry (see [23, 24, 25, 22, 26] for concrete realizations in the
literature). This imposes strict relations between different fermion sectors
and increases model predictivity. The Yukawa sector of each model thus
consists of a set of (non-renormalizable) operators.

We omit the high-energy model building in the neutrino sector. Also, we remain
agnostic about the details of the breaking sector; the only requirement is
that it breaks $\mathrm{SU}(5)$ symmetry to the Standard Model group
$G_{321}:=\mathrm{SU}(3)\times\mathrm{SU}(2)\times\mathrm{U}(1)$ and includes
the Higgs chiral supermultiplets assumed to be present due to their use in the
Yukawa sector.

We are interested in models with as predictive $D=5$ proton decay as possible.
Since this mode of proton decay is mediated by color triplets
$(\mathbf{3},\mathbf{1},-1/3)$ of $G_{321}$, the details of doublet-triplet
(DT) splitting come into the decay prediction. We thus consider as simple a
scenario as possible for DT splitting, where only one triplet pair
$T_{1}\oplus\overline{T}_{1}$ couples to the fermionic sector, and for
simplicity we assume that it is coming from
$\mathbf{5}_{H}\oplus\mathbf{\bar{5}}_{H}$ of $\mathrm{SU}(5)$.

This puts a limitation on which flavor GUT models can be considered in the
context of this paper: the Higgs cannot be contained in the $\mathbf{45}$ or
$\mathbf{\overline{45}}$ of $\mathrm{SU}(5)$, and the non-renormalizable
operators making use of such representations in the lists provided in [27, 28]
are not available, but only those containing $\mathbf{5}_{H}$ or
$\mathbf{\bar{5}}_{H}$.

At the GUT scale $M_{\text{GUT}}$, the GUT models are matched to the softly
broken MSSM as the effective theory, with the addition of the effective $D=5$
proton decay operators in the superpotential. For the soft terms, we assume
constrained MSSM (cMSSM, a.k.a. mSUGRA) boundary conditions at
$M_{\text{GUT}}$ (see e.g. [29]), so they are parametrized by

$\displaystyle m_{0},\quad M_{1/2},\quad A_{0},$ (2.4)

where $m_{0}$ is the universal soft scalar mass, $M_{1/2}$ the universal
gaugino mass and $A_{0}$ the universal $A$-term factor. The choice of cMSSM is
in principle not necessary, but it ensures that our analysis isolates the
effects on nucleon decay from the flavor structure of the theory without added
flavor ambiguities from the SUSY sector.

At the level of the MSSM, the Yukawa sector (including the relevant operators
for proton decay, but ignoring any neutrino terms) consists of

$\displaystyle W_{\text{Yuk}}$
$\displaystyle=\sum_{i,j=1}^{3}\;\;(\mathbf{Y}_{u})_{ij}\;Q_{i}\cdot
H_{u}\,u_{j}^{c}-(\mathbf{Y}_{d})_{ij}\;Q_{i}\cdot
H_{d}\,d^{c}_{j}-(\mathbf{Y}_{e})_{ij}\;L_{i}\cdot H_{d}\,e^{c}_{j}+W_{T},$
(2.5)

where the color indices have been suppressed and the $\mathrm{SU}(2)$
contraction is written as $\Psi\cdot\Phi=\epsilon_{kl}\Psi^{k}\Phi^{l}$, with
$\epsilon_{kl}$ the $2$-index completely anti-symmetric tensor in the
convention $\epsilon_{12}=+1$. Notice that we have written the Yukawa terms in
the left-right (LR) convention.

The triplet part of the Yukawa superpotential consists of

$\displaystyle\begin{split}W_{\text{T}}&=+\tfrac{1}{2}(\tilde{\mathbf{Y}}_{qq})_{ij}\,\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}\,Q_{i}^{\hat{a}}\cdot
Q_{j}^{\hat{b}}\,T_{1}^{\hat{c}}+(\tilde{\mathbf{Y}}_{eu})_{ij}\,e^{c}_{i}\,u^{c}_{j{\hat{a}}}\,T_{1}^{\hat{a}}\\\
&\quad+(\tilde{\mathbf{Y}}_{ql})_{ij}\,Q_{i}^{\hat{a}}\cdot
L_{j}\,\overline{T}_{1{\hat{a}}}+(\tilde{\mathbf{Y}}_{ud})_{ij}\,\epsilon^{{\hat{a}}{\hat{b}}{\hat{c}}}\,u^{c}_{i{\hat{a}}}\,d^{c}_{j{\hat{b}}}\,\overline{T}_{1{\hat{c}}}\\\
&\quad+\sum_{I,J}(\mathbf{M}_{T})_{IJ}\,T_{I}^{\hat{a}}\,\overline{T}_{J{\hat{a}}},\end{split}$
(2.6)

where ${\hat{a}},{\hat{b}},{\hat{c}}$ are $\mathrm{SU}(3)$ indices, while
$\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}$ and
$\epsilon^{{\hat{a}}{\hat{b}}{\hat{c}}}$ are Levi-Civita tensors with
$\epsilon_{123}=\epsilon^{123}=1$. The $\mathrm{SU}(2)$ contractions are again
labeled with a $\cdot$ dot. We refer to the $3\times 3$ matrices
$\tilde{\mathbf{Y}}_{qq}$, $\tilde{\mathbf{Y}}_{eu}$,
$\tilde{\mathbf{Y}}_{ql}$ and $\tilde{\mathbf{Y}}_{ud}$ as quasi-Yukawa
matrices in the remainder of this text.

This ansatz is specific to the models under consideration, since only the
(anti-)triplets of $\mathbf{5}_{H}\oplus\mathbf{\bar{5}}_{H}$ couple to the SM
fermions; we thus omit any triplet indices $I,J$ in the quasi-Yukawas. The
general case with an arbitrary number of triplets is considered in Appendix
B.1.1; the quasi-Yukawas written above would have $I,J=1$, while all other
quasi-Yukawas for $I,J\neq 1$ would be zero.

Integrating out the triplets and adapting the result from Appendix B.1.1 to
our specific case, we obtain the following effective $D=5$ proton decay
operators:

$\displaystyle W_{T}$
$\displaystyle\to-\tfrac{1}{2}\mathbf{C}^{ijkl}_{5L}\,\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}\,(Q_{i}^{\hat{a}}\cdot
L_{j})(Q_{k}^{\hat{b}}\cdot
Q_{l}^{\hat{c}})-\mathbf{C}^{ijkl}_{5R}\,\epsilon^{{\hat{a}}{\hat{b}}{\hat{c}}}\,u^{c}_{i{\hat{a}}}\,d^{c}_{j{\hat{b}}}\,e^{c}_{k}\,u^{c}_{l{\hat{c}}}\,,$
(2.7)

where the $\mathbf{C}_{5}$ coefficients for our models equal

$\displaystyle\mathbf{C}^{ijkl}_{5L}$
$\displaystyle=(1/M^{\text{eff}}_{T})\,(\tilde{\mathbf{Y}}_{ql})_{ij}\,(\tilde{\mathbf{Y}}_{qq})_{kl},$
(2.8) $\displaystyle\mathbf{C}^{ijkl}_{5R}$
$\displaystyle=(1/M^{\text{eff}}_{T})\,(\tilde{\mathbf{Y}}_{ud})_{ij}\,(\tilde{\mathbf{Y}}_{eu})_{kl},$
(2.9)

with

$\displaystyle M^{\text{eff}}_{T}:=((\mathbf{M}_{T}^{-1})_{11})^{-1}$ (2.10)

denoting the “triplet effective mass”. We provide a more detailed discussion
on the form of the $2\times 2$ matrix $\mathbf{M}_{T}$ and the effective
triplet mass in Appendix A. It is clear, however, that in this scheme the
effective triplet mass does not necessarily correspond to any physical scale
of particles, and as such can have a value above the Planck scale without
losing consistency.

In the $\mathrm{SU}(5)$ context, the Yukawa and quasi-Yukawa matrices come
from the following GUT operators:

$\displaystyle\mathbf{10}_{F}\,\mathbf{10}_{F}\,\mathbf{X}$
$\displaystyle\quad\supset\quad\mathbf{Y}_{u},\tilde{\mathbf{Y}}_{qq},\tilde{\mathbf{Y}}_{eu},$
(2.11) $\displaystyle\mathbf{10}_{F}\,\mathbf{\bar{5}}_{F}\,\mathbf{X}$
$\displaystyle\quad\supset\quad\mathbf{Y}_{d},\mathbf{Y}_{e},\tilde{\mathbf{Y}}_{ql},\tilde{\mathbf{Y}}_{ud},$
(2.12)

where $\mathbf{X}$ is a placeholder for one or multiple Higgs fields in a
Yukawa operator, possibly different for different matrix entries. The matrices
arising from the same operator are related via $\mathrm{SU}(5)$ CG
coefficients, which depend on the combination of Higgs fields $\mathbf{X}$
(entry by entry).

###  2.2 Specific models: their Yukawa and quasi-Yukawa sectors

We consider two models in this paper, referred to as model 1 and model 2. They
in fact make use of the Yukawa textures of the models constructed in section 4
of [21]. We are interested only in the Yukawa sector here, so we remain
agnostic about the details of doublet-triplet splitting (see Appendix A for
some possibilities).

The Yukawa textures of the two models are:

$\displaystyle\text{model 1}:$ $\displaystyle\mathbf{Y}_{\mathbf{10}}$
$\displaystyle=\begin{pmatrix}\ast&\ast&0\\\ \ast&\ast&\ast\\\ 0&\ast&\ast\\\
\end{pmatrix},$ $\displaystyle\mathbf{Y}_{\mathbf{\bar{5}}}$
$\displaystyle=\begin{pmatrix}0&\ast&0\\\ \ast&\ast&0\\\ 0&0&\ast\\\
\end{pmatrix},$ (2.13) $\displaystyle\text{model 2}:$
$\displaystyle\mathbf{Y}_{\mathbf{10}}$
$\displaystyle=\begin{pmatrix}\ast&\ast&\ast\\\ \ast&\ast&\ast\\\
\ast&\ast&\ast\\\ \end{pmatrix},$ $\displaystyle\mathbf{Y}_{\mathbf{\bar{5}}}$
$\displaystyle=\begin{pmatrix}\ast&0&0\\\ 0&\ast&0\\\ 0&0&\ast\\\
\end{pmatrix},$ (2.14)

where $\ast$ denotes a non-vanishing complex entry and
$\mathbf{Y}_{\mathbf{10}}$ are always symmetric:
$\mathbf{Y}_{\mathbf{10}}=\mathbf{Y}_{\mathbf{10}}^{\textsf{T}}$. Note that in
order to consider the two models on the same footing, we do not assume any
spontaneous CP violation mechanism555The use of a spontaneous CP violation
mechanism [30, 31] in model 1 could force the
$(\mathbf{Y}_{\mathbf{10}})_{12}$ entry imaginary and the other Yuakwa
parameters real [32, 24]; this would automatically give a correct prediction
for the CKM CP-violating phase. of model 1, so the parameters are taken as
complex.

The non-zero entries in these matrices are populated by the operators
specified in Table 1 (chosen from Table 1 of [27]). The $\mathrm{SU}(5)$
representation in the index position indicates the type of contraction used
between the pieces (analogous to integrating out the mediator in the index
from a renormalizable interaction of that mediator with the fields in the
parentheses). The vacuum expectation value (VEV) of the $\mathbf{24}_{H}$ is
at the GUT scale.

Table 1: Table of operators used in each Yukawa entry for models 1 and 2. entry | model 1 operator | model 2 operator
---|---|---
$(\mathbf{Y}_{\mathbf{10}})_{ij}$ |  $\mathbf{10}_{Fi}\mathbf{10}_{Fj}\mathbf{5}_{H}$ | $\mathbf{10}_{Fi}\mathbf{10}_{Fj}\mathbf{5}_{H}$
$(\mathbf{Y}_{\mathbf{\bar{5}}})_{11}$ | / | $(\mathbf{10}_{F{1}})_{\mathbf{\overline{45}\otimes\mathbf{\overline{45}}}}\;(\mathbf{\bar{5}}_{H}\mathbf{24}_{H})_{\mathbf{45}}\;(\mathbf{\bar{5}}_{F{1}}\mathbf{24}_{H})_{\mathbf{45}}$
$(\mathbf{Y}_{\mathbf{\bar{5}}})_{12}$ |  $(\mathbf{10}_{F{1}})_{\mathbf{\overline{45}\otimes\mathbf{\overline{45}}}}\;(\mathbf{\bar{5}}_{H}\mathbf{24}_{H})_{\mathbf{45}}\;(\mathbf{\bar{5}}_{F{2}}\mathbf{24}_{H})_{\mathbf{45}}$ | /
$(\mathbf{Y}_{\mathbf{\bar{5}}})_{21}$ |  $(\mathbf{10}_{F2}\mathbf{24}_{H})_{\mathbf{\overline{10}}}\;(\mathbf{\bar{5}}_{F1}\mathbf{\bar{5}}_{H})_{\mathbf{10}}$ | /
$(\mathbf{Y}_{\mathbf{\bar{5}}})_{22}$ |  $(\mathbf{10}_{F2}\mathbf{24}_{H})_{\mathbf{\overline{10}}}\;(\mathbf{\bar{5}}_{F2}\mathbf{\bar{5}}_{H})_{\mathbf{10}}$ | $(\mathbf{10}_{F2}\mathbf{24}_{H})_{\mathbf{\overline{10}}}\;(\mathbf{\bar{5}}_{F2}\mathbf{\bar{5}}_{H})_{\mathbf{10}}$
$(\mathbf{Y}_{\mathbf{\bar{5}}})_{33}$ |  $(\mathbf{10}_{F{3}}\mathbf{\bar{5}}_{H})_{\mathbf{\bar{5}}}\;(\mathbf{10}_{F{3}}\mathbf{24}_{H})_{\mathbf{5}}$ | $(\mathbf{10}_{F{3}}\mathbf{\bar{5}}_{H})_{\mathbf{\bar{5}}}\;(\mathbf{10}_{F{3}}\mathbf{24}_{H})_{\mathbf{5}}$

This set of operators gives the following CG factor relations for the Yukawas
and quasi-Yukawas, which we computed explicitly (the reader can also cross-
check with [27]):

  * •

model 1:

$\displaystyle\mathbf{Y}_{e}$
$\displaystyle=\begin{pmatrix}\phantom{+}0&\phantom{+}6&\phantom{+}0\\\
-\tfrac{1}{2}&\phantom{+}6&\phantom{+}0\\\
\phantom{+}0&\phantom{+}0&-\tfrac{3}{2}\end{pmatrix}\cdot\mathbf{Y}_{d}^{\textsf{T}}\,,$
$\displaystyle\tilde{\mathbf{Y}}_{ql}$
$\displaystyle=\begin{pmatrix}\phantom{+}0&-1&\phantom{+}0\\\
-1&-1&\phantom{+}0\\\
\phantom{+}0&\phantom{+}0&\phantom{+}\tfrac{3}{2}\end{pmatrix}\cdot\mathbf{Y}_{d}\,,$
$\displaystyle\tilde{\mathbf{Y}}_{ud}$
$\displaystyle=\begin{pmatrix}\phantom{+}0&\phantom{+}\tfrac{2}{3}&\phantom{+}0\\\
-4&-4&\phantom{+}0\\\
\phantom{+}0&\phantom{+}0&\phantom{+}1\end{pmatrix}\cdot\mathbf{Y}_{d}\,,$
(2.15) $\displaystyle\tilde{\mathbf{Y}}_{qq}$
$\displaystyle=-\mathbf{Y}_{u}\,,$ $\displaystyle\tilde{\mathbf{Y}}_{eu}$
$\displaystyle=+\mathbf{Y}_{u}\,.$ (2.16)

  * •

model 2:

$\displaystyle\mathbf{Y}_{e}$
$\displaystyle=\begin{pmatrix}-\tfrac{1}{2}&\phantom{+}0&\phantom{+}0\\\
\phantom{+}0&\phantom{+}6&\phantom{+}0\\\
\phantom{+}0&\phantom{+}0&-\tfrac{3}{2}\end{pmatrix}\cdot\mathbf{Y}_{d}^{\textsf{T}}\,,$
$\displaystyle\tilde{\mathbf{Y}}_{ql}$
$\displaystyle=\begin{pmatrix}-1&\phantom{+}0&\phantom{+}0\\\
\phantom{+}0&-1&\phantom{+}0\\\
\phantom{+}0&\phantom{+}0&\phantom{+}\tfrac{3}{2}\end{pmatrix}\cdot\mathbf{Y}_{d}\,,$
$\displaystyle\tilde{\mathbf{Y}}_{ud}$
$\displaystyle=\begin{pmatrix}\phantom{+}\tfrac{2}{3}&\phantom{+}0&\phantom{+}0\\\
\phantom{+}0&-4&\phantom{+}0\\\
\phantom{+}0&\phantom{+}0&\phantom{+}1\end{pmatrix}\cdot\mathbf{Y}_{d}\,,$
(2.17) $\displaystyle\tilde{\mathbf{Y}}_{qq}$
$\displaystyle=-\mathbf{Y}_{u}\,,$ $\displaystyle\tilde{\mathbf{Y}}_{eu}$
$\displaystyle=+\mathbf{Y}_{u}\,.$ (2.18)

The dot $\cdot$ indicates entry-wise multiplication of matrices.

##  3 Considerations and procedure for a numerical analysis

###  3.1 Implementation of the Yukawa sector

Following the considerations in Section 2, we implement models 1 and 2 at the
GUT scale in the framework of the MSSM. The Yukawa sector is determined by
specifying the two GUT Yukawa matrices $\mathbf{Y}_{\mathbf{\bar{5}}}$ and
$\mathbf{Y}_{\mathbf{10}}$, which is equivalent to specifying the up- and
down-quark matrices $\mathbf{Y}_{u}$ and $\mathbf{Y}_{d}$, since all the other
Yukawa and quasi-Yukawa parameters are related to them via CG factors, see Eq.
(2.15)-(2.16) and (2.17)-(2.18) for model 1 and 2, respectively.

Given the textures of the two models, $\mathbf{Y}_{u}$ is a complex symmetric
matrix, with the added constraint
$(\mathbf{Y}_{u})_{13}=(\mathbf{Y}_{u})_{31}=0$ on its entries in the case of
model 1. Any symmetric complex $3\times 3$ matrix $\mathbf{Y}$ can be
decomposed, using Takagi decomposition, as

$\displaystyle\mathbf{Y}=\mathbf{U}\;\operatorname{diag}(y_{1},y_{2},y_{3})\;\mathbf{U}^{\textsf{T}},$
(3.1)

where the singular values $y_{i}\geq 0$, and $\mathbf{U}$ is a $3\times 3$
unitary matrix. For the parametrization of $\mathbf{U}$, we shall make use of
the definitions

$\displaystyle\mathbf{U}_{12}(\theta,\delta)$
$\displaystyle:=\begin{pmatrix}\phantom{-}\cos\theta\,\phantom{e^{i\delta}}&\sin\theta\,e^{-i\delta}&0\\\
-\sin\theta\,e^{i\delta}&\cos\theta&0\\\ 0&0&1\end{pmatrix},$ (3.2)
$\displaystyle\mathbf{U}_{23}(\theta,\delta)$
$\displaystyle:=\begin{pmatrix}1&0&0\\\
0&\phantom{+}\cos\theta\,\phantom{e^{i\delta}}&\sin\theta\,e^{-i\delta}\\\
0&-\sin\theta\,e^{i\delta}&\cos\theta\end{pmatrix},$ (3.3)
$\displaystyle\mathbf{U}_{13}(\theta,\delta)$
$\displaystyle:=\begin{pmatrix}\phantom{-}\cos\theta\,\phantom{e^{i\delta}}&0&\sin\theta\,e^{-i\delta}\\\
0&1&0\\\ -\sin\theta\,e^{i\delta}&0&\cos\theta\\\ \end{pmatrix},$ (3.4)

where $\mathbf{U}_{12}$, $\mathbf{U}_{23}$ and $\mathbf{U}_{13}$ represent
unitary matrices of $2$-D complex rotations, each in principle depending upon
an angle and a phase. The most general $3\times 3$ unitary matrix can be
explicitly written by using the following family of real parameterizations:

$\displaystyle\mathbf{U}$
$\displaystyle=\operatorname{diag}(e^{i\eta_{1}},e^{i\eta_{2}},e^{i\eta_{3}})\;\mathbf{U}_{23}(\theta^{uL}_{23},\delta_{0})\;\mathbf{U}_{13}(\theta^{uL}_{13},\delta)\;\mathbf{U}_{12}(\theta^{uL}_{12},0)\;\operatorname{diag}(e^{i\phi_{1}/2},e^{i\phi_{2}/2},1),$
(3.5)

where the family parameter $\delta_{0}$ is fixed, and the $9$ free parameters
describing the unitary rotation consist of angles $\theta_{ij}$ and phases
$\delta$, $\phi_{i}$ and $\eta_{i}$:

$\displaystyle\theta_{12},\quad\theta_{13},\quad\theta_{23},\quad\delta,\quad\phi_{1},\quad\phi_{2},\quad\eta_{1},\quad\eta_{2},\quad\eta_{3}.$
(3.6)

The first and last factor of the product in Eq. (3.5) are the two diagonal
matrices containing the “outer phases” $\eta_{i}$ and the “inner phases”
$\phi_{i}$ (with inner/outer referring to their location in the product of Eq.
3.1). There are only $2$ $\phi$-phases, since a third such phase can be
absorbed as a common phase factor into the $\eta$ phases. Furthermore, the
$1/2$ factors for the $\phi$-phases are there for later convenience; we shall
view them as phases of the singular values in the Takagi decomposition rather
than as part of the unitary rotation.

Due to greater convenience, we use different parameterizations in the two
models. In model 1, we set $\delta_{0}=\pi/2$ and impose the additional
(complex) condition $(\mathbf{Y}_{u})_{13}=0$ by solving for the real
variables $\theta_{13}$ and $\delta$. For model 2 we set $\delta_{0}=0$. The
advantages of this choice of $\delta_{0}$ will be discussed later for each
model separately.

The use of Takagi decomposition for the parametrization of $\mathbf{Y}_{u}$
has the advantage that it provides direct control over the singular values
(Yukawa couplings) and mixing angles in the up-quark Yukawa matrix. Also, note
that we can redefine each fermionic representation $\mathbf{10}_{Fi}$ and
$\bar{\mathbf{5}}_{Fi}$ with an arbitrary complex phase, implying that $6$
phases in the Yukawa sector are unphysical. In our convenient parameterization
of $\mathbf{Y}_{u}$, a phase redefinition of $\mathbf{10}_{Fi}$ allows us to
fix the outer phases without loss of generality:

$\displaystyle\eta_{1}=\eta_{2}=\eta_{3}$ $\displaystyle=0.$ (3.7)

We shall use the phase redefinitions of $\bar{\mathbf{5}}_{Fi}$ to remove the
phase of one complex entry in each column of $\mathbf{Y}_{d}$. In model 2,
this implies that $\mathbf{Y}_{d}$ is diagonal with (positive) real values,
while in model 1 a complex phase remains in one $2$nd column entry; we choose
that entry to be $(\mathbf{Y}_{d})_{12}$.

Putting all the above considerations together, the final form of the
parameterizations of the Yukawa matrices $\mathbf{Y}_{u}$ and
$\mathbf{Y}_{d}$, which uniquely determine the entire Yukawa sector, are the
following:

  * •

Model 1:
According to the textures in Eq. (2.13), the parametrization of the down- and
up-quark Yukawa matrix in the context of the MSSM at the GUT scale is given by

$\displaystyle\mathbf{Y}_{d}$
$\displaystyle=\begin{pmatrix}0&y^{d}_{12}\,e^{i\varphi}&0\\\
y^{d}_{21}&y^{d}_{22}&0\\\ 0&0&y^{d}_{3}\end{pmatrix},$ (3.8)
$\displaystyle\mathbf{Y}_{u}$
$\displaystyle=\mathbf{U}_{1}\;\text{diag}(y^{u}_{1}\,e^{i\phi_{1}},y^{u}_{2}\,e^{i\phi_{2}},y^{u}_{3})\;\mathbf{U}_{1}^{\textsf{T}},$
(3.9)

with

$\displaystyle\mathbf{U}_{1}$
$\displaystyle=\mathbf{U}_{23}(\theta^{uL}_{23},\pi/2)\;\mathbf{U}_{13}(\theta^{uL}_{13},\delta)\;\mathbf{U}_{12}(\theta^{uL}_{12},0),$
(3.10)

where the values of $\theta^{uL}_{13}$ and $\delta$ depend on the other
parameters of the up-sector; they are determined from demanding
$(\mathbf{Y}_{u})_{13}=0$:666We omit an overall $\cos\theta^{uL}_{13}$ factor
in the $(\mathbf{Y}_{u})_{13}=0$ condition. Such a factor allows for an
alternative solution $\theta^{uL}_{13}=\pi/2$, which also sets
$(\mathbf{Y}_{u})_{23}=0$, so we reject it.

$\displaystyle(C_{1}e^{i\delta}+C_{2}e^{-i\delta})\sin\theta^{uL}_{13}+C_{3}$
$\displaystyle=0,$ (3.11)

with

$\displaystyle\begin{split}C_{1}&:=y^{u}_{3}\,\cos\theta^{uL}_{23},\\\
C_{2}&:=(y^{u}_{2}e^{i\phi_{2}}\sin^{2}\theta^{uL}_{12}-y^{u}_{1}e^{i\phi_{1}}\cos^{2}\theta^{uL}_{12})\,\cos\theta^{uL}_{23},\\\
C_{3}&:=i\,(y^{u}_{1}e^{i\phi_{1}}+y^{u}_{2}e^{i\phi_{2}})\,\sin\theta^{uL}_{12}\,\cos\theta^{uL}_{12}\,\sin\theta^{uL}_{23}.\\\
\end{split}$ (3.12)

Once the numeric values of the free parameters are specified, the $C$
coefficients in Eq. (3.12) become complex numbers, allowing for the complex
equation (3.11) to be solved numerically. Given the hierarchical structure of
masses and mixings in the quark sector, this equation can always be solved for
$\delta$ and $\theta^{uL}_{13}$ in a numerically stable fashion. The values
for the angle are very small: $\theta^{13}\lesssim 10^{-5}$, which suggests a
possible alternative which we mention here in passing: we could omit solving
Eq. (3.11) and simply take $\theta^{uL}_{13}=0$ ($\delta$ dependence then
disappears). This would be a relatively good approximation yielding
$|(\mathbf{Y}_{u})_{13}|\lesssim 10^{-5}\approx 0$.

The CKM matrix comes from diagonalizing both $\mathbf{Y}_{u}$ and
$\mathbf{Y}_{d}$. It turns out that for this particular texture, choosing the
$\delta_{0}=\pi/2$ parametrization for the up-sector then already provides a
good CKM phase for $\varphi\approx 0$ in $\mathbf{Y}_{d}$.

Also, note that the complex phase $\phi$ here is in the
$(\mathbf{Y}_{d})_{12}$ entry, in contrast to the class of models in [24]
where a phase is in the $(\mathbf{Y}_{d})_{21}$ entry. Despite the same
texture, the difference comes from removing the phases in $\mathbf{Y}_{d}$
using $\bar{\mathbf{5}}_{Fi}$ in our case (removes one phase per column),
while [24] prefers to remove phases via $\mathbf{10}_{Fi}$ redefinitions
(removes one phase per row), so that they can make use of the
$\bar{\mathbf{5}}_{Fi}$ freedom in their implementation of the neutrino
sector, which we do not consider in the present model.

The charged lepton Yukawa matrix $\mathbf{Y}_{e}$ and the quasi Yukawa
matrices $\tilde{\mathbf{Y}}_{qq}$, $\tilde{\mathbf{Y}}_{eu}$,
$\tilde{\mathbf{Y}}_{ql}$ and $\tilde{\mathbf{Y}}_{ud}$ for model 1 are
determined from Eq. (2.15) and (2.16).

  * •

Model 2:
Corresponding to the textures in Eq. (2.14), the parametrization of the down-
and up-quark Yukawa matrix in the context of the MSSM at the GUT scale is
given by

$\displaystyle\mathbf{Y}_{d}$
$\displaystyle=\text{diag}(y^{d}_{1},y^{d}_{2},y^{d}_{3}),$ (3.13)
$\displaystyle\mathbf{Y}_{u}$
$\displaystyle=\mathbf{U}_{2}\;\text{diag}(y^{u}_{1}\,e^{i\phi_{1}},y^{u}_{2}\,e^{i\phi_{2}},y^{u}_{3})\;\mathbf{U}_{2}^{\textsf{T}},$
(3.14)

where

$\displaystyle\mathbf{U}_{2}$
$\displaystyle=\mathbf{U}_{23}(\theta^{uL}_{23},0)\;\mathbf{U}_{13}(\theta^{uL}_{13},\delta)\;\mathbf{U}_{12}(\theta^{uL}_{12},0),$
(3.15)

while $\mathbf{Y}_{e}$, $\tilde{\mathbf{Y}}_{qq}$, $\tilde{\mathbf{Y}}_{eu}$,
$\tilde{\mathbf{Y}}_{ql}$ and $\tilde{\mathbf{Y}}_{ud}$ are then determined
from Eq. (2.17) and (2.18).

In model 2, the CKM matrix is controlled entirely by the up-sector: choosing
the $\delta_{0}=0$ parametrization in the up-sector thus provides the standard
parametrization of the CKM matrix and the $3$ angles and phases in
$\mathbf{U}_{2}$ are directly those of the CKM values at the GUT scale.

###  3.2 Model parameters

Given the implementation/parameterization of the Yukawa sectors of the two
models in the previous subsection, we can now turn our attention to
comprehensively analyzing the parameters and observables. The two models
contain the following input parameters at the GUT scale:

  * •

Model 1 contains the input parameters

$\displaystyle y^{d}_{12},\ y^{d}_{21},\ y^{d}_{22},\ y^{d}_{3},\ y^{u}_{1},\
y^{u}_{2},\ y^{u}_{3},\ \theta^{uL}_{12},\ \theta^{uL}_{23},\ \varphi,\
\phi_{1},\ \phi_{2};\quad m_{0},\ M_{1/2},\ A_{0},\ \tan\beta,$ (3.16)

where the Yukawa sector parameters $y^{u,d}$, $\theta^{uL}$, $\varphi$ and
$\phi$ are used for the constructions in Eqs. (3.8)–(3.12).

  * •

Model 2 contains the following input parameters at the GUT scale:

$\displaystyle y^{d}_{1},\ y^{d}_{2},\ y^{d}_{3},\ y^{u}_{1},\ y^{u}_{2},\
y^{u}_{3},\ \theta^{uL}_{12},\ \theta^{uL}_{23},\ \theta^{uL}_{13},\ \delta,\
\phi_{1},\ \phi_{2};\quad m_{0},\ M_{1/2},\ A_{0},\ \tan\beta,$ (3.17)

where $y^{u,d}$, $\theta^{uL}$, $\delta$ and $\phi$ are used in constructions
of Eqs. (3.13)–(3.15).

Both models consist of $16$ input parameters, $12$ of which construct the
Yukawa sector, and $4$ parameters are intrinsic to the cMSSM: $m_{0}$,
$M_{1/2}$, $A_{0}$ specify the mSUGRA boundary conditions at the GUT scale,
and $\tan\beta$ is the ratio of the MSSM Higgs VEVs.

Note that the two phases $\phi_{1}$ and $\phi_{2}$ in both models do not
change the SM prediction of fermion masses and mixing angles, since they get
removed in the construction of the CKM matrix, while they do have an impact on
the quasi-Yukawa sector and hence proton decay predictions. These phases are
typically referred to as “GUT phases” in the literature [33, 34].

All initial parameters are real. Taking into account that singular values in a
Takagi decomposition are non-negative, and all quark masses are non-vanishing,
the parameters in the Yukawa sector are allowed to be in the following ranges
(for both models):

$\displaystyle y^{d}_{12},y^{d}_{21},y^{d}_{22},\
y^{d}_{1},y^{d}_{2},y^{d}_{3},\ y^{u}_{1},y^{u}_{2},y^{u}_{3}$
$\displaystyle>0,$ (3.18)
$\displaystyle\theta^{u}_{12},\theta^{u}_{23},\theta^{u}_{13}$
$\displaystyle\in[0,\pi/2],$ (3.19) $\displaystyle\delta,\phi_{1},\phi_{2}$
$\displaystyle\in[0,2\pi),$ (3.20) $\displaystyle\varphi$
$\displaystyle\in[-\pi,\pi).$ (3.21)

The cMSSM parameters can in principle also be varied in a fit. In this paper,
however, our goal is an analysis which disentangles the flavor effects from
the SUSY spectrum effects on proton decay. For this reason, we fix the cMSSM
parameters to common example values in both models, where the SUSY threshold
corrections allow for a good fit to (charged) fermion masses and the Higgs
mass. We thus consider only a particular cMSSM benchmark point in our
analyses.

To simplify the fitting procedure in the two models, we fix the gauge
couplings at the GUT scale $M_{\text{GUT}}=2\cdot 10^{16}\,\mathrm{GeV}$:
$g_{3}=0.698$, $g_{2}=0.697$ and $g_{1}=0.704$, which is consistent with the
experimental values at low energies (assuming MSSM running above the SUSY
scale of around $3\,\mathrm{TeV}$). Furthermore, since the effective triplet
mass scale $M^{\text{eff}}_{T}$ only affects the nucleon decay widths by an
overall scaling, i.e. $\Gamma\propto 1/(M^{\text{eff}}_{T})^{2}$, we can
choose a fixed value $M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}$ in the
following analyses, and can easily relate the results to any other value.
Finally, in order for the RGE of the ratio $y_{\mu}/y_{s}$ of the 2nd family
to start at $6$ at the GUT scale (based on our two models) and run to the
experimental value at the scale $M_{Z}$, we choose the sign of the $\mu$
coupling (connecting $H_{u}$ and $H_{d}$) in the MSSM to be
$\operatorname{sign}\mu=+1$, such that the dominant gluino loop diagram in
SUSY thresholds drives $y_{s}$ in the correct direction [27, 35].

###  3.3 Observables

Both models are fit to the same $14$ SM observables at low energies:

$\displaystyle y_{u},\ y_{c},\ y_{t},\ y_{d},\ y_{s},\ y_{b},\
\theta^{\text{CKM}}_{12},\ \theta^{\text{CKM}}_{23},\
\theta^{\text{CKM}}_{13},\ \delta^{\text{CKM}},\ y_{e},\ y_{\mu},\ y_{\tau},\
m_{h},$ (3.22)

corresponding to the masses of fermions in the up, down and charged lepton
sector, the CKM parameters of the quark sector, and the SM Higgs mass $m_{h}$.

Although there are $16$ input parameters in both models, $2$ GUT phases
$\phi_{1}$ and $\phi_{2}$ have no influence on the listed observables, but
only on nucleon decay discussed later, so the fit has effectively only $14$
relevant parameters. Since that number equals the number of listed
observables, a reasonable expectation is that a good fit for both models can
be obtained, with no flat directions in the $\chi^{2}$.

The values of the SM Yukawa observables at $M_{Z}$ and the error bars are
taken from [36], where for the Yukawa couplings of the charged leptons we take
a standard deviation of $1\,\%$ of the central value, which roughly
corresponds to the precision of the running. The experimental value of $m_{h}$
is given in PDG (2018) [37], where we take a standard deviation of
$3\,\mathrm{GeV}$, namely $125.6\pm 3\,\mathrm{GeV}$, due to the theoretical
uncertainties present in the calculation (for which we use FeynHiggs, see
Section 3.4).

Other observables of interest predicted by the model are the nucleon decay
rates. We use our extended SusyTC code SusyTCProton to compute $13$ different
decay channels, $8$ of which are proton decay channels and $5$ are neutron
decay channels (we sum over the neutrino flavor states, since the decay
experiments do not distinguish them); they are listed in Table 2 along with
their experimental bounds (we took the most strict bound we found across
multiple experiments). We shall not consider these observables in a $\chi^{2}$
fit, since the decay rates themselves are only predicted up to an overall
factor, unless the effective triplet mass parameter $M^{\text{eff}}_{T}$ is
also known. Taking a high enough $M^{\text{eff}}_{T}$ can always avoid all
proton decay bounds for $D=5$ contributions (we remind the reader that this
parameter in our setup is not directly connected to any physical scale of
triplets in the theory, so it can be taken above the Planck scale). Despite
the lack of predictivity in nucleon rates themselves, our setup does allow for
unambiguous predictions of branching ratios and, more generally, ratios of
decay rates in different channels, e.g. a ratio between a proton and neutron
decay channel.

Table 2:  The $90\,\%$ confidence level exclusion bounds on proton and neutron decay from the Super-Kamiokande experiment for various ($B-L$)-preserving decay channels with a meson and a lepton in the final state. For each decay channel the corresponding lifetime bound $\tau/\mathcal{B}$ and the translated partial decay width $\Gamma_{\text{partial}}$ are listed, where $\mathcal{B}$ is the branching ratio for the decay channel. The data is taken from Figure $5$-$3$ in [38] and further updated with other sources [39, 40, 8, 41, 42, 43].  decay channel |  $\tau/\mathcal{B}$ $[\text{year}]$ |  $\Gamma_{\text{partial}}$ $[\mathrm{GeV}]$
---|---|---
Proton |  |
$p\rightarrow\pi^{0}\,e^{+}$ | $>2.4\cdot 10^{34}$ | $<8.7\cdot 10^{-67}$
$p\rightarrow\pi^{0}\,\mu^{+}$ | $>1.6\cdot 10^{34}$ | $<1.3\cdot 10^{-66}$
$p\rightarrow\eta^{0}\,e^{+}$ | $>4.1\cdot 10^{33}$ | $<5.1\cdot 10^{-66}$
$p\rightarrow\eta^{0}\,\mu^{+}$ | $>1.2\cdot 10^{33}$ | $<1.7\cdot 10^{-65}$
$p\rightarrow K^{0}\,e^{+}$ | $>1.1\cdot 10^{33}$ | $<1.9\cdot 10^{-65}$
$p\rightarrow K^{0}\,\mu^{+}$ | $>1.6\cdot 10^{34}$ | $<1.3\cdot 10^{-66}$
$p\rightarrow\pi^{+}\,\bar{\nu}$ | $>2.8\cdot 10^{32}$ | $<7.4\cdot 10^{-65}$
$p\rightarrow K^{+}\,\bar{\nu}$ | $>6.6\cdot 10^{33}$ | $<3.2\cdot 10^{-66}$
Neutron |  |
$n\rightarrow\pi^{-}\,e^{+}$ | $>2.1\cdot 10^{33}$ | $<1.0\cdot 10^{-65}$
$n\rightarrow\pi^{-}\,\mu^{+}$ | $>9.9\cdot 10^{32}$ | $<2.1\cdot 10^{-65}$
$n\rightarrow\pi^{0}\,\bar{\nu}$ | $>9.9\cdot 10^{32}$ | $<2.1\cdot 10^{-65}$
$n\rightarrow\eta^{0}\,\bar{\nu}$ | $>5.6\cdot 10^{32}$ | $<3.7\cdot 10^{-65}$
$n\rightarrow K^{0}\,\bar{\nu}$ | $>1.2\cdot 10^{32}$ | $<1.7\cdot 10^{-64}$

###  3.4 Calculational procedure

For given values of the parameters in model 1 and 2, listed in Eq. (3.16) and
(3.17) respectively, the Yukawa matrices and the soft terms are implemented at
the GUT scale $M_{\text{GUT}}=2\cdot 10^{16}\,\mathrm{GeV}$ according to
Section 3.1. Moreover, the dimension $5$ operators which are relevant for
nucleon decay are calculated at $M_{\text{GUT}}$ as well, using the formulas
from Eq. (2.8) and (2.9).

The RG running in the MSSM is performed by means of the Mathematica package
SusyTCProton (cf. Appendix C.1), which is an extension of the Mathematica
package SusyTC [13] used within REAP [14]. The running of the Yukawa matrices,
gauge couplings and soft terms is computed at $2$-loop, whereas the running of
the dimension $5$ operators is computed at $1$-loop.

The SUSY scale $M_{\text{SUSY}}$ is determined dynamically in SusyTCProton by
the geometric mean of the two up-type squark masses which have the largest
mixing with the stop flavour eigenstates $\tilde{t}_{1}$ and $\tilde{t}_{2}$,
namely $M_{\text{SUSY}}=\sqrt{m_{\tilde{t}_{1}}m_{\tilde{t}_{2}}}$. The values
of the Yukawa matrices, gauge couplings and soft terms at $M_{\text{SUSY}}$
are used to calculate the partial decay widths of the proton and neutron
listed in Table 2 by means of the Mathematica package ProtonDecay (cf.
Appendix C.2). In addition, the electroweak (EW) Higgs mass is calculated by
using FeynHiggs 2.16.0 [44, 45, 46, 47, 48, 49, 50, 51].

In SusyTCProton, the MSSM is matched to the SM by taking the SUSY threshold
corrections into account. Furthermore, the running from the SUSY scale to the
$Z$-boson mass scale $M_{Z}=91.2\,\mathrm{GeV}$ is calculated at $2$-loop. At
$M_{Z}$, the observables in the quark and charged lepton sector, namely the
Yukawa couplings and the CKM parameters, are computed.

For the statistical model analyses we use the following procedures:

  * •

As a measure of goodness of fit we define the $\chi^{2}$ function in the usual
way:

$\displaystyle\chi^{2}(\vec{x})=\sum_{i}\frac{\big{(}f_{i}(\vec{x})-y_{i}\big{)}^{2}}{\sigma_{i}^{2}}\equiv\sum_{i}\chi^{2}_{i},$
(3.23)

where the vector $\vec{x}$ consists of the input parameters of the models from
either Eq. (3.16) or (3.17), and the index $i$ goes over all observables. The
$y_{i}$ denote the central values from the (experimental) data, $\sigma_{i}$
are their corresponding standard deviation errors, while $f_{i}(\vec{x})$ are
the predictions for the $i$-th observable given the parameter point $\vec{x}$.
Some observables may be equipped with asymmetric errors $\sigma_{i+}$ and
$\sigma_{i-}$ when $f_{i}(\vec{x})>y_{i}$ or $f_{i}(\vec{x})<y_{i}$,
respectively. Each $\chi^{2}_{i}$ term in the $\chi^{2}$ function is referred
to as a pull associated to a particular observable; it will also be convenient
to refer to a particular pull by replacing the index $i$ with the label of the
$i$-th observable.

  * •

The best-fit point for a given model setup is calculated by using a
differential evolution algorithm for global minimization using custom
Mathematica code.

  * •

To calculate the posterior density of the parameters and observables for a
certain model setup, we perform an MCMC analysis by using an adaptive
Metropolis-Hastings algorithm [52]. All prior probability distributions for
the parameters are chosen to be flat. Furthermore, we calculate $12$
independent chains with $5\cdot 10^{4}$ valid points (beyond the “burn-in
phase”), which gives $6\cdot 10^{5}$ points in total for each dataset. All
chains are started at the best-fit point of the corresponding model setup.

##  4 Results for nucleon decay from SusyTCProton

In this section, we present the results of our numerical analyses for the two
models we considered. Our main interest are the predicted nucleon decay rates
for the points with a good fit to the low energy data.

Expressions for nucleon decay rates mediated by $D=5$ operators contain many
contributions, all of which depend on both the masses of the SUSY sparticles
that are used to dress the $D=5$ operators, as well as flavor effects both
from the SM Yukawa sector and the sfermion sector. More details on the
expressions used are given in Appendix B and references specified there.

The SusyTCProton.m and ProtonDecay.m packages developed for such a nucleon
decay calculation, see Appendices C.1 and C.2, can handle all these effects.
Our main interest in this paper are SUSY GUTs, and as such we are interested
in the effect the Yukawa textures are having on proton decay. For this reason,
we would like to isolate this effect from other contributions: taking cMSSM
boundary conditions already insures no extra flavor violation in the SUSY
sector, while the effects of variable sparticle masses can be removed by
fixing the soft parameters.

We therefore choose in our analysis to fix the cMSSM parameters and consider
only a benchmark point. We take the following values:

$\displaystyle\begin{split}\tan\beta&=\phantom{+}28.4,\\\
m_{0}&=\phantom{+0}4000\,\mathrm{GeV},\\\
M_{1/2}&=\phantom{+0}1700\,\mathrm{GeV},\\\ A_{0}&=-10000\,\mathrm{GeV}.\\\
\end{split}$ (4.1)

We use these values in both models, and they allow us to obtain good fits for
both777Both models have by design the same CG coefficient between the down and
charged lepton sectors in the $3$-$3$ and $2$-$2$ entries at the GUT scale,
and thus require essentially the same SUSY threshold corrections. A good fit
to the fermion data is thus possible with the same values of soft parameters.,
thus being able to compare them on an equal footing and isolate the flavor
effects. We have checked that EW vacuum stablility for this parameter point by
using Vevacious 1.2.03 [53], SPheno 4.0.3 [16, 17], and the predefined model
from SARAH 4.14.1 [54, 55] with possible charge breaking via stau VEVs. The
output shows that the produced EW vacuum is stable.

Fixing the cMSSM parmeters leaves $12$ free parameters to vary for a fit, with
10 having an effect on the fit, while varying the GUT phases $\phi_{1}$ and
$\phi_{2}$ will only effect the nucleon decay rates.

###  4.1 Best-fit points and nucleon decay

Minimizing the $\chi^{2}$ with respect to the flavor parameters, we obtain a
best-fit point for each of the models. The parameter values for the best-fit
points are shown in Table 3, while the $\chi^{2}$ results and various dominant
pulls are shown in Table 4. Most of the tension in the fit is coming from the
pull in the $s$-quark mass.

Table 3: The flavor input parameters for the best-fit point of models 1 and 2.
The GUT phases $\phi_{1}$ and $\phi_{2}$ have no influence on the predictions.
Model 1
---
$y^{u}_{1}$ | $y^{u}_{2}$ | $y^{u}_{3}$ | $\theta^{uL}_{12}$ | $\theta^{uL}_{23}$ | $\varphi$ | $y^{d}_{12}$ | $y^{d}_{21}$ | $y^{d}_{22}$ | $y^{d}_{3}$
2.84 | 1.39 | 4.90 | 8.66 | 3.53 | $-4.49$ | 4.84 | 5.59 | 2.11 | 1.54
$\times 10^{-6}$ | $\times 10^{-3}$ | $\times 10^{-1}$ | $\times 10^{-2}$ | $\times 10^{-2}$ | $\times 10^{-3}$ | $\times 10^{-4}$ | $\times 10^{-4}$ | $\times 10^{-3}$ | $\times 10^{-1}$
Model 2
$y^{u}_{1}$ | $y^{u}_{2}$ | $y^{u}_{3}$ | $\theta^{uL}_{12}$ | $\theta^{uL}_{23}$ | $\theta^{uL}_{13}$ | $\delta$ | $y^{d}_{1}$ | $y^{d}_{2}$ | $y^{d}_{3}$
2.85 | 1.39 | 4.90 | 2.27 | 3.45 | 7.40 | 5.88 | 1.24 | 2.19 | 1.54
$\times 10^{-6}$ | $\times 10^{-3}$ | $\times 10^{-1}$ | $\times 10^{-1}$ | $\times 10^{-2}$ | $\times 10^{-3}$ |  | $\times 10^{-4}$ | $\times 10^{-3}$ | $\times 10^{-1}$
Table 4: The total $\chi^{2}$ and the dominant pulls $\chi^{2}_{i}$, as defined in Eq. (3.23), for the best-fit points of models 1 and 2 of Table 3. The subdominant pulls $\chi^{2}_{i}<0.01$ are not shown.  model |  $\chi^{2}$ |  $\chi^{2}_{y_{d}}$ |  $\chi^{2}_{y_{s}}$ |  $\chi^{2}_{y_{b}}$ |  $\chi^{2}_{\theta^{CKM}_{23}}$ |  $\chi^{2}_{\theta^{CKM}_{13}}$ |  $m_{h}$
---|---|---|---|---|---|---|---
1 | 4.55 |  0.20 |  3.29 |  0.11 |  0.08 |  0.18 |  0.67
2 | 5.92 |  0.01 |  4.79 |  0.11 |  0.14 |  0.18 |  0.67

We can now use these best-fit points to obtain predictions of proton and
neutron decay rates in various decay channels. We postpone a more
comprehensive comparison of all the channels to Section 4.2, and will only
focus here on some selected and most interesting aspects of the results,
namely the dependence of the rates and their ratios on the GUT phases
$\phi_{1}$ and $\phi_{2}$, since these phases are not determined by the fit
and we thus have an embedded ambiguity of their values in the models. We plot
the results in this section as contour plots in the $\phi_{1}$-$\phi_{2}$
plane, and take a fixed value for the effective triplet mass
$M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}$ where needed to obtain specific
numerical results. Each of these plots is equipped with it’s own scale legend
to the right. Note that despite the same color range for contour regions, the
values represented by these colors change between different plots. Also, some
of the plots provide information about the location and value of the minimum
and maximum, which are added to the plots in red color. This information of
the minimum and maximum is important, since it indicates that the predictions
robustly lie in a compact region, and the decay rates cannot be made to vanish
by a suitable choice of $\phi_{1}$ and $\phi_{2}$.

To draw the plots, we computed the values on an equidistant 2D-grid with
increments $\pi/30$ for both the $\phi_{1}$ and $\phi_{2}$ directions, while a
separate minimization/maximization routine was run to obtain the extrema in
the $\phi_{1}$-$\phi_{2}$ plane.

Figure 1: Contour plots for the GUT phase dependence of the dominant proton
decay channel $p\to K^{+}\bar{\nu}$ (top) and neutron decay channel $n\to
K^{0}\bar{\nu}$ (bottom) for the best-fit points of model 1 (left) and model 2
(right). We assumed $M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}$. Figure 2:
Contour plots for the GUT phase dependence of the $p\to K^{0}e^{+}$ proton
decay channel for the best-fit points of model 1 (left) and model 2 (right).
We assumed $M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}$.

Figure 3: Contour plots for the GUT phase dependence of the decay rate ratios
$\Gamma_{p\to K^{+}\bar{\nu}}/\Gamma_{p\to\pi^{+}\bar{\nu}}$ (top) and
$\Gamma_{p\to K^{+}\bar{\nu}}/\Gamma_{n\to K^{0}\bar{\nu}}$ (bottom) for the
best-fit points of model 1 (left) and model 2 (right). The results are
independent of the effective triplet mass $M^{\text{eff}}_{T}$.

We make the following observations of the results:

  * •

The dominant decay channels of both the proton and neutron are those with
neutrinos in the final state; furthermore, the decay channel with a kaon
always dominates over that of the pion. We plot the GUT phase dependence of
the dominant proton and neutron decay rate for each model in Figure 1.

We see that the rates for each process in a given model can differ by a factor
of up to $\sim 3$, depending on the unknown GUT phases, but can never vanish
completely. Observing the decay rate values in each panel, we see that the
dominant neutron decay channel dominates over the proton decay channel in the
case of model 1, while the opposite is true in model 2. Also, for model 1 the
$\phi_{1}$-$\phi_{2}$ regions minimizing the proton decay rates roughly
correspond to those minimizing the neutron decay rate; similarly, the proton
and neutron decay rates are simultaneously maximal in a similar region. For
model 2, on the other hand, minimizing the dominant proton decay roughly
maximizes the neutron decay rate, and vice versa.

  * •

We selected to show another decay channel, $p\to K^{0}e^{+}$, in Figure 2. The
decay rates of this channel are highly suppressed compared to the dominant
channels, and the plots show a much simpler GUT phase dependence. For model 2,
where the charged lepton Yukawa is diagonal, the rate effectively depends only
on the phase $\phi_{1}$; in model 1, where the charged lepton mass matrix is
not diagonal, the relative changes of the rates are smaller, but the
$\phi_{2}$ dependence is not negligible.

  * •

A result independent of the effective triplet mass can be obtained by
considering ratios of decay rates. We choose to present two such ratios in the
plots: the decay ratio $\Gamma_{p\to
K^{+}\bar{\nu}}/\Gamma_{p\to\pi^{+}\bar{\nu}}$ represents the ratio between
the two most dominant proton decay channels, while the ratio $\Gamma_{p\to
K^{+}\bar{\nu}}/\Gamma_{n\to K^{0}\bar{\nu}}$ is between the dominant proton
and neutron channels. The GUT phase dependence of the ratios is shown in
Figure 3.

We observe that the plots for ratios have a similar feature to those for the
dominant rates: the ratios are simultaneously minimized or maximized for model
1, but minimizing one ratio roughly maximizes the other in model 2. This
implies that the process $p\to\pi^{+}\bar{\nu}$, which is not shown in Figure
1, has a similar minimizing/maximizing region to the dominant proton decay
process; the switch of the regions happens for model 2 only when migrating to
neutron decays.

Another important observation is that the predicted range for a given ratio
does not overlap when comparing the two models. This implies that
experimentally determining such a ratio of nucleon decay channels would enable
us to discriminate between the two flavor GUT models, despite the unknown GUT
phases and despite the two models having the same soft SUSY-breaking
parameters. This difference in the predicted ratios is thus purely a flavor
effect, as could be expected; we provide additional discussion of using ratios
for discrimination in Section 4.2.

###  4.2 MCMC results for nucleon decay

The nucleon decay results of the previous section involved varying only the
two GUT phases $\phi_{1}$ and $\phi_{2}$, while we keep the other flavor
parameters confined to the best-fit point in each model.

We now analyze the two models by varying all flavor parameters of Eq. (3.16)
and (3.17), such that the fit is still close to the best-fit point; the cMSSM
prameters are still fixed, however, to their benchmark values of Eq. (4.1).
For this analysis we use the Markov Chain Monte Carlo (MCMC) method, which
allows us to obtain estimates of posterior probability densities on the
parameter space. This results in obtaining highest posterior density (HPD)
intervals for the observables of interest. For more technical details on our
MCMC procedure see Section 3.4 (and references there-in).

The plots in this section show HPD intervals of nucleon decay rates or ratios
predicted by the two models. The HPD intervals are colored red (green) for
model 1 (2). We show intervals corresponding to the $1$-$\sigma$ and
$2$-$\sigma$ HPD regions as darker and lighter colored, respectively, and
overlayed one on top of the other.

We make the following observations of the results:

  * •

Our SusyTCProton and ProtonDecay packages are able to compute $13$ different
decay channels (specified in Table 2). The MCMC results for the HPD intervals
in each channel are shown in Figure 4. Note that we assumed the effective
triplet mass to be $M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}$; the decay rate
results for a different $M^{\text{eff}}_{T}$ can be easily obtained by noting
that $\Gamma\propto(M^{\text{eff}}_{T})^{-2}$, implying that we linearly
translate the intervals down by $2$ on the log scale for each order of
magnitude that $M^{\text{eff}}_{T}$ increases over the reference value. The
experimental bounds for the various channels are shown as blue bars; choosing
a different $M^{\text{eff}}_{T}$ can control whether we violate any of the
experimental bounds or not with the predictions. For the reference value of
$M^{\text{eff}}_{T}$, the experimental bound on $p\to K^{+}\bar{\nu}$ would
exclude model 2 but not model 1.

We see that the dominant channels are those with a neutrino in the final
state, while those with a charged lepton have decay rates at least $4$ orders
of magnitude below the dominant ones. The dominant proton decay channel is
$p\to K^{+}\bar{\nu}$, while the dominant neutron decay channel is $n\to
K^{0}\bar{\nu}$. Which of the two is bigger differs between the two models.
This result is consistent with the dominant decay channels when only the GUT
phases were varied in Section 4.1.

Also, we notice that when a charged lepton is a decay product, the predictions
for the rates are much sharper when the charged lepton is $e^{+}$ (the
$2$-$\sigma$ ranges vary $\sim 20\,\%$ around the average), while the
prediction for $\mu^{+}$ in the final state typically spans $1$ to $2$ orders
of magnitude. Given the decay rates of these processes, however, they are at
least $4$ orders of magnitude below the dominant channels and thus
experimentally not within reach.

Comparing the HPD intervals between the two models process by process, we
notice that for a given $M^{\text{eff}}_{T}$ most channels are predicted with
a similar decay rate, except for $p\to K^{+}\bar{\nu}$; this process is thus
most suitable for discrimination between the two models, assuming another
channel has been measured as well. Luckily, the decay $p\to K^{+}\bar{\nu}$
also has either the biggest or 2nd biggest decay rate, so it would likely be
among the first experimentally measured channels if one of the considered
models is realized in nature.

  * •

The quantities independent of $M^{\text{eff}}_{T}$ are the ratios of various
decay rates. Note that a ratio for each parameter point in the MCMC dataset is
computed separately, and only then is statistics performed on these values;
this is an important point, since there could be additional correlations
between various decay channels in a model not apparent when viewing only the
HPD intervals of decay rates, i.e. when considering the numerator and
denominator in the ratio separately.

In Figure 5 we show selected ratios which best discriminate between the two
models. They effectively have no overlap of predicted $2$-$\sigma$ HPD
intervals when comparing the two models, thus measuring even one such ratio
has sufficient discriminating power. All ratios between the $13$ channels not
shown in the Figure all have some overlap in the $2$-$\sigma$ intervals. More
broadly, given any specific model, predictions for a set of ratios effectively
provides a fingerprint of that model, in principle allowing for it to be
excluded if multiple channels of decays were measured experimentally.

The most interesting ratios experimentally are those involving the dominant
decay channels, which would be the first ones measured. Also, a very big or
small prediction of the ratio (compared to a $1:1$ ratio) implies a very big
discrepancy in the rates of the numerator and denominator, again making it
unlikely that both processes could be measured in the near future. Given these
considerations, the most experimentally relevant ratios are the first three,
which compare the high rate processes with a neutrino in the final state.

All but one of the ratios in Figure 5 include the $p\to K^{+}\bar{\nu}$ decay
as one of the two processes, confirming our comment for Figure 4 that this
process is a good discriminator between the two models. In terms of the
sharpness of the prediction of ratios, the value for the fourth ratio in the
figure has the smallest allowed range: the $2$-$\sigma$ HPD intervals for
$\Gamma_{p\to K^{0}e^{+}}/\Gamma_{p\to\eta^{0}e^{+}}$ are $[5.7,6.9]$ and
$[7.5,9.3]$ for models 1 and 2, respectively. The decay channels in that ratio
have very small rates and are experimentally out of reach in the considered
scenarios; assuming an experimental measurement of such a ratio (or even just
one channel and having bounds on the other) provides, however, an excellent
way to exclude the model under consideration.

A final note on ratios: some ratios are predicted exactly the same with a
value $(1/\sqrt{2})^{\pm 1}$ for every point due to the form of the low-energy
effective Lagrangian in chiral perturbation theory. In particular, these are
the $3$ ratios where only a replacement of $\pi^{0}$ with $\pi^{\pm}$ takes
place in the final state:

$\displaystyle\frac{\Gamma_{p\to\pi^{0}e^{+}}}{\Gamma_{n\to\pi^{-}e^{+}}},\quad\frac{\Gamma_{p\to\pi^{0}\mu^{+}}}{\Gamma_{n\to\pi^{-}\mu^{+}}},\quad\frac{\Gamma_{p\to\pi^{+}\bar{\nu}}}{\Gamma_{n\to\pi^{0}\bar{\nu}}}.$
(4.2)

These ratios are fixed and have no dependence on any high-energy or SM
parameters.

Figure 4: The $1$-$\sigma$ (dark) and $2$-$\sigma$ (light) HPD regions
predicted via MCMC for nucleon decay rates in both models using a benchmark
cMSSM point and a fixed effective triplet mass
$M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}.$ The blue line segments indicate
current experimental bounds.  Figure 5: The $1$-$\sigma$ (dark) and
$2$-$\sigma$ (light) HPD regions predicted via MCMC for selected ratios of
decay rates of two channels; we assume the same cMSSM benchmark point in both
models.

##  5 Conclusions

While the observation of nucleon decay would be a smoking gun of Grand Unified
Theories (GUTs) in general, the ratios between the decay rates of the various
channels carry rich information about the specific GUT model realization. We
investigated in this paper these nucleon decay fingerprints of two flavor SUSY
GUT models based on the group $\mathrm{SU}(5)$. The two models differ in the
texture of their Yukawa matrices, and the types of GUT operators generating
the SM fermion masses, cf. Section 2. Both admit essentially equally good fits
to the low energy Yukawa data of SM masses and mixing angles (with a
$\chi^{2}\approx 5$, cf. Section 4.1), even given the same sparticle spectrum.
Yet their nucleon decay computation shows that the partial decay rates of
various channels are good observables to discriminate between them, despite
the decay ambiguity introduced by $2$ unknown GUT phases.

We investigated the dependence of the decay rate on the GUT phases at the
best-fit point for the flavor parameters of the model in Section 4.1. The
nucleon decay fingerprinting has been performed by computing the decay rates
in $8$ different proton decay channels and $5$ different nucleon decay
channel, cf. Table 2. We obtained their predicted ranges by varying the flavor
parameters in a way still compatible with a low $\chi^{2}$ of the low-energy
observables using the MCMC algorithm. The main results of the analysis are
presented in Section 4.2, see especially Figures 4 and 5. The general features
are those expected in a SUSY GUT: the $p\to K^{+}\bar{\nu}$ and $n\to
K^{0}\bar{\nu}$ channels are the dominant modes of proton and neutron decay,
respectively, and the rates with kaons in the final states are bigger than
those with pions when the other decay product is the same. Other features
found in the fingerprints are that the neutrino decay channels dominate over
those with charged leptons in the final state and that the channels with
$e^{+}$ in the final state are predicted much more sharply than those with
$\mu^{+}$.

The decay rates invariably depend on the scale of the heavy color triplets
that mediate $D=5$ nucleon decay. We considered a doublet-triplet setup where
the overall scale of the rates depends only on the triplet effective mass
parameter $M^{\text{eff}}_{T}$. Robust and $M^{\text{eff}}_{T}$-independent
predictions for both models are obtained by considering ratios of decay rates
between different channels. The key connection of nucleon decay to flavor is
that the couplings of the SM fermions to the $D=5$ nucleon decay mediating
triplets are related to the Yukawa coupling of the theory, since they arise
from the same GUT operators, which is a generic feature of any GUT
realization. The results for the considered models show that the impact of
their different flavor structure is already sufficient to discriminate between
them. In particular, the dominant proton decay mode $p\to K^{+}\bar{\nu}$
turns out to differ relative to other dominant modes with neutrinos in the
final states such as $p\to\pi^{+}\bar{\nu}$ and $n\to\pi^{0}\bar{\nu}$.

The computations of nucleon decay rates for the $13$ channels in this paper
were performed by use of the Mathematica software modules SusyTCProton and
ProtonDecay, which we developed for the purposes of this paper and are rolling
them out for public use. The SusyTCProton is an extension of the SusyTC
package within the REAP framework, and allows to compute the RG running of
$D=5$ operators from the GUT scale to the $Z$ scale alongside the running of
all other parameters in a softly broken MSSM, performing the MSSM to SM
matching and dressing the $D=5$ operators at the SUSY scale. The ProtonDecay
module admits as input either the SusyTCProton output or an SLHA file amended
with $D=5$ operator values at the SUSY scale, and computes the resulting
nucleon decay rates in the various channels. The module is able to include
$D=6$ contributions to nucleon decay if the operator values are provided at
some high scale where the SM is taken as a valid effective theory, so the
package can be used for computing nucleon decay in non-SUSY GUTs as well. We
list the documentation for the use of these software packages (a “User’s
guide”) in Appendix C.

When used in tandem, the above software packages allow for the computation of
the decay rates in any softly broken SUSY GUT if one is able to input the
parameter values at a scale where the MSSM becomes valid (in our case this was
immediately below the GUT scale). It computes the rates in all generality,
including all the flavor effects from both the Yukawa and SUSY sectors, and
considers all diagrams in these decays following [12]. This allows for the
fingerprinting of nucleon decay within such a model, opening up a new avenue
for the model to be tested in experiments. The fingerprint is most robust when
considered at the level of ratios of decay rates for various channels.
Experimentally, the measurement of decay events in just one channel either for
the proton or neutron, alongside experimental bounds in the other channels,
already allows for strategies for cross-checking with the fingerprint, since
the model prediction might suggest that decay events in some other channels
should have been observed already as well. A situation with decay events
measured in multiple channels of course allows for even better tests of the
model. These fingerprinting techniques can be readily used either to
discriminate between otherwise indistinguishable models from the point of view
of current experiments, or to possibly rule them out.

This paper represents a thorough analysis of nucleon decay in two flavor GUT
toy models, demonstrating the efficacy of using their nucleon decay
fingerprints to distinguish them. Alongside this analysis, we also roll out
new software tools which can be employed by the community for further
investigation of nucleon decay in models of various types, thus adding to the
nucleon decay fingerprint collection.

## Acknowledgements

This work has been supported by the Swiss National Science Foundation. The
authors would like to thank Liliana Velasco-Sevilla for discussion.

##  Appendix A Details on doublet-triplet splitting

DT-splitting [56] needs to be addressed in any GUT theory. In the SUSY
setting, the main issue is the following dichotomy:

  * [D]:

One mass-eigenpair of a doublet $D\sim(\mathbf{1},\mathbf{2},+1/2)$ and anti-
doublet $\overline{D}\sim(\mathbf{1},\mathbf{2},-1/2)$ is identified with
$H_{u}\oplus H_{d}$ of the MSSM, and their mass should be negligible compared
to the GUT scale in order for EW symmetry breaking to occur. MSSM extensions
with multiple doublets are of course possible, but we consider the MSSM as the
effective theory below $M_{\text{GUT}}$ in this paper.

  * [T]:

The triplet states $T\sim(\mathbf{3},\mathbf{1},-1/3)$ and anti-triplet states
$\overline{T}\sim(\mathbf{\bar{3}},\mathbf{1},+1/3)$ can mediate dimension $5$
proton decay, so they should be heavy (at the GUT scale).

The two conditions are difficult to fulfill simultaneously due to irreducible
representations of he unified group always containing a triplet alongside a
doublet. Even in $\mathrm{SU}(5)$, which is the smallest possible unified
group containing the SM and a subgroup of all other (true) unified group
candidates, both doublets and triples are present in the representations
$\mathbf{5}$ and $\mathbf{45}$ (and their conjugates). The only (low-
dimensional) exception is the representation $\mathbf{50}$ of
$\mathrm{SU}(5)$, which contains a triplet but no doublet; this feature is
used in the missing partner mechanisms (see e.g. [57]) to achieve DT-
splitting.

Since $D=5$ proton decay is mediated by heavy triplets, it unavoidably depends
on the details of DT splitting, which in turn depends on the breaking sector
of the theories. We wish to be as model independent in the breaking sector as
possible for our nucleon decay analysis in this paper; for this reason we
assumed the simplest possible scenario where only one triplet couples to the
fermionic representations. As discussed in Section 2.1, this implies that only
the $1$-$1$ element of the inverse mass matrix $\mathbf{M}_{T}^{-1}$ is
important for proton decay.

This choice represents a minimal parametrization of DT-splitting with only one
parameter: the effective triplet mass $M^{\text{eff}}_{T}$, cf. Eq. (2.10).
Due to Eq. (2.8)-(2.9), the amplitude for all $D=5$ nucleon decay channels is
proportional to $(M^{\text{eff}}_{T})^{-1}$, and the decay rates are
proportional to $(M^{\text{eff}}_{T})^{-2}$. Raising the effective triplet
mass can suppress the overall proton decay while keeping the decay rate ratios
fixed (ignoring $D=6$ contributions).

In this appendix section we present two concrete realizations of such a
scenario; despite multiple possibilities for the underlying dynamics, only the
effective mass $M^{\text{eff}}_{T}$ is the relevant parameter. Note that while
this parameter has mass dimension $1$, it does not (necessarily) represent a
physical mass of a particle or a characteristic energy scale at which a
process takes place, and as such can take trans-Planckian values. Indeed, for
our numerical analysis in Section 4 we took for the effective triplet mass the
benchmark value $M^{\text{eff}}_{T}=10^{19}\,\mathrm{GeV}$, so the proposed
mechanisms should be able to accommodate such a high parameter value.

We present two possibilities how to achieve a setup admitting high
$M^{\text{eff}}_{T}$: a modification of the ordinary missing partner mechanism
is discussed in Appendix A.1, while the double missing partner mechanism is
described in Appendix A.2.

###  A.1 The modified missing partner mechanism

####  A.1.1 The missing partner mechanism — the usual setup

In the usual setup, the missing partner mechanism in $\mathrm{SU}(5)$ [57, 58]
requires that the color triplets $T$ (or anti-triplets $\overline{T}$) come
from the following $\mathrm{SU}(5)$ representations:

$\displaystyle\mathbf{5}_{H}\oplus\mathbf{\bar{5}}_{H}\oplus\mathbf{50}_{H}\oplus\mathbf{\overline{50}}_{H}.$
(A.1)

The Higgs doublets $(\mathbf{1},\mathbf{2},\pm 1/2)$ on the other hand are
present only in the $\mathbf{5}_{H}\oplus\mathbf{\bar{5}}_{H}$ part; as
mentioned earlier in Appendix A, the representations $\mathbf{50}$
($\overline{\mathbf{50}}$) of $\mathrm{SU}(5)$ remarkably contains only a
triplets (anti-triplet), but not an (anti-)doublet.

The doublet and anti-doublet directly correspond to the MSSM doublets $H_{u}$
and $H_{d}$, respectively:

$\displaystyle H_{u}$ $\displaystyle\subseteq\mathbf{5}_{H},$ $\displaystyle
H_{d}$ $\displaystyle\subseteq\mathbf{\bar{5}}_{H}.$ (A.2)

The (anti-)triplets in Eq. (A.1) on the other hand are labeled by

$\displaystyle T_{1}$ $\displaystyle\subseteq\mathbf{5}_{H},$
$\displaystyle\overline{T}_{1}$ $\displaystyle\subseteq\mathbf{\bar{5}}_{H},$
(A.3) $\displaystyle T_{2}$ $\displaystyle\subseteq\mathbf{50}_{H},$
$\displaystyle\overline{T}_{2}$
$\displaystyle\subseteq\mathbf{\overline{50}}_{H}.$ (A.4)

We assume then an $\mathrm{SU}(5)$-level superpotential of the form

$\displaystyle W_{DT}$
$\displaystyle=\lambda\;\mathbf{5}_{H}\cdot\mathbf{\overline{50}}_{H}\cdot\mathbf{75}_{H}+\bar{\lambda}\;\mathbf{\overline{5}}_{H}\cdot\mathbf{50}_{H}\cdot\mathbf{75}_{H}+m_{50}\;\mathbf{50}_{H}\cdot\mathbf{\overline{50}}_{H},$
(A.5)

where the $G_{321}$ singlet component in $\mathbf{75}_{H}$ obtains a GUT scale
VEV $v=\langle\mathbf{75}\rangle$.888Alternatively, one could avoid the use of
the representation $\mathbf{75}$ by utilizing a symmetric product of a pair of
$\mathbf{24}$s instead. Crucially, the explicit mass term
$\mathbf{5}_{H}\cdot\mathbf{\overline{5}}_{H}$ is missing (this can be
achieved for example by imposing discrete symmetries). The superpotential of
Eq. (A.5) yields the mass terms

$\displaystyle\mathbf{W}\supset
H_{u}\;\mathbf{M}_{D}\;H_{d}+\begin{pmatrix}T_{1}&T_{2}\\\
\end{pmatrix}\;\mathbf{M}_{T}\;\begin{pmatrix}\overline{T}_{1}\\\
\overline{T}_{2}\\\ \end{pmatrix}+\cdots,$ (A.6)

where the doublet and triplet mass matrices $\mathbf{M}_{D}$ and
$\mathbf{M}_{T}$ are given by

$\displaystyle\mathbf{M}_{D}$ $\displaystyle=0,$ $\displaystyle\mathbf{M}_{T}$
$\displaystyle=\begin{pmatrix}0&\lambda v\\\ \bar{\lambda}v&m_{50}\\\
\end{pmatrix}.$ (A.7)

In general we have $\lambda,\bar{\lambda},v,m\in\mathbb{C}$, but the phases of
$\lambda v$, $\bar{\lambda}v$ and $m$ can be absorbed into the overall phase
definition of the irreps $\mathbf{\overline{50}}_{H}$, $\mathbf{75}_{H}$ and
$\mathbf{50}_{H}$ in Eq. (A.5), respectively. This implies that all entries in
$\mathbf{M}_{T}$ can be taken to be real numbers.

This situation in Eq. (A.7) clearly achieves DT-splitting, since the doublet
remains effectively massless, while all triplets are roughly at the GUT scale
due to the presence of the VEV $v$. Furthermore, using Eq. (2.10), such an
$\mathbf{M}_{T}$ gives the effective triplet mass to be

$\displaystyle M^{\text{eff}}_{T}$
$\displaystyle=-\lambda\bar{\lambda}\,v^{2}/m_{50}.$ (A.8)

####  A.1.2 The perturbativity problem and modifying the usual setup

Clearly, $M^{\text{eff}}_{T}$ from Eq. (A.8) can be enhanced by taking a small
value for the parameter $m_{50}$, i.e. $m_{50}\ll M_{\text{GUT}}$. Although
the scale of the triplet masses is not impacted, since those are controlled by
the VEV $v$, the remainder of $\mathbf{50}\oplus\overline{\mathbf{50}}$
(consisting of particles of $1\times 1$ mass matrices, since they are not
found in $\mathbf{5}_{H}\oplus\overline{\mathbf{5}}_{H}$) is then moved to the
scale of $m_{50}$.

This turns out to present a problem, since introducing a lot of matter below
$M_{\text{GUT}}$ may increase the gauge coupling value to such a degree at
$M_{\text{GUT}}$ that the Landau pole of the coupling is unacceptably close to
the unification scale. Let’s investigate this situation:

  1. 1.

The solution to the 1-loop RG running of the gauge coupling $g$ between two
scales $\mu_{r}$ and $\mu^{\prime}_{r}$ is described by

$\displaystyle\alpha^{-1}(\mu_{r}^{\prime})=\alpha^{-1}(\mu_{r})+\tfrac{\beta}{2\pi}\log(\mu_{r}^{\prime}/\mu_{r}),$
(A.9)

where $\alpha=g^{2}/(4\pi)$ and $\beta$ depends on the model. We remind the
reader that in SUSY theories $\beta=3D_{2}(V)-D_{2}(C)$, where $D_{2}$ is the
Dynkin index, with $V$ and $C$ referring to all vector and chiral
supermultiplets in the theory with masses below a given renormalization scale,
respectively. Adding a lot of matter (chiral supermultiplets) makes $\beta$
smaller (more negative), causing $g$ to blow up sooner when running to higher
scales.

  2. 2.

Suppose we change the theory by moving some representations from the scale
$\mu_{0}$ to $\mu_{0}/c$ ($c>1$), thus changing the $\beta$-function between
$\mu_{0}/c$ and $\mu_{0}$ from $\beta$ to $\beta^{\prime}$. The new coupling
$\alpha^{\prime-1}$ and the old coupling $\alpha^{-1}$ should be the same at
$\mu_{0}/c$, and using Eq. (A.9) we can relate them at the scale $\mu_{0}$ via

$\displaystyle\alpha^{\prime-1}(\mu_{0})$
$\displaystyle=\alpha^{-1}(\mu_{0})+\tfrac{\beta^{\prime}-\beta}{2\pi}\log c,$
(A.10)

which depends on the difference of the two beta functions
$\Delta\beta=\beta^{\prime}-\beta$. For simplicity, we consider the unified
coupling and a setup where an entire $\mathbf{50}\oplus\overline{\mathbf{50}}$
is moved from $\mu_{0}=M_{\text{GUT}}$ to $\mu_{0}/c$: since it is a complete
$\mathrm{SU}(5)$ representation, it does not change the unification scale
$\mu_{0}$, but only the value of the unified coupling from
$\alpha^{-1}(\mu_{0})$ to $\alpha^{\prime-1}(\mu_{0})$.

  3. 3.

The scale of the Landau pole $\mu_{L}$ is defined implicitly by
$\alpha^{-1}(\mu_{L})=0$, i.e. when $g$ has a pole. We define $C$ as the
factor above the GUT scale where the Landau pole occurs, i.e.
$\mu_{L}=C\mu_{0}$. The Landau pole condition in the new theory
$\alpha^{\prime-1}(C\mu_{0})=0$ can be related to $\alpha^{\prime-1}(\mu)$ by
using the solution in Eq. (A.9) (assuming $\beta_{f}$ for the beta function
above $\mu_{0}$ and no new thresholds); when combined with Eq. (A.10) this
relates $c$, $C$ and the original gauge coupling $\alpha^{-1}$ at the GUT
scale:

$\displaystyle 2\pi\alpha^{-1}(\mu_{0})+\beta_{f}\log C+\Delta\beta\;\log c$
$\displaystyle=0.$ (A.11)

  4. 4.

We now turn to our specific case of $\mathrm{SU}(5)$ and assume the full
theory has at least the following chiral supermultiplets in the fermionic and
Higgs sectors (no messengers, no flavons, we assume GUT breaking only with the
$\mathbf{24}$ and $\mathbf{75}$):

$\displaystyle
3\times(\mathbf{10}_{F}\oplus\mathbf{\bar{5}}_{F})\oplus(\mathbf{50}_{H}\oplus\overline{\mathbf{50}}_{H}\oplus\mathbf{5}_{H}\oplus\overline{\mathbf{5}}_{H})\oplus\mathbf{75}_{H}\oplus\mathbf{24}_{H}.$
(A.12)

This yields $\beta_{f}=-57$ (see [59, 60, 61] for tables of Dynkin indices),
and assuming $g(\mu_{0})\approx 0.7$ and $\Delta\beta=-35$ (we move the
$\mathbf{50}_{H}\oplus\overline{\mathbf{50}}$ a factor $c$ below
$M_{\text{GUT}}$), we obtain from Eq. (A.11) the approximate value

$\displaystyle c\approx 10^{2}/C^{1.63}.$ (A.13)

Even for a nearby Landau pole a factor only $C=10$ above $M_{\text{GUT}}$, we
can move the $\mathbf{50}$s only a factor $c\approx 2.35$ below the GUT scale.

The above calculation implies that to avoid a Landau pole too close to
$M_{\text{GUT}}$, the $\mathbf{50}$ cannot really be much lighter than the GUT
scale. A small $m_{50}$ alone thus cannot be used to enhance
$M^{\text{eff}}_{T}$ by a large factor in the usual missing partner mechanism.

To avoid problems with perturbativity, we need to modify the missing partner
mechanism so that the $2$-$2$ element of $\mathbf{M}_{T}$ has a richer
structure. We can for example introduce representations
$\mathbf{24}^{\prime}_{H}$ and $\mathbf{75}^{\prime}_{H}$ (possibly those that
we already have of this dimension, depends on the full model), and two new
terms to the superpotential $W_{DT}$ in Eq. (A.5):

$\displaystyle W^{\prime}_{DT}$
$\displaystyle=W_{DT}+\lambda_{24}\;\mathbf{50}_{H}\cdot\langle\mathbf{24}^{\prime}_{H}\rangle\cdot\mathbf{\overline{50}}_{H}+\lambda_{75}\;\mathbf{50}_{H}\cdot\langle\mathbf{75}^{\prime}_{H}\rangle\cdot\mathbf{\overline{50}}_{H}.$
(A.14)

This generates the triplet mass matrix of the new form

$\displaystyle\mathbf{M}_{T}$ $\displaystyle=\begin{pmatrix}0&\lambda v\\\
\bar{\lambda}v&\tilde{m}_{50}\\\ \end{pmatrix},$ (A.15)

where

$\displaystyle\tilde{m}_{50}$
$\displaystyle=m_{50}+\lambda_{24}\langle\mathbf{24}^{\prime}_{H}\rangle+\lambda_{75}\langle\mathbf{75}^{\prime}_{H}\rangle.$
(A.16)

The expression for $\tilde{m}_{50}$ can now have all terms at the GUT scale,
but the couplings $\lambda_{24}$ and $\lambda_{75}$ can be tuned such that
$\tilde{m}_{50}$ itself for the triplets is small, while keeping all SM
representations in the $\mathbf{50}_{H}\oplus\mathbf{\overline{50}}_{H}$ other
than the triplet heavy. This can be seen from the (relative) CG factors in
front of the $3$ operators generating their masses that we computed explicitly
in Table 5. We can see from this table that the CG factor pairs for the
operators $(\lambda_{24},\lambda_{75})$ in each row are unique; since the mass
of any SM irreducible representation in the
$\mathbf{50}_{H}\oplus\mathbf{\overline{50}}_{H}$ has the form of Eq. (A.16),
but with CG factors added to the $\lambda_{24}$ and $\lambda_{75}$ terms, the
uniqueness of the CG pairs implies that any one SM irreducible representations
can be tuned to be light while keeping the others heavy; we choose to tune it
to a low value for the triplets.

In this modified missing partner mechanism, we achieved what we wanted: it is
possible to perform DT-splitting, as well as tune $M^{\text{eff}}_{T}$ to
arbitrarily high values. All physical masses of particles in the Higgs sector
are around the GUT scale (except the MSSM doublet pair), while precision of
the tuning in the parameter $\tilde{m}_{50}$ controls how high the effective
triplet mass $M^{\text{eff}}_{T}$ is:

$\displaystyle M^{\text{eff}}_{T}$
$\displaystyle=-\lambda\bar{\lambda}\,v^{2}/\tilde{m}_{50}.$ (A.17)

For $\lambda\bar{\lambda}v^{2}\approx(10^{16}\,\mathrm{GeV})^{2}$, and
effective triplet mass $|M^{\text{eff}}_{T}|=10^{19}\,\mathrm{GeV}$ used as a
benchmark value in this paper, we require $\tilde{m}_{50}\approx
10^{-3}M_{\text{GUT}}$, i.e. a relatively small tuning of one part in
$10^{3}$, which retains the advantage of the missing partner mechanism over a
brute force fine-tuning of a doublet from $M_{\text{GUT}}$ to the EW scale.

Table 5: Clebsch factors (up to an overall VEV normalization factor) for all different SM representations in the $\mathbf{50}$ (or $\overline{\mathbf{50}}$) of $\mathrm{SU}(5)$ in the three different operators generating their masses. $G_{321}$ representation | $\mathbf{50}\cdot\overline{\mathbf{50}}$ | $\mathbf{50}\cdot\mathbf{24}\cdot\overline{\mathbf{50}}$ | $\mathbf{50}\cdot\mathbf{75}\cdot\overline{\mathbf{50}}$
---|---|---|---
$(\mathbf{1},\mathbf{1},-2)$ | $1$ | $-12$ | $-3$
$(\mathbf{3},\mathbf{1},-\frac{1}{3})$ | $1$ | $-2$ | $-1$
$(\mathbf{\bar{3}},\mathbf{2},-\tfrac{7}{6})$ | $1$ | $-7$ | $-1$
$(\mathbf{\bar{6}},\mathbf{3},-\tfrac{1}{3})$ | $1$ | $-2$ | $1$
$(\mathbf{6},\mathbf{1},+\tfrac{4}{3})$ | $1$ | $8$ | $-1$
$(\mathbf{8},\mathbf{2},+\tfrac{1}{2})$ | $1$ | $3$ | $0$

###  A.2 Double missing partner mechanism

A straightforward alternative to the modified missing partner mechanism
discussed in Appendix A.1 is the double missing partner mechanism, cf. [21].
This also solves the DT splitting problem and provides only one effective
triplet mass that can be trans-Planckian.

In this scenario, the representations of the doublets/triplets are doubled:
the DT sector consists of

$\displaystyle(\mathbf{5}_{H}\oplus\mathbf{\bar{5}}_{H}\oplus\mathbf{50}_{H}\oplus\mathbf{\overline{50}}_{H})\oplus(\mathbf{5}^{\prime}_{H}\oplus\mathbf{\bar{5}}^{\prime}_{H}\oplus\mathbf{50}^{\prime}_{H}\oplus\mathbf{\overline{50}}^{\prime}_{H}),$
(A.18)

where the $G_{321}$ doublets $(\mathbf{1},\mathbf{2},\pm 1/2)$ are contained
in

$\displaystyle D_{1}$ $\displaystyle\subseteq\mathbf{5}_{H},$
$\displaystyle\bar{D}_{1}$ $\displaystyle\subseteq\mathbf{\bar{5}}_{H},$
(A.19) $\displaystyle D^{\prime}_{1}$
$\displaystyle\subseteq\mathbf{5}^{\prime}_{H},$
$\displaystyle\bar{D}^{\prime}_{1}$
$\displaystyle\subseteq\mathbf{\bar{5}}^{\prime}_{H},$ (A.20)

and the triplets $T_{I}\sim(\mathbf{3},\mathbf{1},-1/3)$ and anti-triplets
$\overline{T}_{I}\sim(\mathbf{\bar{3}},\mathbf{1},+1/3)$ are in

$\displaystyle T_{1}$ $\displaystyle\subseteq\mathbf{5}_{H},$
$\displaystyle\overline{T}_{1}$ $\displaystyle\subseteq\mathbf{\bar{5}}_{H},$
$\displaystyle T^{\prime}_{1}$
$\displaystyle\subseteq\mathbf{5}^{\prime}_{H},$
$\displaystyle\overline{T}^{\prime}_{1}$
$\displaystyle\subseteq\mathbf{\bar{5}}^{\prime}_{H},$ $\displaystyle T_{2}$
$\displaystyle\subseteq\mathbf{50}_{H},$ $\displaystyle\overline{T}_{2}$
$\displaystyle\subseteq\mathbf{\overline{50}}_{H},$ $\displaystyle
T^{\prime}_{2}$ $\displaystyle\subseteq\mathbf{50}^{\prime}_{H},$
$\displaystyle\overline{T}^{\prime}_{2}$
$\displaystyle\subseteq\mathbf{\overline{50}}^{\prime}_{H}.$ (A.21)

We then arrange the presence of such operators that the mass terms

$\displaystyle\mathbf{W}\supset\begin{pmatrix}D_{1}\\\ D^{\prime}_{1}\\\
\end{pmatrix}^{\\!{\textsf{T}}}\mathbf{M}_{D}\begin{pmatrix}\overline{D}_{1}\\\
\overline{D}^{\prime}_{1}\\\ \end{pmatrix}+\begin{pmatrix}T_{1}\\\
T^{\prime}_{1}\\\ T_{2}\\\ T^{\prime}_{2}\\\
\end{pmatrix}^{\\!{\textsf{T}}}\mathbf{M}_{T}\begin{pmatrix}\overline{T}_{1}\\\
\overline{T}^{\prime}_{1}\\\ \overline{T}_{2}\\\ \overline{T}^{\prime}_{2}\\\
\end{pmatrix}+\cdots,$ (A.22)

have the following specific form for the mass matrices:

$\displaystyle\mathbf{M}_{D}$ $\displaystyle=\begin{pmatrix}0&0\\\
0&\mu^{\prime}\\\ \end{pmatrix},$ $\displaystyle\mathbf{M}_{T}$
$\displaystyle=\begin{pmatrix}0&0&0&\alpha_{1}V\\\
0&\mu^{\prime}&\alpha_{2}V&0\\\ \alpha_{3}V&0&M_{50}&0\\\
0&\alpha_{4}V&0&M^{\prime}_{50}\\\ \end{pmatrix},$ (A.23)

where $\alpha_{i}$ are dimensionless couplings, while $V$ could for example be
a VEV $V=\langle\mathbf{75}\rangle$ at the GUT scale (see [21] for more
details on some realizations). In this set up, we require that only
$\mathbf{5}_{H}\oplus\overline{\mathbf{5}}_{H}$ (and not those with primes)
can couple to the SM fermions in
$\mathbf{10}_{Fi}\oplus\overline{\mathbf{5}}_{Fi}$, which yields the effective
triplet setup with the effective mass parameter equal to

$\displaystyle
M^{\text{eff}}_{T}=((\mathbf{M}_{T}^{-1})_{11})^{-1}=-\alpha_{1}\alpha_{2}\alpha_{3}\alpha_{4}\frac{V^{4}}{\mu^{\prime}\,M_{50}\,M^{\prime}_{50}}.$
(A.24)

To postpone the appearance of the Landau pole from the large (negative)
contributions to the beta function, we can provide the representations
$\mathbf{50}\oplus\mathbf{50}^{\prime}\oplus\overline{\mathbf{50}}\oplus\overline{\mathbf{50}}^{\prime}$
with large mass terms close to the Planck scale, i.e. we take the hierarchy

$\displaystyle\mu^{\prime}\quad\ll\quad V\approx M_{\text{GUT}}\quad\ll\quad
M_{50}\approx M^{\prime}_{50}\quad\lesssim\quad M_{Pl}.$ (A.25)

We obtain from this setup the physical doublet masses to be $0$ (the MSSM
doublet) and $\mu^{\prime}$, while the physical triplet masses are at scales
$\eta^{\pm 1}V$ (two triplets at each scale), where $\eta\equiv V/M_{50}$.
This constitutes a successful DT splitting yielding

$\displaystyle M^{\text{eff}}_{T}\approx-\eta^{2}\frac{V^{2}}{\mu^{\prime}},$
(A.26)

allowing for the effective mass to be large provided that $\mu^{\prime}$ is
small. For $V=2\cdot 10^{16}\,\mathrm{GeV}$, $\eta^{2}\approx 10^{-4}$ (so
that $M_{50}$ and $M^{\prime}_{50}$ are close to the Planck scale) and
$|M^{\text{eff}}_{T}|\gtrsim 10^{19}\,\mathrm{GeV}$, we need
$|\mu^{\prime}|\lesssim 4\cdot 10^{9}\,\mathrm{GeV}$. Note that $\mu^{\prime}$
is associated to the scale of the second doublet (a rather small
representation not irreparably spoiling the unification of gauge couplings),
but not to any other particle, since the $\mathbf{5}$s contain beside the
doublets only the already considered triplets. For a more precise calculation,
such a scenario does require adding an extra doublet to the MSSM at the
intermediate scale $\mu^{\prime}$.

##  Appendix B Documentation for the proton decay calculation

###  B.1 Numerical procedure

####  B.1.1 Definition of the dimension 5 operators

The effective superpotential of the dimension $5$ operators, which are
relevant for proton and neutron decay, reads

$\displaystyle\begin{split}W_{5}&=-\tfrac{1}{2}\mathbf{C}^{ijkl}_{5L}\,\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}\,(Q_{i}^{\hat{a}}\cdot
L_{j})(Q_{k}^{\hat{b}}\cdot Q_{l}^{\hat{c}})\\\
&\quad-\mathbf{C}^{ijkl}_{5R}\,\epsilon^{{\hat{a}}{\hat{b}}{\hat{c}}}\,u^{c}_{i{\hat{a}}}\,d^{c}_{j{\hat{b}}}\,e^{c}_{k}\,u^{c}_{l{\hat{c}}}\,,\end{split}$
(B.1)

where ${\hat{a}},{\hat{b}},{\hat{c}}$ are $\mathrm{SU}(3)_{\text{C}}$ indices
with $\epsilon_{123}=\epsilon^{123}=1$, and the contraction of
$\mathrm{SU}(2)_{\text{L}}$ indices is given by
$\Psi\cdot\Phi=\epsilon_{ab}\,\Psi^{a}\,\Phi^{b}$ with $\epsilon_{12}=1$.

The quasi-Yukawa operators and the mass term for the (anti-)triplets are given
by

$\displaystyle\begin{split}W_{\text{T}}&=\phantom{+}\tfrac{1}{2}(\tilde{\mathbf{Y}}_{qq})_{Iij}\,\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}\,Q_{i}^{\hat{a}}\cdot
Q_{j}^{\hat{b}}\,T_{I}^{\hat{c}}+(\tilde{\mathbf{Y}}_{eu})_{Iij}\,e^{c}_{i}\,u^{c}_{j{\hat{a}}}\,T_{I}^{\hat{a}}\\\
&\quad+(\tilde{\mathbf{Y}}_{ql})_{Jij}\,Q_{i}^{\hat{a}}\cdot
L_{j}\,\overline{T}_{J{\hat{a}}}+(\tilde{\mathbf{Y}}_{ud})_{Jij}\,\epsilon^{{\hat{a}}{\hat{b}}{\hat{c}}}\,u^{c}_{i{\hat{a}}}\,d^{c}_{j{\hat{b}}}\,\overline{T}_{J{\hat{c}}}\\\
&\quad+(\mathbf{M}_{T})_{IJ}\,T_{I}^{\hat{a}}\,\overline{T}_{J{\hat{a}}}\,.\end{split}$
(B.2)

After integrating out the (anti-)triplets, the effective baryon and lepton
number violating operators are obtained:

$\displaystyle\begin{split}W_{5}&=-\tfrac{1}{2}(\mathbf{M}_{T}^{-1})_{IJ}\,(\tilde{\mathbf{Y}}_{ql})_{Iij}\,(\tilde{\mathbf{Y}}_{qq})_{Jkl}\,\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}\,(Q_{i}^{\hat{a}}\cdot
L_{j})(Q_{k}^{\hat{b}}\cdot Q_{l}^{\hat{c}})\\\
&\quad-(\mathbf{M}_{T}^{-1})_{IJ}\,(\tilde{\mathbf{Y}}_{ud})_{Iij}\,(\tilde{\mathbf{Y}}_{eu})_{Jkl}\,\epsilon^{{\hat{a}}{\hat{b}}{\hat{c}}}\,u^{c}_{i{\hat{a}}}\,d^{c}_{j{\hat{b}}}\,e^{c}_{k}\,u^{c}_{l{\hat{c}}}\\\
&\quad+...\,,\end{split}$ (B.3)

where the dots indicate terms which are not relevant for nucleon decay.
According to Eq. (B.1), the dimension $5$ operators at the GUT scale are then
calculated as follows:

$\displaystyle\mathbf{C}^{ijkl}_{5L}$
$\displaystyle=(\mathbf{M}_{T}^{-1})_{IJ}\,(\tilde{\mathbf{Y}}_{ql})_{Iij}\,(\tilde{\mathbf{Y}}_{qq})_{Jkl}\,,$
(B.4) $\displaystyle\mathbf{C}^{ijkl}_{5R}$
$\displaystyle=(\mathbf{M}_{T}^{-1})_{IJ}\,(\tilde{\mathbf{Y}}_{ud})_{Iij}\,(\tilde{\mathbf{Y}}_{eu})_{Jkl}\,.$
(B.5)

####  B.1.2 Running of the dimension 5 operators

The running of the dimension $5$ operators from the GUT scale to the SUSY
scale is calculated in the framework of the MSSM. The $\beta$-functions of the
dimension $5$ operators in the MSSM at $1$-loop in the
$\overline{\mathrm{DR}}$ scheme, by using the left-right convention for the
Yukawa matrices, are given by

$\displaystyle\begin{split}\beta(\mathbf{C}_{5L}^{ijkl})&=\phantom{+}\mathbf{C}_{5L}^{mjkl}\big{(}\mathbf{Y}_{d}\mathbf{Y}_{d}^{\dagger}+\mathbf{Y}_{u}\mathbf{Y}_{u}^{\dagger}\big{)}^{i}{}_{m}+\mathbf{C}_{5L}^{imkl}\big{(}\mathbf{Y}_{e}\mathbf{Y}_{e}^{\dagger}+\mathbf{Y}_{\nu}\mathbf{Y}_{\nu}^{\dagger}\big{)}^{j}{}_{m}\\\
&\quad+\mathbf{C}_{5L}^{ijml}\big{(}\mathbf{Y}_{d}\mathbf{Y}_{d}^{\dagger}+\mathbf{Y}_{u}\mathbf{Y}_{u}^{\dagger}\big{)}^{k}{}_{m}+\mathbf{C}_{5L}^{ijkm}\big{(}\mathbf{Y}_{d}\mathbf{Y}_{d}^{\dagger}+\mathbf{Y}_{u}\mathbf{Y}_{u}^{\dagger}\big{)}^{l}{}_{m}\\\
&\quad+\mathbf{C}_{5L}^{ijkl}\big{(}-\tfrac{2}{5}g_{1}^{2}-6g_{2}^{2}-8g_{3}^{2}\big{)}\,,\end{split}$
(B.6)
$\displaystyle\begin{split}\beta(\mathbf{C}_{5R}^{ijkl})&=\phantom{+}\mathbf{C}_{5R}^{mjkl}\big{(}2\mathbf{Y}_{u}^{\dagger}\mathbf{Y}_{u}\big{)}_{m}{}^{i}+\mathbf{C}_{5R}^{imkl}\big{(}2\mathbf{Y}_{d}^{\dagger}\mathbf{Y}_{d}\big{)}_{m}{}^{j}\\\
&\quad+\mathbf{C}_{5R}^{ijml}\big{(}2\mathbf{Y}_{e}^{\dagger}\mathbf{Y}_{e}\big{)}_{m}{}^{k}+\mathbf{C}_{5R}^{ijkm}\big{(}2\mathbf{Y}_{u}^{\dagger}\mathbf{Y}_{u}\big{)}_{m}{}^{l}\\\
&\quad+\mathbf{C}_{5R}^{ijkl}\big{(}-\tfrac{12}{5}g_{1}^{2}-8g_{3}^{2}\big{)}\,.\end{split}$
(B.7)

The $\beta$-functions are computed by means of the general formulas in [62]
for the $1$-loop anomalous dimension matrix ${\gamma_{m}}^{i}$, and the non-
renormalization theorem of the superpotential [63, 64, 65]. We use the GUT
normalization for the $\mathrm{U}(1)_{\text{Y}}$ gauge coupling $g_{1}$; the
SM normalization can be recovered by $g_{1}^{\text{SM}}=\sqrt{3/5}\,g_{1}$.
The RGEs are then given by

$\displaystyle\frac{\mathrm{d}x}{\mathrm{d}\log\mu}=\frac{1}{16\pi^{2}}\beta(x)\,,$
(B.8)

where $x\in\\{\mathbf{C}_{5L}^{ijkl},\mathbf{C}_{5R}^{ijkl}\\}$ and $\mu$ is
the renormalization scale.

####  B.1.3 Dressing of the dimension 5 operators

The dimension $6$ four-fermion operators are built at $1$-loop level by
dressing the dimension $5$ operators with gluino, chargino and neutralino
exchange diagrams, which takes place at the SUSY scale. An example diagram is
shown in Figure 6. The complete set of formulas to calculate the dressing of
the dimension $5$ operators is stated in [12], in particular Eqs.
(A.1.13)–(A.1.16) and Eqs. (A.2.8)–(A.2.20). Furthermore, in [66] all relevant
diagrams for nucleon decay from dimension $5$ operators can be found.

Figure 6: An example diagram for the proton decay channel $p\rightarrow
K^{+}\bar{\nu}$. The big black circle represents the effective dimension $5$
operator, and the dressing of this operator takes place by using a mass
insertion of the winos $\widetilde{W}^{\pm}$.

Beside the dimension $5$ operators, the calculation of the resulting dimension
$6$ operators involves the rotation matrices of the sfermions, charginos and
neutralinos, as well as the mass eigenvalues of the sparticles. It is
convenient to use the standardized SLHA [15] and SLHA2 [67] conventions to
specify these quantities. In Section B.2 the translation of the conventions
used in [12] into the SLHA conventions is stated.

Apart from integrating out the sparticles at the SUSY scale, we also switch to
the EW symmetry broken phase of the SM. In the following $u_{i}$, $d_{i}$ and
$e_{i}$ $(i\in\\{1,2,3\\})$ represent the mass eigenstates of the up-type
quarks, the down-type quarks and the charged leptons. In contrast, the
$\nu_{i}$ represent the neutrinos in the interaction (flavor) basis. The
additional label $L$ or $R$ indicates whether a state is a left- or a right-
handed Weyl spinor. The dimension $6$ operators $\widetilde{\mathbf{C}}$ (at
the Lagrangian level), which are relevant for the calculation of proton decay,
include only mass eigenstates of the fermions which are lighter than the
nucleons and are given by

$\displaystyle\begin{split}\mathcal{L}_{\not{B}}&=\frac{1}{16\pi^{2}}\epsilon_{{\hat{a}}{\hat{b}}{\hat{c}}}\Big{(}\widetilde{\mathbf{C}}_{LL}(udue)^{ik}(u_{L}^{\hat{a}}d_{Li}^{\hat{b}})(u_{L}^{\hat{c}}e_{Lk}^{\phantom{{\hat{b}}}})+\widetilde{\mathbf{C}}_{RL}(udue)^{ik}(u_{R}^{\hat{a}}d_{Ri}^{\hat{b}})(u_{L}^{\hat{c}}e_{Lk}^{\phantom{{\hat{b}}}})\\\
&\qquad\qquad\;\;+\widetilde{\mathbf{C}}_{LR}(udue)^{ik}(u_{L}^{\hat{a}}d_{Li}^{\hat{b}})(u_{R}^{\hat{c}}e_{Rk}^{\phantom{{\hat{b}}}})+\widetilde{\mathbf{C}}_{RR}(udue)^{ik}(u_{R}^{\hat{a}}d_{Ri}^{\hat{b}})(u_{R}^{\hat{c}}e_{Rk}^{\phantom{{\hat{b}}}})\\\
&\qquad\qquad\;\;+\widetilde{\mathbf{C}}_{LL}(udd\nu)^{ijk}(u_{L}^{\hat{a}}d_{Li}^{\hat{b}})(d_{Lj}^{\hat{c}}\nu_{Lk}^{\phantom{{\hat{b}}}})+\widetilde{\mathbf{C}}_{RL}(udd\nu)^{ijk}(u_{R}^{\hat{a}}d_{Ri}^{\hat{b}})(d_{Lj}^{\hat{c}}\nu_{Lk}^{\phantom{{\hat{b}}}})\\\
&\qquad\qquad\;\;+\tfrac{1}{2}\widetilde{\mathbf{C}}_{RL}(ddu\nu)^{ijk}(d_{Ri}^{\hat{a}}d_{Rj}^{\hat{b}})(u_{L}^{\hat{c}}\nu_{Lk}^{\phantom{{\hat{b}}}})\Big{)}+h.c.\,.\end{split}$
(B.9)

Since only the lightest up-type quark $u_{1}$ is lighter than the proton and
the neutron, only this state is present in Eq. (B.9) and the label $i$ is
neglected.999Note that in the calculation of the decay widths, the index $i$
in $d_{i}$ and $e_{i}$ is restricted to the values $1,2$, because the heaviest
down-type quark and charged lepton $(i=3)$ are heavier than the two nucleons.
In contrast, all three states $\nu_{i}$ are present.

####  B.1.4 Running of the dimension 6 operators

The running of the dimension $6$ operators from the SUSY scale to the mass
scale of the nucleons is calculated in the framework of the SM. Following
[68], the main contribution to the running of the dimension $6$ operators
comes from the strong gauge coupling $g_{3}$. At $1$-loop the
$\beta$-functions in the $\overline{\mathrm{MS}}$ scheme of the dimension $6$
operators $\widetilde{\mathbf{C}}$ from Eq. (B.9), written in the general form
$\widetilde{\mathbf{C}}^{ijkl}$, and of the gauge coupling $g_{3}$ are given
by

$\displaystyle\beta(\widetilde{\mathbf{C}}^{ijkl})$
$\displaystyle=-4\,g_{3}^{2}\,\widetilde{\mathbf{C}}^{ijkl}\,,$ (B.10)
$\displaystyle\beta(g_{3})$
$\displaystyle=\Big{(}-11+\frac{2}{3}N_{F}\Big{)}g_{3}^{3}\,,$ (B.11)

where $N_{F}$ is the number of quarks whose masses are below the
renormalization scale $\mu$. The RGEs have the form as specified in Eq. (B.8).
The running of the dimension $6$ operators just corresponds to an overall
scaling factor, which is also referred to as the long range effect on the
effective operators $\widetilde{\mathbf{C}}$.

####  B.1.5 Partial decay widths

The formula for the calculation of the partial decay widths of a nucleon
$B_{i}$ with a meson $M_{j}$ and a lepton $l_{k}$ in the final state is taken
from [12] and is given by

$\displaystyle\Gamma(B_{i}\rightarrow M_{j}\,l_{k})$
$\displaystyle=\frac{m_{i}}{32\pi}\,\Big{(}1-\frac{m^{2}_{j}}{m^{2}_{i}}\Big{)}^{2}\,\frac{1}{f_{\pi}^{2}}\,\Big{(}\big{|}A_{L}^{ijk}\big{|}^{2}+\big{|}A_{R}^{ijk}\big{|}^{2}\Big{)}\,,$
(B.12)

where $m_{i}$ and $m_{j}$ are the masses of the nucleon and meson,
respectively. Moreover, $f_{\pi}=0.130\,\mathrm{GeV}$ (cf. PDG [37]) is the
pion decay constant.

The amplitudes $A_{L}^{ijk}$ and $A_{R}^{ijk}$ for different decay channels
are calculated by using the dimension $6$ operators at the nucleon mass scale,
and are specified in Table 1 in [12]. In that table, the nucleon mass $m_{N}$
is taken as the average of the proton mass $m_{p}$ and the neutron mass
$m_{n}$, and the baryon mass $m_{B^{\prime}}$ is given by the average of
$m_{\Lambda}$ and $m_{\Sigma}$. The additional constants with the
corresponding (approximate) values are the constants from hyperon decay
$F\approx 0.463$, $D\approx 0.804$ (cf. [69, 7]), and the parameters for the
proton decay matrix element $\alpha\approx 0.0090\,\mathrm{GeV}^{3}$,
$\beta\approx 0.0096\,\mathrm{GeV}^{3}$ (cf. [70, 7]).

The following decay channels for the proton $p$ and the neutron $n$ are
considered in the calculation of the partial decay widths:

$\displaystyle\bullet\text{ proton: }$ $\displaystyle
p\rightarrow\pi^{0}\,e_{k}^{+},\quad p\rightarrow\eta^{0}\,e_{k}^{+},\quad
p\rightarrow K^{0}\,e_{k}^{+},\quad p\rightarrow\pi^{+}\,\bar{\nu}_{k},\quad
p\rightarrow K^{+}\,\bar{\nu}_{k},$ (B.13) $\displaystyle\bullet\text{
neutron: }$ $\displaystyle n\rightarrow\pi^{-}\,e_{k}^{+},\quad
n\rightarrow\pi^{0}\,\bar{\nu}_{k},\quad
n\rightarrow\eta^{0}\,\bar{\nu}_{k},\quad n\rightarrow K^{0}\,\bar{\nu}_{k}.$

As discussed above, $k\in\\{1,2\\}$ for $e_{k}^{+}$, since the heaviest
charged lepton $\tau$ corresponding to $k=3$ is heavier than the nucleons, and
$k\in\\{1,2,3\\}$ for $\bar{\nu}_{k}$.

###  B.2 Matching with SLHA conventions

In the following the relation between the conventions in [12] (later referred
to as GN), which are used in the dressing of the dimension $5$ operators, and
the standardized set of conventions in SLHA [15] and SLHA2 [67] are specified.

  * •

The Yukawa matrices are related as follows:

GN | SLHA | our notation
---|---|---
$f_{U}$ | $Y_{U}$ | = | $f_{U}$
$f_{D}$ | $Y_{D}$ | = | $-f_{D}$
$f_{L}$ | $Y_{E}$ | = | $f_{L}^{\textsf{T}}$

where the Yukawa matrices are given in GN: Eq. (A.1.1) and SLHA: Eq. (3). In
particular, GN uses the LR convention in the up- and down-quark sector
($f_{U}$ and $f_{D}$), and the right-left (RL) convention in the charged
lepton sector $(f_{L})$. On the other hand, SLHA uses the LR convention in all
sectors.

  * •

Both GN and SLHA use the SM normalization for the gauge coupling $g_{1}$ of
$\mathrm{U}(1)_{\text{Y}}$, and they also use the same convention for the
$\mu$-term.

  * •

In the EW symmetry broken phase the mass matrices of the sfermions are written
in a basis where the sfermions are aligned with their SM superpartners.

    * –

GN uses the basis where the mass matrices $M_{U}=M_{U}^{\text{diag}}$ (up-
quarks) and $M_{L}=M_{L}^{\text{diag}}$ (charged leptons) are diagonal, and
the mass matrices $M_{D}=V_{\text{CKM}}^{*}\,M_{D}^{\text{diag}}$ (down-
quarks) and
$M_{\nu}=V_{\text{PMNS}}^{*}\,M_{\nu}^{\text{diag}}\,V_{\text{PMNS}}^{\textsf{T}}$
(neutrinos) are rotated by the CKM and PMNS matrix, respectively ($u_{L}$ and
$d_{L}$ as well as $e_{L}$ and $\nu_{L}$ are assumed to form doublets under
rotations).

    * –

SLHA uses the basis where all fermion mass matrices are diagonal, which
corresponds to the super-CKM/PMNS basis.

The relation between the sfermion, the chargino and the neutralino mass
matrices in GN and SLHA are then the following:

GN | SLHA |  |  |
---|---|---|---|---
$\mathcal{M}^{2}_{\widetilde{u}}$ | $(\mathcal{M}^{2}_{\widetilde{u}})_{\text{sCKM}}$ | $=$ | $(\mathcal{M}^{2}_{\widetilde{u}})^{\textsf{T}}$ |
$\mathcal{M}^{2}_{\widetilde{d}}$ | $(\mathcal{M}^{2}_{\widetilde{d}})_{\text{sCKM}}$ | $=$ | $U_{d}(\mathcal{M}^{2}_{\widetilde{d}})^{\textsf{T}}U_{d}^{\dagger}$ | with $U_{d}=\left(\begin{smallmatrix}V^{\dagger}_{\text{CKM}}&0\\\ 0&\mathbb{1}\end{smallmatrix}\right)$
$\mathcal{M}^{2}_{\widetilde{e}}$ | $(\mathcal{M}^{2}_{\widetilde{e}})_{\text{sPMNS}}$ | $=$ | $\mathcal{M}^{2}_{\widetilde{e}}$ |
$\mathcal{M}^{2}_{\widetilde{\nu}}$ | $(\mathcal{M}^{2}_{\widetilde{\nu}})_{\text{sPMNS}}$ | $=$ | $U_{\nu}(\mathcal{M}^{2}_{\widetilde{\nu}})^{\textsf{T}}U_{\nu}^{\dagger}$ | with $U_{\nu}=V_{\text{PMNS}}^{\dagger}$
$\mathcal{M}_{C}$ | $\mathcal{M}_{\widetilde{\psi}^{+}}$ | $=$ | $-\mathcal{M}_{C}$ |
$\mathcal{M}_{N}$ | $\mathcal{M}_{\widetilde{\psi}^{0}}$ | $=$ | $\mathcal{M}_{N}$ |

The mass matrices of the sfermions are given in GN: Eqs. (A.1.6)–(A.1.9) and
SLHA2: Eq. (11), (12), (24) and (25). In addition, the mass matrices of the
charginos and neutralinos are stated in GN: Eq. (A.1.11) and SLHA: Eq. (21)
and (22).

The mass matrices of the squarks and charged sleptons $\widetilde{f}$
$(f\in\\{u,d,e\\})$ are $6\times 6$-dimensional, written in the basis
$(\widetilde{f}_{L1},\widetilde{f}_{L2},\widetilde{f}_{L3},\widetilde{f}_{R1},\widetilde{f}_{R2},\widetilde{f}_{R3})$,
where $\widetilde{f}_{L,Ri}$ is the superpartner of $f_{L,Ri}$, i.e.

$\displaystyle\mathcal{M}^{2}_{\widetilde{f}}$
$\displaystyle=\begin{pmatrix}\mathcal{M}^{2}_{LL}&\mathcal{M}^{2}_{LR}\\\
\mathcal{M}^{2}_{RL}&\mathcal{M}^{2}_{RR}\end{pmatrix},$ (B.14)

where each entry represents a $3\times 3$-block. In contrast, the mass matrix
$\mathcal{M}^{2}_{\widetilde{\nu}}\equiv\mathcal{M}^{2}_{LL}$ of the
sneutrinos is only $3\times 3$-dimensional, since only left-handed neutrinos
are present.

Moreover, the $2\times 2$-dimensional mass matrix $\mathcal{M}_{C}$ of the
charginos is written in the basis
$(\widetilde{W}^{\pm},\widetilde{H}^{\pm}_{u})$, and for the $4\times
4$-dimensional mass matrix $\mathcal{M}_{N}$ of the neutralinos the basis
$(\widetilde{B},\widetilde{W}^{0},\widetilde{H}^{0}_{d},\widetilde{H}^{0}_{u})$
is used.

  * •

The definitions and relations of the sparticle rotation matrices in GN and
SLHA are given by:

GN | SLHA |  | Relation
---|---|---|---
$\widetilde{U}_{U}(\mathcal{M}^{2}_{\widetilde{u}})^{\textsf{T}}\widetilde{U}_{U}^{\dagger}=\text{diag}$ | $R_{u}(\mathcal{M}^{2}_{\widetilde{u}})_{\text{sCKM}}R_{u}^{\dagger}=\text{diag}$ | $\Rightarrow$ | $\widetilde{U}_{U}=R_{u}$
$\widetilde{U}_{D}(\mathcal{M}^{2}_{\widetilde{d}})^{\textsf{T}}\widetilde{U}_{D}^{\dagger}=\text{diag}$ | $R_{d}(\mathcal{M}^{2}_{\widetilde{d}})_{\text{sCKM}}R_{d}^{\dagger}=\text{diag}$ | $\Rightarrow$ | $\widetilde{U}_{D}=R_{d}U_{d}$
$\widetilde{U}_{L}^{\dagger}\mathcal{M}^{2}_{\widetilde{e}}\widetilde{U}_{L}=\text{diag}$ | $R_{e}(\mathcal{M}^{2}_{\widetilde{e}})_{\text{sPMNS}}R_{e}^{\dagger}=\text{diag}$ | $\Rightarrow$ | $\widetilde{U}_{L}=R_{e}^{\dagger}$
$\widetilde{U}_{N}^{\dagger}\mathcal{M}^{2}_{\widetilde{\nu}}\widetilde{U}_{N}=\text{diag}$ | $R_{\nu}(\mathcal{M}^{2}_{\widetilde{\nu}})_{\text{sPMNS}}R_{\nu}^{\dagger}=\text{diag}$ | $\Rightarrow$ | $\widetilde{U}_{N}=(R_{\nu}U_{\nu})^{\dagger}$
$-U_{-}^{\dagger}\mathcal{M}_{C}U_{+}=\text{diag}$ | $U(\mathcal{M}_{C})_{\text{SLHA}}V^{\textsf{T}}=\text{diag}$ | $\Rightarrow$ |  $U_{-}=U^{\dagger}$, $U_{+}=V^{\textsf{T}}$
$U_{N}^{\textsf{T}}\mathcal{M}_{N}U_{N}=\text{diag}$ | $N^{*}(\mathcal{M}_{N})_{\text{SLHA}}N^{\dagger}=\text{diag}$ | $\Rightarrow$ | $U_{N}=N^{\dagger}$

where “$\operatorname{diag}$” indicates a diagonal matrix with real, positive
entries. The definitions for the sfermion rotation matrices are taken from GN:
Eq. (A.1.10) and SLHA2: Eqs. (28)–(31), whereas the definitions for the
chargino and neutralino rotation matrices are stated in GN: Eq. (A.1.12) and
SLHA: Eq. (12) and (15).

##  Appendix C User’s guide to Mathematica packages SusyTCProton and
ProtonDecay

This section provides a guide for the usage of the Mathematica packages
SusyTCProton.m and ProtonDecay.m. Visit the web page
http://particlesandcosmology.unibas.ch/downloads/protondecay.html to obtain
the packages. The zip file ProtonDecay.zip there contains the README.md file
with installation instructions, the .m package file and a directory containing
some example notebooks for evaluation. For SusyTCProton, follow the link on
the web page specified above. Mathematica package requirements: REAP for
SusyTCProton, none for ProtonDecay.

The documentation below uses the typewriter font to specify functions and
variables. In particular, the dimension $5$ operators $\mathbf{C}^{ijkl}_{5L}$
and $\mathbf{C}^{ijkl}_{5R}$ from Eq. (B.1) are labeled as CL5qlqq and
CR5udeu, respectively.

###  C.1 SusyTCProton.m

The package SusyTCProton.m has the same functionality as the package SusyTC.m
from [13], with the following extensions/modifications:

  * •

RGEGetSolution:
This function has the same functionality as described in [13], where in
addition the values of the dimension $5$ operators RGECL5qlqq and RGECR5udeu
can be output at a given energy scale.

  * •

STCGetProtonDecayOperators[]:
Returns a list of replacement rules for the dimension $5$ operators in the
flavour basis and in the SCKM basis (flavor indices in all sectors are rotated
so that the Yukawa matrices are diagonal, see e.g. [13] for SCKM definition):

    * –

"FlavourBasis[CL5qlqq]"

    * –

"FlavourBasis[CR5udeu]"

    * –

"SCKMBasis[CL5qlqq]"

    * –

"SCKMBasis[CR5udeu]"

  * •

STCGetSUSYSpectrum[], STCGetSCKMValues[], STCGetOneLoopValues[]:
These functions are the same as described in [13].

  * •

STCGetInternalValues[]:
Returns a list of replacement rules of the quantities which are internally
used for the calculation of the threshold corrections and the sparticle
spectrum, i.e. the combined results from STCGetSCKMValues[],
STCGetSUSYSpectrum[] and STCGetProtonDecayOperators[], with $1$-loop corrected
parameters from STCGetOneLoopValues[] replacing tree-level ones if available.

###  C.2 ProtonDecay.m

The package ProtonDecay.m contains functions to calculate $D=5$ and $D=6$
nucleon decay according to the considerations in Section B.

  * •

RunDim5OperatorsToMSUSY[{Y,MR,g,C5},scales]:
Calculates the running of the dimension $5$ operators from the GUT scale to
the SUSY scale at $1$-loop, according to Eq. (B.6) and (B.7), and the use of
beta functions of the MSSM superpotential parameters therein (see e.g. [13]
for their RGE). The input is specified at the GUT scale and has the following
form:

    * –

Y = {Yu,Yd,Ye,Ynu}: Yukawa matrices using left-right (SLHA) convention,

    * –

MR: right-handed neutrino mass matrix,

    * –

g = {g1,g2,g3}: gauge couplings (GUT normalization for g1 is used),

    * –

C5 = {C5Lqlqq,C5Rudeu}: dimension $5$ operators in the SCKM basis,

    * –

scales = {MGUT,MSUSY}: GUT scale and SUSY scale.

The output is given at the SUSY scale and has the form {Y,MR,g,C5}, where the
same definitions as for the input are used.

  * •

ProtonDecayWidths:
Calculates the partial decay widths of the proton and neutron for the decay
channels listed in Eq. (B.13). To calculate $D=5$ nucleon decay (optionally
combined with $D=6$ nucleon decay), the function accepts two different inputs
at the SUSY scale, where the input of the dimension $6$ operators

C6 = {CLLudue,CRLudue,CLRudue,CRRudue,CLLuddn,CRLuddn,CRLddun} (C.1)

is optional (cf. Eq. (B.9) for the definitions, where $\widetilde{\mathbf{C}}$
are tensors of dimension $3\times 3$ or $3\times 3\times 3$, i.e.
$i,(j),k\in\\{1,2,3\\}$):

    * –

ProtonDecayWidths[STCInternalValues,C6(optional)], where

      * $\bullet$

STCInternalValues: output of the function STCGetInternalValues[] from the
package SusyTCProton,

      * $\bullet$

C6: dimension $6$ operators at $M_{\text{SUSY}}$.

    * –

ProtonDecayWidths[pathSLHA2File,C5,C6(optional)], where

      * $\bullet$

pathSLHA2File: path to a SLHA2 output file,

      * $\bullet$

C5 = {C5Lqlqq,C5Rudeu}: dimension $5$ operators,

      * $\bullet$

C6: dimension $6$ operators at $M_{\text{SUSY}}$.

If only $D=6$ nucleon decay is calculated, the function accepts the following
input at an input mass scale $M_{\text{In}}$ (typically the EW scale):

    * –

ProtonDecayWidths[MIn,g3,C6], where

      * $\bullet$

MIn: input mass scale $M_{\text{In}}$,

      * $\bullet$

g3: gauge coupling of $\mathrm{SU}(3)_{\text{C}}$ at $M_{\text{In}}$,

      * $\bullet$

C6: dimension $6$ operators at $M_{\text{In}}$.

The function has several option values which appear in the calculation of the
partial decay widths (see Sec. B.1.5):

    * –

Options[ProtonDecayWidths], where

      * $\bullet$

"f$\pi$" $\rightarrow$ 0.13: pion decay constant (in GeV),

      * $\bullet$

"F" $\rightarrow$ 0.463: constant from hyperon decay,

      * $\bullet$

"D" $\rightarrow$ 0.804: —"—,

      * $\bullet$

"$\alpha$p" $\rightarrow$ 0.0090: parameter for the proton decay matrix
element,

      * $\bullet$

"$\beta$p" $\rightarrow$ 0.0096: —"—.

The output of the function consists of a list with elements of the form
{label,decayWidths}, where label specifies the decay channel (e.g.
$p\rightarrow\pi^{0}\,e^{+}_{k}$). Furthermore, decayWidths is the list of the
corresponding decay widths in GeV at the nucleon mass scale, where
$k\in\\{1,2\\}$ if the lepton in the final state is a charged lepton, and
$k\in\\{1,2,3\\}$ if the lepton is a neutrino. Note that the charged leptons
are considered in the mass eigenbasis, whereas the neutrinos are considered in
the interaction (flavor) basis.

## References

  * [1] H. Georgi and S. L. Glashow, Phys. Rev. Lett. 32 (1974), 438-441
  * [2] H. Fritzsch and P. Minkowski, Annals Phys. 93 (1975), 193-266
  * [3] H. Georgi, AIP Conf. Proc. 23 (1975), 575-582
  * [4] N. Sakai and T. Yanagida, Nucl. Phys. B 197 (1982), 533
  * [5] S. Weinberg, Phys. Rev. D 26 (1982), 287
  * [6] G. Senjanovic, AIP Conf. Proc. 1200 (2010) no.1, 131-141 [arXiv:0912.5375 [hep-ph]].
  * [7] P. Nath and P. Fileviez Perez, Phys. Rept. 441 (2007), 191-317 [arXiv:hep-ph/0601023 [hep-ph]].
  * [8] S. Mine [Super Kamiokande], J. Phys. Conf. Ser. 718 (2016) no.6, 062044
  * [9] R. Acciarri et al. [DUNE], [arXiv:1512.06148 [physics.ins-det]].
  * [10] K. Abe et al. [Hyper-Kamiokande], [arXiv:1805.04163 [physics.ins-det]].
  * [11] F. An et al. [JUNO], J. Phys. G 43 (2016) no.3, 030401 [arXiv:1507.05613 [physics.ins-det]].
  * [12] T. Goto and T. Nihei, Phys. Rev. D 59 (1999), 115009 [arXiv:hep-ph/9808255 [hep-ph]].
  * [13] S. Antusch and C. Sluka, JHEP 07 (2016), 108 [arXiv:1512.06727 [hep-ph]].
  * [14] S. Antusch, J. Kersten, M. Lindner, M. Ratz and M. A. Schmidt, JHEP 03 (2005), 024 [arXiv:hep-ph/0501272 [hep-ph]].
  * [15] P. Z. Skands, B. C. Allanach, H. Baer, C. Balazs, G. Belanger, F. Boudjema, A. Djouadi, R. Godbole, J. Guasch and S. Heinemeyer, et al. JHEP 07 (2004), 036 [arXiv:hep-ph/0311123 [hep-ph]].
  * [16] W. Porod, Comput. Phys. Commun. 153 (2003), 275-315 [arXiv:hep-ph/0301101 [hep-ph]].
  * [17] W. Porod and F. Staub, Comput. Phys. Commun. 183 (2012), 2458-2469 [arXiv:1104.1573 [hep-ph]].
  * [18] B. C. Allanach, Comput. Phys. Commun. 143 (2002), 305-331 [arXiv:hep-ph/0104145 [hep-ph]].
  * [19] D. Chowdhury, R. Garani and S. K. Vempati, Comput. Phys. Commun. 184 (2013), 899-918 [arXiv:1109.3551 [hep-ph]].
  * [20] A. Djouadi, J. L. Kneur and G. Moultaka, Comput. Phys. Commun. 176 (2007), 426-455 [arXiv:hep-ph/0211331 [hep-ph]].
  * [21] S. Antusch, I. de Medeiros Varzielas, V. Maurer, C. Sluka and M. Spinrath, JHEP 09 (2014), 141 [arXiv:1405.6962 [hep-ph]].
  * [22] S. Antusch, C. Gross, V. Maurer and C. Sluka, Nucl. Phys. B 877 (2013), 772-791 [arXiv:1305.6612 [hep-ph]].
  * [23] S. Antusch, C. Hohl and V. Susič, JHEP 02 (2020), 086 [arXiv:1911.12807 [hep-ph]].
  * [24] S. Antusch, C. Hohl, C. K. Khosa and V. Susic, JHEP 12 (2018), 025 [arXiv:1808.09364 [hep-ph]].
  * [25] S. Antusch and C. Hohl, JHEP 10 (2017), 155 [arXiv:1706.04274 [hep-ph]].
  * [26] S. Antusch, C. Gross, V. Maurer and C. Sluka, Nucl. Phys. B 866 (2013), 255-269 [arXiv:1205.1051 [hep-ph]].
  * [27] S. Antusch and M. Spinrath, Phys. Rev. D 79 (2009), 095004 [arXiv:0902.4644 [hep-ph]].
  * [28] S. Antusch, S. F. King and M. Spinrath, Phys. Rev. D 89 (2014) no.5, 055027 [arXiv:1311.0877 [hep-ph]].
  * [29] S. P. Martin, Adv. Ser. Direct. High Energy Phys. 21 (2010), 1-153 [arXiv:hep-ph/9709356 [hep-ph]].
  * [30] S. Antusch, S. F. King, C. Luhn and M. Spinrath, Nucl. Phys. B 850 (2011), 477-504 [arXiv:1103.5930 [hep-ph]].
  * [31] S. Antusch, S. F. King and M. Spinrath, Phys. Rev. D 87 (2013) no.9, 096018 [arXiv:1301.6764 [hep-ph]].
  * [32] S. Antusch, S. F. King, M. Malinsky and M. Spinrath, Phys. Rev. D 81 (2010), 033008 [arXiv:0910.5127 [hep-ph]].
  * [33] J. R. Ellis, M. K. Gaillard and D. V. Nanopoulos, Phys. Lett. B 88 (1979), 320-324
  * [34] J. Ellis, J. L. Evans, N. Nagata, K. A. Olive and L. Velasco-Sevilla, Eur. Phys. J. C 80 (2020) no.4, 332 [arXiv:1912.04888 [hep-ph]].
  * [35] S. Antusch and M. Spinrath, Phys. Rev. D 78 (2008), 075020 [arXiv:0804.0717 [hep-ph]].
  * [36] S. Antusch and V. Maurer, JHEP 11 (2013), 115 [arXiv:1306.6879 [hep-ph]].
  * [37] M. Tanabashi et al. [Particle Data Group], Phys. Rev. D 98 (2018) no.3, 030001
  * [38] R. Brock, C. K. Jung, C. E. M. Wagner, K. S. Babu, R. Dermisek, B. Dutta, P. F. Perez, K. S. Ganezer, I. Gogoladze and Y. Kamyshkov, et al. “Proton Decay,” Workshop on Fundamental Physics at the Intensity Frontier, 111-130.
  * [39] A. Takenaka et al. [Super-Kamiokande], [arXiv:2010.16098 [hep-ex]].
  * [40] B. Bajc, J. Hisano, T. Kuwahara and Y. Omura, Nucl. Phys. B 910 (2016), 1-22 [arXiv:1603.03568 [hep-ph]].
  * [41] K. Abe et al. [Super-Kamiokande], Phys. Rev. D 95 (2017) no.1, 012004 [arXiv:1610.03597 [hep-ex]].
  * [42] K. Abe et al. [Super-Kamiokande], Phys. Rev. D 90 (2014) no.7, 072005 [arXiv:1408.1195 [hep-ex]].
  * [43] C. Regis et al. [Super-Kamiokande], Phys. Rev. D 86 (2012), 012006 [arXiv:1205.6538 [hep-ex]].
  * [44] H. Bahl, T. Hahn, S. Heinemeyer, W. Hollik, S. Paßehr, H. Rzehak and G. Weiglein, Comput. Phys. Commun. 249 (2020), 107099 [arXiv:1811.09073 [hep-ph]].
  * [45] H. Bahl, S. Heinemeyer, W. Hollik and G. Weiglein, Eur. Phys. J. C 78 (2018) no.1, 57 [arXiv:1706.00346 [hep-ph]].
  * [46] H. Bahl and W. Hollik, Eur. Phys. J. C 76 (2016) no.9, 499 [arXiv:1608.01880 [hep-ph]].
  * [47] T. Hahn, S. Heinemeyer, W. Hollik, H. Rzehak and G. Weiglein, Phys. Rev. Lett. 112 (2014) no.14, 141801 [arXiv:1312.4937 [hep-ph]].
  * [48] M. Frank, T. Hahn, S. Heinemeyer, W. Hollik, H. Rzehak and G. Weiglein, JHEP 02 (2007), 047 [arXiv:hep-ph/0611326 [hep-ph]].
  * [49] G. Degrassi, S. Heinemeyer, W. Hollik, P. Slavich and G. Weiglein, Eur. Phys. J. C 28 (2003), 133-143 [arXiv:hep-ph/0212020 [hep-ph]].
  * [50] S. Heinemeyer, W. Hollik and G. Weiglein, Eur. Phys. J. C 9 (1999), 343-366 [arXiv:hep-ph/9812472 [hep-ph]].
  * [51] S. Heinemeyer, W. Hollik and G. Weiglein, Comput. Phys. Commun. 124 (2000), 76-89 [arXiv:hep-ph/9812320 [hep-ph]].
  * [52] Gareth O. Roberts and Jeffrey S. Rosenthal, “Examples of Adaptive MCMC,” Journal of Computational and Graphical Statistics, 18:2 (2009), 349-367.
  * [53] J. E. Camargo-Molina, B. O’Leary, W. Porod and F. Staub, Eur. Phys. J. C 73 (2013) no.10, 2588 [arXiv:1307.1477 [hep-ph]].
  * [54] F. Staub, [arXiv:0806.0538 [hep-ph]].
  * [55] F. Staub, Comput. Phys. Commun. 185 (2014), 1773-1790 [arXiv:1309.7223 [hep-ph]].
  * [56] S. Dimopoulos and H. Georgi, Nucl. Phys. B 193 (1981), 150-162
  * [57] A. Masiero, D. V. Nanopoulos, K. Tamvakis and T. Yanagida, Phys. Lett. B 115 (1982), 380-384
  * [58] B. Grinstein, Nucl. Phys. B 206 (1982), 387
  * [59] R. Slansky, Phys. Rept. 79 (1981), 1-128
  * [60] R. Feger and T. W. Kephart, Comput. Phys. Commun. 192 (2015), 166-195 [arXiv:1206.6379 [math-ph]].
  * [61] R. Feger, T. W. Kephart and R. J. Saskowski, Comput. Phys. Commun. 257 (2020), 107490 [arXiv:1912.10969 [hep-th]].
  * [62] S. P. Martin and M. T. Vaughn, Phys. Rev. D 50 (1994), 2282 [erratum: Phys. Rev. D 78 (2008), 039903] [arXiv:hep-ph/9311340 [hep-ph]].
  * [63] J. Wess and B. Zumino, Phys. Lett. B 49 (1974), 52
  * [64] J. Iliopoulos and B. Zumino, Nucl. Phys. B 76 (1974), 310
  * [65] M. T. Grisaru, W. Siegel and M. Rocek, Nucl. Phys. B 159 (1979), 429
  * [66] B. Bajc, P. Fileviez Perez and G. Senjanovic, Phys. Rev. D 66 (2002), 075005 [arXiv:hep-ph/0204311 [hep-ph]].
  * [67] B. C. Allanach, C. Balazs, G. Belanger, M. Bernhardt, F. Boudjema, D. Choudhury, K. Desch, U. Ellwanger, P. Gambino and R. Godbole, et al. Comput. Phys. Commun. 180 (2009), 8-25 [arXiv:0801.0045 [hep-ph]].
  * [68] T. Nihei and J. Arafune, Prog. Theor. Phys. 93 (1995), 665-669 [arXiv:hep-ph/9412325 [hep-ph]].
  * [69] N. Cabibbo, E. C. Swallow and R. Winston, Ann. Rev. Nucl. Part. Sci. 53 (2003), 39-75 [arXiv:hep-ph/0307298 [hep-ph]].
  * [70] N. Tsutsui et al. [CP-PACS and JLQCD], Phys. Rev. D 70 (2004), 111501 [arXiv:hep-lat/0402026 [hep-lat]].