Probing baryogenesis with neutron-antineutron oscillations

In the near future, the Deep Underground Neutrino Experiment and the European Spallation Source aim to reach unprecedented sensitivity in the search for neutron-antineutron ($n\text{-}\bar{n}$) oscillations, whose observation would directly imply $|\Delta B| = 2$ violation and hence might hint towards a close link to the mechanism behind the observed baryon asymmetry of the Universe. In this work, we explore the consequences of such a discovery for baryogenesis first within a model-independent effective field theory approach. We then refine our analysis by including a source of CP violation and different hierarchies between the scales of new physics using a simplified model. We analyse the implication for baryogenesis in different scenarios and confront our results with complementary experimental constraints from dinucleon decay, LHC, and meson oscillations. We find that for a small mass hierarchy between the new degrees of freedom, an observable rate for $n\text{-}\bar{n}$ oscillation would imply that the washout processes are too strong to generate any sizeable baryon asymmetry, even if the CP violation is maximal. On the other hand, for a large hierarchy between the new degrees of freedom, our analysis shows that successful baryogenesis can occur over a large part of the parameter space, opening the window to be probed by current and future colliders and upcoming $n\text{-}\bar{n}$ oscillation searches.


Introduction
Baryon number (B) is an accidental global symmetry in the Standard Model (SM), which explains the empirically severe experimental limits from the non-observation of the proton decay. Baryon number conservation is violated in the SM only at finite temperature through non-perturbative instanton effects. However, B is expected to be violated in many well motivated ultraviolet (UV) extensions, e.g. grand unified theories (GUTs) naturally violate B, given the fact that quarks and (anti)leptons are often placed in the same representation(s) of the GUT gauge group. Finally, one of the big questions of particle physics is the observed baryon asymmetry in the Universe: the overabundance of baryons over antibaryons, quantified by the measured baryon-to-photon number density ratio [1] η obs B = (6.20 ± 0.15) × 10 −10 .
A theoretical explanation of the dynamical generation of such a baryon asymmetry requires the three Sakharov conditions [2] -(i) B − L violation (where L is the lepton number), (ii) C and CP violation and (iii) a departure from thermal equilibrium -to be fulfilled, where the first condition can be induced via B violation. Proton decay modes, e.g. p → e + π 0 , mediated via dimension-six operators, attributed to |∆B| = 1 and |∆(B − L)| = 0 can directly probe very high scales, O(10 16 GeV). The severity of the experimental limits on the non-observation of the proton decay might lead to the naïve expectation that dimension-nine |∆B| = 2 operators mediating neutron-antineutron (n-n) oscillations must be even more suppressed when compared with single-nucleon decay modes. However, this is true only when a single heavy new physics (NP) mass scale is involved. In fact, the presence of more than one new scale beyond the SM might suppress single-nucleon decay, while mediating n-n oscillations at a level comparable to the current experimental limits. Intriguingly, |∆B| = 2 observables such as n-n oscillations or dinucleon decays can be intimately connected to the baryon asymmetry of the Universe as these processes violate B − L by two units (|∆(B − L)| = 2). Therefore, |∆B| = 2 processes like n-n oscillations can be used to verify and probe baryogenesis mechanisms through B (and B − L) violation directly 1 .
During the recent years, on the one hand, lattice-QCD calculations have been improved tremendously in computing the QCD matrix elements that connect amplitudes of B-violating interactions to n-n oscillations [17], on the other hand, the current experimental sensitivities and future prospects for the observation of n-n oscillations have improved significantly. Currently, the most stringent constraint from the bound n-n oscillation is due to the Super-Kamiokande experiment [18], which provides a limit on the n-n oscillation lifetime τ SK n-n ≥ 4.7 × 10 8 s. The current best limit from the free n-n oscillation is due to the ILL experiment [19] τ ILL n-n ≥ 0.86 × 10 8 s. Experimental sensitivities from both, free and bound n-n oscillation times, are expected to be improved significantly in future experiments. The DUNE experiment [20] is expected to achieve a sensitivity of τ DUNE n-n ≥ 7 × 10 8 s using 40 Ar nuclei, while NNBAR [21] will exploit the Large Beam Port of the ESS facility to search for free n-n oscillations and is expected to achieve an impressive sensitivity of τ NNBAR n-n ≥ 3 × 10 9 s. Given the current stringent limits and expected future improvements in the experimental sensitivities for n-n oscillations, a detailed study of the phenomenological implications of such searches for baryogenesis mechanisms is timely and of high theoretical importance, since n-n oscillations are among very few observables which provide an opportunity to directly probe baryogenesis mechanisms and to distinguish underlying NP scenarios in synergy with direct searches at the high-energy frontier.
In this work, we explore the phenomenological possibility of probing baryogenesis using n-n oscillations and other complementary observables at the high-energy and high-intensity frontiers. We commence with an effective field theory (EFT) framework for n-n oscillations and study the impact of the current and future experimental limits of the n-n oscillation lifetime on the viability of realising a successful baryogenesis mechanism. Taking into account the latest lattice-QCD computations of the QCD matrix elements and relevant renormalisation group (RG) running effects, we 1 Note that in some early baryogenesis realisations in GUT theories [3][4][5][6] containing baryon number violation, B − L is conserved while B + L is violated. In such scenarios any asymmetry generated below the unification scale gets washed out by electroweak sphaleron interactions [7] and therefore a successful baryogenesis mechanism cannot be realised. Some alternatives include electroweak baryogenesis (which does not work within the SM alone, but can potentially work in some SM extensions, see e.g. Ref. [8]), leptogenesis [9] connected to the seesaw mechanism [10][11][12][13][14][15][16] of neutrino masses, as well as high-scale baryogenesis. In this work we will primarily be interested in the latter scenario.
first present a general framework for estimating the washout processes to derive model-independent limits on the viable scale for baryogenesis. In order to accommodate the possibility of different hierarchies of NP within the effective operator and an additional new source of CP violation, we explore then a simplified set-up to perform a comprehensive phenomenological study of the viability of baryogenesis above the electroweak scale in the context of an observable n-n oscillation lifetime. In particular, we focus on the B (and B − L) violating trilinear scalar coupling topology for n-n oscillation originally proposed in Ref. [22]. In order to make our analysis as general as possible, we consider a minimal simplified extension of the SM with diquark scalar fields coupling to SM quark fields and a B (and B − L) violating trilinear scalar coupling involving diquark scalar fields. Interestingly, a split scenario featuring some diquarks at TeV scale and some around GUT scale leads to the interesting possibility of realising a high-scale baryogenesis mechanism that can be readily embedded in many well motivated UV realisations. Moreover, it can be probed using the synergy between n-n oscillations and direct searches at the colliders, where the scalar diquarks are subject to extensive searches. The most stringent current constraints on the mass of the diquarks are already at the level of a few TeV for order unity couplings and are expected to be improved significantly in future searches. However, phenomenologically, diquark masses can generally also lie within the range from a few TeV to the GUT scale 2 . While there are a few instances of relevant studies for baryogenesis for some comparable scenarios in the literature [27][28][29], in this work we present for the first time a detailed and consistent prescription for the derivation of applicable Boltzmann equations and study the different relevant cases of phenomenological interest. We include in a comprehensive manner all experimental and theoretical constraints relevant for constraining the parameter space comprising n-n oscillations, meson oscillations, dinucleon decay, LHC constraints, and limits from a colour preserving vacuum. For instance, in contrast to n-n oscillations, the dinucleon decay is particularly relevant for TeV-scale masses of diquarks due to the stringent constraints on the couplings of diquarks to first generation quarks from meson oscillations.
Our findings suggest that the complementarity between n-n oscillations and LHC searches for diquarks can probe the baryogenesis mechanism extensively, ruling out the possibility of successful baryogenesis in some scenarios. The current best limit from the Super-Kamiokande experiment on the n-n oscillation lifetime together with the latest CMS limits exclude a large part of the viable parameter space for successful high-scale baryogenesis where one of the diquarks features a GUT scale mass while another one lies in the collider-accessible TeV scale range. However, we demonstrate that in case the future searches for n-n oscillations observe a signal, high-scale baryogenesis still remains a viable option to generate the correct observed baryon asymmetry of the Universe. On the other hand, for a scenario with all scalar diquark masses being in a similar mass range ( 10 5 TeV), the corresponding washout processes prove to be too strong to create a sizeable baryon asymmetry. Therefore, if the relevant scalar diquark fields are to be discovered by the LHC or future collider searches to lie within a few tens of TeV, then the possibility of a viable baryogenesis scenario featuring no additional particles at a high scale is completely ruled out. Baryogenesis scenarios with all diquarks having masses 10 5 TeV can still work successfully, however, remain unfortunately largely inaccessible by current and future experiments.
The paper is organised as follows: We introduce in Sec. 2, an effective field theory (EFT) framework for the operators mediating n-n oscillations, which we then use to set limits on the corresponding Wilson coefficients and consequently on NP mass scales and couplings in the subsequent sections. In Sec. 3, we first present a model-independent analysis of the washout of a pre-existing baryon asymmetry due to the effective operators mediating n-n oscillations without introducing any new sources of CP violation. To explore the possibility of different hierarchies of NP within the effective operators and the impact of a new source of CP violation, we present two possible topologies for realising n-n oscillations and study one of them using a simplified model set-up to perform a comprehensive phenomenological analysis. To conclude this section, we present the Boltzmann equation framework that we used to study the evolution of the baryon asymmetry including a detailed discussion of the CP violating decays and all relevant washout processes. In Sec. 4, we discuss all phenomenologically relevant constraints on such a set-up including the limits from LHC and future collider searches, neutral meson oscillations, dinucleon decay, colour preserving vacua and comment on some other observables. In Sec. 5, we combine all experimental constraints with the predictions for the final baryon asymmetry for two distinct scenarios (classified by the hierarchies of the NP scales involved) to present the parameter space for successful baryogenesis and discuss the implications for current and future experiments. Finally, in Sec. 6 we summarise and make concluding remarks.

Neutron-antineutron oscillations in effective field theory
Neutron-antineutron (n-n) oscillations violate baryon number by two units (|∆B| = 2), and therefore must be induced by some NP beyond the SM in case of an observation. Given any new physics model (e.g. the simplified model that we will consider in Sec 3), it is convenient to match the new physics operators mediating n-n oscillations to the effective operators involving only light fields at the scale where the heavy NP is integrated out. The effect of the heavy new physics then can be encoded in the Wilson coefficients of the effective operators, while the operators can be rundown to the scale of n-n oscillation and be identified with hadronic matrix elements, which are available from the lattice QCD computations. Therefore, we proceed to discuss below the relevant EFT formalisms which are of particular interest and provide an independent set of operators (and their relation to other commonly used operator bases in the literature) and hadronic matrix elements, which we then subsequently use for our simplified model.
At the QCD scale, the effective Lagrangian for n-n oscillations (after integrating out the heavy degrees of freedom) consists out of six SM quark fields with associated Wilson coefficients. These operators correspond to scattering processes violating baryon number at the temperature of interest for baryogenesis (up to the effects of RGE running of the Wilson coefficients). Therefore, the relevant effective operators (and the associated Wilson coefficients) are directly correlated with the washout processes which provide the possibility of probing the effectiveness of a given baryogenesis mechanism using the current and expected future experimental limits on the n-n oscillation lifetime. Since the n-n oscillation operators will be important for our study, we briefly summarise the formalism for the effective operator bases and the relevant RG running effects that we will later use for the analysis to obtain the corresponding n-n oscillation rates. In subsection 2.1, we first survey one of the commonly used SU (3) c × U (1) EM invariant EFT bases, and comment on possible connections of this basis with a SMEFT formalism. In subsection 2.2, we introduce the operator basis that we use in the rest of this work. This basis is also SU (3) c × U (1) EM invariant and additionally obeys a chiral SU (2) L ⊗ SU (2) R symmetry, which is commonly used in the literature for lattice QCD computations of the relevant hadronic matrix elements. We also provide the relations between the operator basis in subsection 2.1 and the one in subsection 2.2, and provide a list of independent hadronic matrix elements necessary to compute the n-n oscillation rate for a given NP model. Finally, in subsection 2.3, we provide the prescription to compute the n-n oscillation rate including the RG running effects for the hadronic matrix elements.
Note that an effective Lagrangian for n-n oscillations relevant at scales below the electroweak symmetry breaking must ensure that the relevant six-quark operators preserve SU (3) c × U (1) EM . A complete basis of such six-quark operators of the form uudddd, relevant for n-n oscillations, can be constructed as follows [30][31][32][33][34]: Hereby χ = {L, R} indicates the chirality with P L,R = 1 2 (1∓γ 5 ) being the chiral projection operators. The contraction of the spinor indices are implicitly assumed in the parentheses, C is the chargeconjugation operator and the quark colour tensors are defined as where { } denotes index symmetrisation and [ ] denotes index antisymmetrisation. Note that the operators involving the diquark invariants of the vector form (q T CP χ γ µ q) or tensor form (qCP χ σ µν q) are not independent and can be expressed as linear combinations of the operators given in Eq. (2) by performing Fierz transformations on a relevant subset of four fermions. Each of the operators in Eq. (2) leads to eight distinct operators when all possible combinations of chiralities are considered, leading to 24 operators in total. However, imposing the relations due to antisymmetrisation each of the operators in Eq. (2) leads to only six distinct operators making the total number of operators 18. Out of these 18 possible operators, four can further be eliminated due to the relation where χ, χ ∈ [L, R]. If the NP fields mediating the n-n oscillations are much heavier than the electroweak symmetry breaking scale, then the relevant effective six-quark operators in the effective Lagrangian (valid above the electroweak symmetry breaking scale), after integrating out the heavy NP degrees of freedom, must be SM gauge group invariant.
In passing we note that, since at the energy scale of n-n oscillations the whole unbroken SM gauge group need not to be respected a SU (3) c × U (1) EM invariant basis is the more appropriate choice of EFT. In case the requirement of invariance under SU (2) L × U (1) Y is imposed in addition to the SU (3) c × U (1) EM (e.g. in the case of a SMEFT [35][36][37][38][39] formulation), then a subset of only four independent operators survive [33,34], e.g.

SU
Most of the recent robust calculations for the hadronic matrix elements are performed using lattice-QCD simulations including nonzero quark masses and matched to massless chiral perturbation theory in which the chiral symmetry SU (2) L ⊗ SU (2) R is approximately preserved [40]. Since the relevant latest hadronic matrix elements are readily available in an SU (3) C × U (1) EM invariant chiral basis with a SU (2) L ⊗ SU (2) R symmetry [17], we find it convenient to work with it for the numerical evaluation of the n-n oscillation rates. In this framework the n-n oscillations can be described by the effective Lagrangian where C i are the Wilson coefficients corresponding to the set of effective operators O i defined as [17]: which are related to the remaining seven independent operators O P i by a parity transformation, accounting for total 14 independent operators. Here ψ corresponds to the isospin doublet ψ = (u, d) T , C corresponds to the charge conjugation operator, τ a denote the Pauli matrices for i = 1, 2, 3 and τ ± = 1 2 (τ 1 ± iτ 2 ). We have dropped the colour subscripts of the fields and the colour tensors T AAS(SSS) for brevity. In Eq. (8), the first equalities provides the relation between the new basis and the SU (3) c × U (1) EM invariant basis defined in Eq. (2).
As an useful remark, we note that many of the NP models (e.g. the simplified model considered in this work in Sec. 3), introduce two additional operatorsÕ 1 andÕ 3 , given by [40] However, these operators are not independent with respect to the complete basis of 14 operators included in Eq. (8). In fact, in dimension D = 4, the operatorsÕ 1 andÕ 3 are equal to O 1 and O 3 , respectively, by Fierz relations. However, such Fierz relations are broken by dimensional regularisation, therefore, in addition to the operators in Eq. (8) one must also include O 1 −Õ 1 and O 3 −Õ 3 as evanescent operators (vanishing for D = 4) for a complete treatment of the EFT at an arbitrary D. Alternatively, one can choose to includeÕ 1 andÕ 3 explicitly as a part of the physical basis of EFT operators.
To compute the n-n oscillation rate, we will be interested in the hadronic matrix elements associated with the operators defined in Eq. (8). To this end, we note that the isospin symmetry approximation further reduces the number of the relevant n-n matrix elements, making the matrix elements associated with three of the operators in Eq. (8) redundant [17], as follows. The hadronic matrix element for O 4 vanishes n|O 4 |n = 0 , in the approximate limit m u = m d (even after including the isospin breaking effects this matrix element is suppressed by powers of (m u − m d )/Λ QCD ). On the other hand, in the presence of isospin symmetry the hadronic matrix elements for O 6 and O 7 are related to that of O 5 by Therefore, for all practical purposes we work with in total four independent hadronic matrix elements to 3 Here we will not try to construct a complete SMEFT invariant basis but we will rather assume that the SMEFT basis is matched to the SU (3) C × U (1) EM invariant EFT introducing the relevant electroweak vacuum expectation values (VEVs). O 1 , O 2 , O 3 , O 5 . Given an explicit model, we first compute the relevant Wilson coefficients corresponding to a given operator in Eq. (8) at the scale where the heavy NP is integrated out and then identify the operator with one of the four independent hadronic matrix elements, which are subject to running effects between the scale where the lattice-QCD nucleon matrix elements are available in the MS scheme and the heavy NP scale where the effective Lagrangian is defined, as discussed in the following subsection.

n-n oscillation lifetime and RG running effects
The transition rate for the n-n oscillations, τ −1 n-n , is related to the Lagrangian in Eq. (7) via Here, M i (µ) is the transition matrix element for operator O i (µ), with M i (µ) ≡ n|O i |n . These can be determined using lattice-QCD techniques as described in [17], which also provides the relevant numerical values for M i (µ) in the MS scheme at µ = 2 GeV. Since we are interested in the case where the NP scale is much higher as compared to the n-n scale, i.e. where the relevant Wilson coefficients are defined at some heavy NP scale, the running of the operators between the different scales needs to be considered. Hence, to evolve the transition nuclear matrix elements from the relevant lattice-QCD scale to some heavy NP scale, it is necessary to perform the RG running of the EFT operators, for which we follow the prescription of [40], as detailed below. The running of the operators from the lattice scale µ 0 to a higher scale µ N P , to first order in the strong coupling constant α S , is described by the following equation [40]: up to O(α 2 s ). Here, q 1 and q 2 are two mass scales, with q 2 < q 1 , and N f is the number of quark flavours with masses above q 2 . Furthermore, the one-and two-loop MS anomalous dimensions γ   i , and the one-loop matching factor r (0) i , for each operator O i taken from [17] and [40].
is given, at 4-loop order, by [41] where L = ln with q α corresponding to the scale (with N f quarks on-shell) at which α S is known. The relevant β-functions in Eqs. (14) and (15) are given by Using Eq. (13), we obtain the matrix elements for each operator at a scale µ NP , by running them from the scale µ 0 = 2 GeV to the corresponding scale of NP µ NP . For a ready reference, we show the the relevant running for the matrix elements as a function of the NP scale in Fig. 1.  Neglecting the next-to-next-to-leading-order perturbative renormalisation effects the n-n transition rate can be expressed in terms of the nuclear matrix element at the scale µ = 2 GeV as [17] τ −1 n-n = 1.52 × 10 18 i=1,2,3,5 where M i (2 GeV) is given in Tab. 1 and U i (µ, 2 GeV) in Eq. 13. The relevant Wilson coefficients are obtained by computing the NP diagrams in a BSM scenario and integrating out the heavy NP degrees of freedom to match O i (µ) and consequently M i (µ) ≡ n|O i (µ)|n . Notice that Eq. (18) is expressed to make all the relative fractions involved dimensionless. It is then straightforward to notice that for C (P ) i ∼ Λ −5 , with Λ being the relevant NP scale, and taking M i (µ) ∼ O(10 −4 ) we obtain a n-n oscillation rate within the reach of current and future generation of n-n oscillation experiments (τ n-n O(few) × 10 8 s) for Λ O(100) TeV. We further notice from Fig. 1 that the the matrix elements run rather moderately as a function of the NP scale keeping the order of magnitude for the matrix elements the same across the relevant mass scales. Therefore, to demonstrate the model-independent implications of an observable rate of n-n oscillations for baryogenesis in the following section, we will take O 1 as a working example with the corresponding matrix element M 1 shown in Fig. 1. As a numerical example, considering the operator O 1 , we take the Wilson coefficient to be C 1 (µ) ∼ Λ −5 1 , expressed in terms of a NP scale Λ 1 , such that for the most stringent current experimental limit from the Super-Kamiokande experiment [18] τ SK n-n ≥ 4.7×10 8 s, we obtain Λ 1 ≥ 7.04 × 10 5 GeV. For the NNBAR experiment with an expected improved future sensitivity by an order of magnitude τ NNBAR n-n ≥ 3 × 10 9 s we find the limit Λ 1 ≥ 1.02 × 10 6 GeV.

Implications for baryogenesis
In the following, we want to study first the model-independent consequences of a possible observation of n-n oscillations on baryogenesis models. In the second part, we then study in particular two baryogenesis scenarios of our simplified model set-up that can lead to successful baryogenesis and confront the resulting parameter space with current and future constraints.
In order to generate a baryon asymmetry the three Sakharov conditions (i) B violation, (ii) CP violation and (iii) departure from thermal equilibrium have to be fulfilled. While the first condition is in principle fulfilled within the SM via the so-called sphaleron processes, the latter two are not sufficiently realised: The known CP violation in the SM is not enough and as the Higgs is too heavy in order to establish a first order phase transition, either new physics has to alter the potential or another mechanism is required such as out-of-equilibrium decays of heavy particles.
We will study in the following the consequences of baryon number violating operators relevant for n-n oscillations on baryogenesis. Hereby, it is important to keep in mind that above the electroweak phase transition, electroweak sphalerons, violating B + L, are highly active such that onlẏ is conserved. This means that when we study baryon number violating interactions, it is convenient to consider the change in where η X ≡ n X /n eq γ describes the number density of quantity X normalised to the photon density n γ , and whereη new ∆B tracks the yield of the new baryon-number-violating interactions. Below around T ∼ 10 12 GeV, when the electroweak sphalerons reach chemical equilibrium, one can easily use the known relation [42] η ∆(B−L) = 79 28 in order to solve for the final baryon asymmetry with the collision termη new ∆B arising from new baryon number violating interactions 3.1 Model-independent implications of n-n oscillation on baryogenesis models In Sec. 2, we have introduced the |∆B| = 2 effective operators that are relevant for n-n oscillations.
When the NP mediators are much heavier than the external quarks in the effective operators, then at any temperature below their mass scale, the operators relevant for n-n oscillations correspond to potential washout processes for any baryon asymmetry generated at a comparably higher scale. This mechanism can be either due to some CP-violating decay of any of the unstable mediators or due to a completely disconnected mechanism. Hence, an observed rate of n-n oscillations at experiments directly indicate washout effects in a model-independent way. In order to estimate their corresponding washout effect on baryogenesis scenarios, we will use the generalised Boltzmann equation formalism as described in Ref. [43]. Hereby, we consistently account for the running of relevant Wilson coefficients of the n-n operators between the scale of n-n oscillations and the scale at which the washout effects are of relevance.
The generic form of a Boltzmann equation for a particle species X is given by (see, e.g. Refs. [44][45][46]) where [Xa · · · ↔ ij · · · ] = n X n a · · · n eq X n eq a · · · γ eq (Xa · · · → ij · · · ) − n i n j · · · n eq i n eq j · · · γ eq (ij · · · → Xa · · · ) , where γ eq is the scattering density in thermal equilibrium. The Hubble rate H is given by with the effective number of degrees of freedom g * ∼ 107 in the SM and the photon density We can use this formalism to describe the evolution of baryon number over time. With the baryon number density per comoving photon defined as where the sum over (u, d) indicates the number of generations in thermal equilibrium. Note that we have used the left-handed fields along with the CP conjugates of the right handed fields, which are left handed antiparticles, (ψ c ) L = (ψ R ) c . In the 2-component Weyl spinor notation, the 4component Dirac spinors are then given e.g. as u = (u L ,ū c ) T . We summarise the notation for with g q = 3 being the number of degrees of freedom of the quarks. The chemical potential of a particle species a is related to the chemical potential of its antiparticle via µ a = −µā. When the SM Yukawa interactions and the sphalerons are in equilibrium, all relevant chemical potentials can be expressed in terms of a single chemical potential µ u L [42], with µ = e,µ,τ µ e L . Given the relations among the chemical potentials, we can express the baryon number density in terms of the chemical potential of a single species (which we choose to be u L ) In the following, we want to consider the n-n oscillation operator O 1 , which corresponds to O 1 in Eq. (8). We want to address the question, what the observation of n-n oscillations would imply for the washout of a baryon asymmetry that might have been created at a higher scale. We can write down the Boltzmann equation for the baryon-to-photon density by differentiating Eq. (27) In the presence of the operator O 1 the evolution of the abundance forū c ,d c , and their antiparticles can be obtained by writing the corresponding Boltzmann equation for a given species abundance; e.g.
where we have assumed three generations of fermions and a universal chemical potential among the three quark generations. The ellipsis denote other possible permutations of 3 ↔ 3 and 2 ↔ 4 processes. In order to arrive at the last line, we used the relation derived in Eq. (30). Similarly, one can write down the evolution forū c , u c and d c . On the other hand, in the absence of any B − L violating interactions involving u L and d L we can make the simplifying approximation d dz (η u L − ηū L ) 0 and d dz (η d L − ηd L ) 0. The evolution of baryon number density per comoving photon η ∆B is then given by Super-K DUNE NNBAR Figure 2: Left: Ratio between the width of the interaction induced by operator O 1 and the Hubble rate H, as a function of temperature T . When this fraction is greater than 1, as indicated by shaded areas, the interaction is assumed to provide a strong washout of baryon asymmetry. Right: Temperature at which the fraction Γ W /H falls below 1 (purple), and a more accurate limit on the out-of-equilibrium temperature coming from Eq. (36) (green). In blue and orange shaded areas, the experimental reach of the LHC and different n-n experiments are shown, respectively.
Following the prescription of [43,47,48] for the computation of thermal rate γ eq , the total washout effect from the operator O 1 can be expressed as where we have only included the 3 ↔ 3 scatterings since the 2 ↔ 4 scatterings are phase space suppressed. Therefore, the washout process corresponding to the n-n operator O 1 with an interaction rate Γ O 1 W can be regarded to be roughly in equilibrium if where c = π 2 c/(ζ(3) 3.3 √ g * ) ≈ 0.3 c and the Planck scale Λ Pl = 1.2 × 10 19 GeV. The Hubble parameter H and the equilibrium photon density n γ are defined in Eqs. (25) and (26), respectively. A more accurate limit on the out-of-equilibrium temperature of the washout process for a successful baryogenesis scenario can be obtained by integrating Eq. (34) to be given bŷ where d rec ≈ 1/27 is the dilution factor due to entropy conservation when the Universe cools down from the temperature of baryon asymmetry generation T = T * to the recombination temperature where g s is a function of temperature that maintains the relation s = (2π 2 /45)g s T 3 , where s is the entropy density. Furthermore, v is the vacuum expectation value of the SM Higgs, and T is the out-of-equilibrium temperature in Eq. (35) is given by In Fig. 2 left panel, we show the washout parameter Γ W /H as a function of temperature for the operator O 1 for different values of the EFT scale Λ corresponding to the integrated out heavy new physics. As discussed in Sec. 2, the most stringent limits from n-n oscillations constrain the NP scale to Λ ≥ 2.4 × 10 5 GeV. If in the future, n-n oscillations would be observed and not involve any CP-violating interactions, they would imply a strong washout down to the scale indicated by the corresponding shaded area in Fig. 2 (left). The implication, on the other hand, becomes more visible in Fig. 2 (right), where we show the out-of-equilibrium temperature for the washout processes corresponding to O 1 as a function of the the EFT scale Λ using Eq. (35) and Eq. (36). An observation of n-n oscillations around the scale of Λ ≈ 10 6 GeV, would imply a strong washout down toT ≈ 1.4 × 10 5 GeV. Under the assumption of a pre-existing, generated asymmetry at a high scale, this would imply such a strong washout, that an asymmetry must be generated below this scale and above the current exclusion limits from the LHC. Hence, such a discovery would hint towards new physics possibly observable at future colliders, e.g. a 100 TeV collider. Taking a completely agnostic approach towards the origin or flavour of the pre-existing asymmetry, as n-n oscillations strictly involve first generation quarks only, the conclusions regarding the washout effects derived in this section are strictly applicable only in the case of a pre-existing baron asymmetry in the first generation quarks. In order to ensure the complete washout of a pre-existing asymmetry in all flavours, a complementary measurement of new physics arising from the operator in question involving second or third generation quarks (e.g. LHC searches, meson oscillations etc.) is needed besides an observation of n-n oscillation, in the absence of flavour transitions. In order to study such a situation and its interplay with other experimental constraints, we will explore a simplified model (low-scale scenario) later in more detail. However, this is a conservative assumption as one would expect washout also in other flavours via spectator effects, as for instance discussed for leptogenesis in [49].

A comprehensive Boltzmann equation formalism for baryogenesis in models featuring n-n oscillations
While in the previous analysis, we have assumed that the new operator does not include any source of CP violation and only contributes to the washout of the baryon asymmetry, we want to refine our analysis by analysing a simplified set-up that allows for (a) including a source of CP violation in the NP operator and (b) different mass hierarchies within this operator. The latter is in particular important, as the the EFT approach presented in the previous subsection provides only an estimate of the washout within the validity of the EFT, i.e. when the masses of all new degrees of freedom are below its cutoff scale. For these purposes, we are interested in the most generic topologies that such an operator could lead to. At tree level, the realisations for the short-range n-n operators given in Eq. (8) can be classified in two possible topologies, which are shown in Fig. 3. In general, for topology I the internal mediators between vertices x 1 − x 2 and x 3 − x 4 can be either vector or scalar fields, while the particle between x 2 − x 3 must be a fermionic field. On the other hand, for topology II, all the internal mediators can be either scalars or vectors, with all possible combinations. Topology I has been explored in the context of many specific model realisations in Refs. [50][51][52][53][54][55][56][57][58][59][60][61][62][63][64][65][66][67][68] and recently, has been extensively discussed in the context of baryogenesis using a minimal simplified model setup in [69]. The main focus of the remaining of this work will be on topology II, which has been proposed originally in Ref. [22] and has been realised in many UV complete and TeV scale models [23-25, 27-29, 70-73] 4 . While there are a few instances of studies for this scenario in the context of baryogenesis [27][28][29], a phenomenologically comprehensive and in depth exploration of the viable parameter space for baryogenesis still remained desirable. To this end, in this work, we explore both high-and low-scale (pre-electroweak) baryogenesis in the context of an observable n-n oscillation lifetime while taking into account constraints from various complementary observables at the highenergy and high-intensity frontiers. Having already discussed the general EFT approach, in the following sections we consider a very general minimal set-up extending the SM with diquark scalar fields coupling to SM quark fields and a B (and B − L) violating trilinear scalar coupling involving only the diquark scalar fields. This simplified set-up not only allows us to develop an in-depth prescription for the Boltzmann equation formalism but also makes the analysis very general and directly applicable to TeV scale and UV complete model realisations of topology II.

The trilinear topology of n-n oscillations: a simplified model with scalar diquarks
We consider the following Lagrangian including scalar diquarks given by where ξ is a neutral complex scalar field, whose real part acquires a VEV and can be written as Hereby, χ corresponds to the relevant goldstone mode associated with the symmetry broken by the VEV and is absorbed by the associated gauge boson. Before the B − L symmetry is broken by ξ acquiring a VEV, the new fields X dd , X uu , and X ud can be assigned consistently a baryon number B = −2/3 (and lepton number L = 0), while ξ carries B − L = 2. However, once B − L is broken, the trilinear term will violate baryon number in units of |∆B| = 2, as can be seen in the n-n oscillation diagram in Fig. 4. The transformation properties of the scalar diquark fields under the SM gauge group and their associated baryon number charges are summarised in Tab. 3. Note that depending on the colour of the scalar diquarks, the associated Table 3: Transformation properties and charges under the SM gauge group of the relevant scalar diquark fields.
Yukawa couplings are either symmetric or antisymmetric with respect to the flavour indices. If X dd (X uu ) transforms as a colour triplet under SU (3) c then f dd ij (f uu ij ) must be flavour antisymmetric. On the other hand, if X dd (X uu ) transforms as a colour sextet then f dd ij (f uu ij ) must be flavour symmetric.
Another point worth noting is that given the transformations of the scalar diquark field X ud under the SM gauge group it can also couple to a left-handed quark doublet: a colour sextet (triplet) X ud can couple to flavour (anti-)symmetric pair of QQ. However, in the presence of X ud couplings to both left-and right-handed quark pairs simultaneously, large chiral enhancements can be induced for flavour changing neutral current (FCNC) ∆F = 2 operators (e.g. neutral meson mixing) and ∆F = 1 operators (e.g b → sγ), which in spite of being loop suppressed (with X ud in the loop) can lead to overwhelming rates disfavoured by the current stringent constraints from experiments for an X ud mass within the collider reach. Therefore, we assume that the couplings of X ud to a QQ pair is negligible even if it gets generated due to radiative corrections.
Furthermore, if X ud and X uu transform as colour triplets, they can potentially also have leptoquark couplings. Leptoquark couplings when present simultaneously with diquark couplings lead to rapid proton decay for low X ud (X uu ) masses in the absence of any symmetry forbidding one of Figure 4: Diagram for n-n oscillation with the diquark fields X dd and X ud descibed in the text.
the couplings [77][78][79]. In this work, we mainly focus on the case of colour sextet diquarks, implying flavour symmetric couplings f dd ij (f uu ij ) and absence of any leptoquark couplings (thereby avoiding any possibility of rapid proton decay). Interestingly, one can directly associate the Lagrangian in Eq. (38) assuming colour-sextet scalar diquarks with a UV complete realisation as discussed in detail in App. A.
As we will later see in more detail, with respect to realising baryogenesis in this model, two regimes are of particular interest: Such a choice is naturally motivated by UV completions such as SO(10) to obtain gauge coupling unification, as discussed in App. A. Low-scale scenario: We assume m Xuu m X dd > m X ud is maintained but with both m X ud and m X dd not too far from ∼ O(TeV), while again m Xuu ∼ m GUT . Even though such a scenario might require a more complex UV completion in order to address the colour vacuum stability as will be discussed in Sec. 4, this phenomenological scenario is particularly interesting as it allows for the possibility of two diquark states accessible at the LHC and future colliders. Also a number of ongoing low-energy experiments searching for baryon-number-violating processes e.g. dinucleon decay, are particularly relevant for such a scenario as will be discussed in Sec. 4.
For both scenarios above, the Lagrangian relevant for n-n oscillations can be written in an effective form (at a scale below m X dd ) as where we have integrated out the heavy degrees of freedom. The relevant tree-level n-n oscillation operator generated can be identified with O 2 RRR , defined in Eq. (2). The effective Lagrangian for n-n oscillation can be expressed, c.f. Eqs (8) and (9), as where the last approximation follows from Eq. (10) up to the uncertainties. For the numerical evaluation of the transition matrix element for operatorÕ 1 (µ), defined asM 1 (µ) ≡ n|Õ 1 |n we i , forÕ 1 given in Table 1. In a realistic model, such as the one discussed above, in addition to the n-n oscillation operator O 1 , a number of additional baryon number violating two-to-two scattering processes with NP particles (e.g. scalar diquark states) as the external states can give rise to washout processes. These processes become in particular important, when the new degrees of freedom feature a large hiearchy not captured within the validity of the previous EFT analysis. Therefore, taking the simplified model discussed above as a working example, we will explore the different relevant processes for CP violation and baryon asymmetry washout in full detail, accounting for different mass hierarchies of new physics.

Derivation of the Boltzmann equation framework
We assume that the baryon asymmetry is generated by the |∆B| = 2 decay X dd → X * ud X * ud (once ξ acquires a VEV) with the CP violation generated by the interference with loop diagrams involving the baryon number conserving decay mode X dd → d c d c , as shown in Fig. 5. It is straightforward to check that there is no vertex correction due to charge conservation. Only the self-energy diagram is present and it necessitates the introduction of an additional heavier X dd state, which we denote by X dd . To make our analysis as general as possible we consider the minimal relevant interactions of X dd as The CP parameter is defined as [27] ≡ where Γ(X dd → X * ud X * ud ) and Γ(X * dd → X ud X ud ) are the decay widths for the decays X dd → X * ud X * ud and X * dd → X ud X ud , respectively. The total decay width Γ tot ( and r denotes the branching ratio for the decay mode X dd → X * ud X * ud . Note that an alternative baryogenesis mechanism that can occur in this construction is the postsphaleron baryogenesis via the decay of S contained in ξ [23][24][25] with S being a real scalar field that can decay into six quarks (S →ū cdcūcdcdcdc ) and antiquarks (S → u c d c u c d c d c d c ) leading to |∆B| = 2. In such a scenario the relevant loop diagrams inducing the CP violation through interference employs W ± loops 5 , which necessitates flavour violation linking the CP violation directly to CKM mixing. In this scenario, it is therefore important to include loop contributions to n-n oscillations in addition to tree level ones. A detailed discussion of the post-sphaleron baryogenesis mechanism is beyond the scope of the current work and will be presented elsewhere [26].
Particle dynamics in the early Universe can be described using Boltzmann equations [44][45][46], which is well studied in particular for leptogenesis scenarios with right-handed Majorana neutrinos. However, there are some crucial differences in this scenario, which must be taken into account to obtain a consistent description that we will discuss in the following.
For our simplified construction, cf. Eq. (24), the relevant processes for baryogenesis can be classified as where the Feynman diagrams for the processes are shown in Figs. 6 and 7. Now to compare this scenario with respect to the standard leptogenesis scenario let us note some key observations. Since X dd can decay into either a pair of X ud or d c (after the B − L breaking), one should include the change of the number density for both of these species in the Boltzmann equation for the baryon asymmetry. However, baryon number is only violated in the interaction X dd ↔ X * ud X * ud and X * dd ↔ X ud X ud . This can seen from the assignments in Tab. 3. Note, that in our analysis, we do not consider the possibility of a three-body decay involving this interaction term under the assumption that m ξ m X dd , leading to an early decoupling of such a decay mode. Due to this baryon number assignment, the decay modes X dd ↔ d c d c and X * dd ↔d cdc do not violate B, cf. Eq. (38). Therefore, in view of the thermal decoupling of X dd (X * dd ), the decay of a pair of X dd and X * dd can be considered as the analog of the heavy out-of-equilibrium decay of the Majorana fermion in the standard leptogenesis scenario. However, in contrast to the standard scenario, it is important to note that in presence of a CP asymmetry between X dd ↔ X * ud X * ud and X * dd ↔ X ud X ud indirectly a CP asymmetry between X dd ↔ d c d c and X * dd ↔d cdc is generated. This is caused under the assumption of CPT invariance such that the total decay widths for X dd and X * dd must be the same, as illustrated in Fig. 8. We note that being coloured particles, X dd and X * dd can also be produced via two gluons gg → X dd X * dd , however such a mode is highly phase space suppressed at a temperature below m X dd (while the inverse process is doubly Boltzmann suppressed by m X dd /T ), which is the interesting temperature for the out-of-equilibrium decay of X ( * ) dd and its implications for baryogenesis. The process gg → X dd X * dd generally dominates for z < 1, as is also demonstrated in [29]. However, both the decay and the scattering washout channels dominate over any effect of gg → X dd X * dd for z > 1, making this gauge scattering mode inconsequential for the final baryon asymmetry and mainly relevant for bringing the X dd density to equilibrium at the onset of baryogenesis. The reason for this is that the rate of gg → X dd X * dd falls off much earlier than the washout via scatterings, due to the different dependence on the number density of X dd for T < m X dd . The gauge scattering features two X dd external legs, while the washout processes we consider have at most a single X dd external leg. Therefore, we do not include gauge scatterings in our analysis 6 .
In the existing literature, both X dd ↔ X * ud X * ud and X dd ↔ d c d c modes have been taken into account to some extent, either while calculating the CP violation (see e.g. [27] and the references therein), or by including the Boltzmann equations for both X ud and d c [29], however, a comprehensive formalism for the relevant Boltzmann equation for the baryon asymmetry is lacking. Here, we present a prescription for obtaining a consistent equation for the evolution of baryon asymmetry.
Based on the Boltzmann equations introduced in Sec. 3.1, in particular Eqs. (23) and (24) and the definition of the baryon asymmetry in Eq. (31), we can evaluate the final baryon asymmetry at the electroweak scale. After X ( * ) dd goes out of equilibrium, both X dd ↔ X * ud X * ud and X dd ↔ d c d c modes (together with their CP conjugate modes) can generate a B−L asymmetry and their interplay dictates the final baryon asymmetry. First, we discuss the necessary equations that describe the evolution of the out-of-equilibrium decay. The number density of X dd (X * dd ) per photon density at a given temperature can be obtained by solving the relevant Boltzmann equation given by where we note that the additional factor of 2 on the left-hand side accounts for the averaging factor since the right-hand side includes both X dd and X * dd mediated processes. The different decay and scattering rates are defined in Eqs. (24) and (43), such that Eq. (44) can further be expressed in terms of the rates in the following form The details of the prescription for obtaining the relevant density of scatterings γ ≡ γ eq is outlined in App. C. We would like to note that even if we initially assume η X dd = η X * dd , an asymmetry can be induced between η X * dd and η X dd through the back reaction of η X * ud η X * ud and η X ud η X ud (d c d c andd cdc ) because of the asymmetry generated between η X ud and η X * ud (d c andd c ). Since during baryogenesis we are interested in the out-of-equilibrium decay of X ( * ) dd , such a secondary asymmetry is relevant for the Boltzmann equation for η X ( * ) dd . However, such a contribution can be easily checked to be proportional to the generated baryon asymmetry η ∆B which is very small compared to η X ( * ) dd for the temperature where the baryon asymmetry is dynamically generated (i.e. before the baryon asymmetry freezes out) c.f. Figs. 15 and 18. Therefore, we find it a good approximation to take η X dd = η X * dd through out the whole regime of baryogenesis. Now, to derive the Boltzmann equation for the baryon asymmetry one must consistently define the net baryon number density per photon density in terms of the number densities of the relevant species carrying SM baryon number that are in thermal equilibrium at the time of baryogenesis, see Eq. (27). Again, we neglect d dz (η u L − ηū L ) and d dz (η d L − ηd L ) as u L and d L do not participate in any baryon-number-violating interactions generating or washing out the baryon asymmetry. Therefore, they can be decoupled from the Boltzmann equation 7 . Note that we include the number density of X ud in our definition of the final baryon asymmetry, which is valid in the regime when X ud is in thermal equilibrium. We have checked that this assumption holds for all scenarios presented in this work. Once X ud goes out of equilibrium one would need to write a new Boltzmann equation with the baryon number expressed in terms of only the SM quarks, taking into account the decay modes of X ud . Since we are interested in the case where m X dd > m X ud m u,d , the baryon asymmetry generation freezes out by z ≡ m X dd /T ∼ O(10). Therefore, for our analysis it will suffice to consider Eq. (46) with the assumption that the asymmetry generated in X ud − X * ud gets redistributed into SM quarks once X ud goes out of equilibrium (at a temperature below the baryon asymmetry generation freeze out), see Fig. 9 for an illustration. Hence, we arrive at the following definition, where the sum over (u, d) is over the number N of generations in thermal equilibrium. Figure 9: A timeline (not to scale) that describes the chronology of events in the baryogenesis mechanism discussed in this section. Before the time of B − L symmetry breaking (T v ), the particle species X dd , X ud , u c and d c , as well as their corresponding antiparticles, are in thermal and chemical equilibrium. When X ( * ) dd has decayed away (T m X dd ), an asymmetry has emerged in the relative number densities of X ud , u c , d c with their antiparticles. In the subsequent decay of X ( * ) ud (at T m X ud ), this asymmetry is translated into a baryon asymmetry of the Universe stored in SM quarks.
Then the Boltzmann equation for the baryon asymmetry per photon density can be obtained by differentiating Eq. (46): Now each term on the right-hand side of Eq. (47) corresponds to the standard Boltzmann equation describing the evolution of the number density of a particle species, (cf. Eqs. (23) and (24)) and is given by Now let us discuss the relevant terms in some detail before employing the chemical potential relations relating all the species in chemical equilibrium. The contributions to the baryon asymmetry

Process
Branching fraction Table 4: Two-body decay modes of X dd and the corresponding branching fractions.
originating from decays can be parameterised in terms of the relevant decay rates as where r corresponds to the average branching fraction for the X dd → X * ud X * ud and X * dd → X ud X ud decay modes and corresponds to average B asymmetry generated in the decay of a X dd (X * dd ), as defined in Eq. (42). We summarise the branching ratios in Tab. 4. We further assume CPT invariance implying the same total decay width for X dd and X * dd , and define γ X dd D is the total decay rate of the X dd and X * dd pair. Using the above parametrisation we obtain Similar to the standard leptogenesis scenario, it is important to note that the s-channel scattering mediated by X dd must be calculated subtracting the CP-violating contribution due to the on-shell contributions in order to avoid double-counting, as this effect is already taken into account by successive (inverse) decays, X * ud X * ud ↔ X dd ↔ d c d c (real intermediate state subtraction). Therefore, we can write with γ on-shell where the rightmost equalities are valid up to O( ). Therefore, after simplifying Eq. (57) we obtain Correspondingly, we obtain for the t-channel scattering mediated by X dd for the scattering processes with X dd in the initial or final state (for s-and t-channel, respectively), and similarly for the quark mediated scatterings, In order to solve the Boltzmann equations, we have to translate the different baryon densities on the right-hand side of Eq. (47) into η ∆B . As a first step, we express the ratio of the number density over the equilibrium density for different species in terms of the chemical potentials, recalling the approximation η a /η eq a ≈ e µa/T ≈ 1 + µa T for a species a with chemical potential µ a . Then we can express all relevant chemical potentials in terms of the chemical potential of a single species (following a similar prescription as in Sec. 3.1) under the assumption that in addition to the SM Yukawa interactions and sphalerons, the X * ud ↔ū cdc interactions are also in thermal equilibrium. We express all appearing chemical potentials for any particle species a in terms of which we summarise in App. B. Finally, we arrive at the baryon asymmetry expressed in terms of one chemical potential µ X * where C B is given by (cf. App. B) with C X * ud = 6 being the colour multiplicity of X * ud . Hence, the combination of decay terms relevant for the evolution of baryon asymmetry in Eq. (54) is given by (cf. Eq. (56)) where xdc is given in Eq. (133). Similarly, the various s-and t-channel scattering contributions in Eqs. (59), (60), (61), (62), and (63) can be expressed as whereγ and xūc and xdc are given in Eq. (133).
We then obtain the final form of the equation governing the evolution of the baryon asymmetry per photon density using Eqs. (68), (69), (70), (71) and (72) with Eq. (54) given by is the on-shell contribution subtracted scattering rate with the on-shell part given by γ on-shell as can be verified from Eqs. (58). Note that after solving Eq. (77), one has to include a factor d γ ≈ 1/27 for the dilution of the photon density in order to obtain the final baryon asymmetry at the recombination epoch T 0 . We included this dilution factor for the later presented numerical analysis to obtain the final η ∆B (T 0 ). For solving the Boltzmann equations for the high-scale and low-scale scenario defined in section 3.2.1 we make the following assumptions which we summarise below.
High-scale scenario: • Given that the SM top Yukawa coupling is in equilibrium for T 10 13 GeV, and bottom as well as tau Yukawa couplings come in equilibrium just below 10 13 GeV, we assume that all these Yukawas interactions as well as gauge interactions are in chemical equilibrium during the whole range of temperature relevant for high scale baryogenesis. This implies that the chemical potential relations due to transfer of the baryon asymmetry from SU (2) L singlet quarks to SU (2) L doublet quarks in equilibrium (c.f. Appendix B) by imposing the relevant chemical potential constraints due to Yukawa interactions. For a first exploration, we restrict ourselves to the single flavour approximation by assuming that the (inverse) decay of X ( * ) dd and X ( * ) ud dominantly occurs along the third generation quarks in flavour space. The respective chemical potential relations then correspond to the case N = 1 in Appendix B. In the future, it would be interesting to perform a detailed analysis including potential flavour correlations.
• Under the assumption that a single flavour of quarks b and t, as well as the gauge interactions, are in equilibrium, the chemical potential relation 2µ Q 3 − µūc 3 − µdc 3 = 0 is enforced, therefore making the chemical potentials of these single flavour quarks no longer affected by the QCD sphaleron constraint 3 i 2µ Q i − µūc i − µdc i = 0. However, the partial equilibriation of Yukawa interactions can lead to a change of the final yield by little (e.g. this effect is ∼ 10% [87] for leptogenesis), which we do not take into account as it has little bearing in our final parameter space conclusions. For simplicity, we assume that the electroweak sphalerons are also in equilibrium together with the QCD sphalerons from the beginning of baryogenesis.
Low-scale scenario: • We assume that the Yukawa couplings for all three generation of quarks as well as the gauge interactions are in equilibrium during the baryon asymmetry generation. The species X * ud , and the process X * ud ↔ū c +d c are also assumed to be in equilibrium for all three generations. We further assume that both QCD sphalerons and electroweak sphalerons are also in equilibrium. In this case the chemical potential relation 2µ Q i − µū i c − µd i c = 0 for i = 1, 2, 3 is enforced, which takes into account the effect of transfer of asymmetry from singlet to doublet quarks. The final chemical potential relations correspond to the case N = 3 in the Appendix B.
• To simplify the analysis we consider the case of flavour universal and diagonal couplings of X dd and X ud to three quark generations, implying that X For more details, we refer to e.g. [88,89] discussing potential flavour effects and their implications, which is beyond the scope of this work.
• To present the numerical results for the low-scale scenario in section 5, we allow for the possibility that X dd and X ud Yukawa couplings f ud and f dd are hierarchical among themselves (e.g. we consider the two benchmark cases f ud = f dd and f ud = 10f dd ), we assume that f ud and f dd are universal across all three generations.
Before discussing in detail the parameter space that leads to the observed baryon asymmetry, we will discuss the relevant phenomenological constraints.

Phenomenological constraints
Diquarks are subject to different phenomenological constraints, both from experimental observables as well as theoretical conditions. The relevant constrains are particularly strong for diquark states that are not much heavier than the electroweak symmetry breaking scale. For instance, the LHC probes already O(TeV) mass scales and provides high precision measurements of the Flavour Changing Neutral Current (FCNC) processes in mesonic observables. Moreover, dinucleon decays, which provide a complimentary probe to n-n oscillations, are also of particular interest with the expected improvements on future experimental sensitivities. Furthermore, considerations of a colour preserving vacuum provides useful constraints on the B − L breaking scale and mass hierarchy between diquark states, when more than one of the diquarks are light. In this section, we provide an overview of all relevant constraints for our study of the baryogenesis parameter space and will comment on additional constraints, which can be relevant in scenarios beyond the studied ones.

Direct LHC searches
Scalar diquarks have been studied extensively in the context of collider searches in the literature, see e.g. [28,[90][91][92][93][94][95][96][97][98][99][100][101][102][103][104]. At the LHC or future colliders they can be produced through the annihilation of a pair of quarks via an s-channel resonance decaying into two quarks producing a dijet final state. In principle, it is also possible to produce a pair of scalar diquarks e.g. through gluon-gluon fusion. For diquark couplings of O(0.1) and higher, the resonant production cross section dominates over the pair production [98]. For smaller couplings, the pair production cross section can potentially dominate over the resonant production (the resonant production contains one power of the diquark coupling, while the gluon fusion pair production only includes gauge couplings) [105]. As limits from the current collider searches for the pair production mode is only available till roughly TeV scale [106], while the resonant production search limits reach O(10) TeV [107], we consider only the resonant production, providing a very good estimate of the current LHC reach.
Since we assume in our simplified framework that X uu is very heavy, it will be beyond the collider reach; however, X ud being around few TeV of mass is actively probed by collider searches. The accessibility of X dd at colliders depends on the respective baryogenesis scenario (with a mass around the GUT scale in the high-scale scenario or around a few TeV for the low-scale scenario).
Depending on the generation of the quarks, we can distinguish between a resonant dijet or a top+jet signature. The partonic differential cross section for the latter (u i d j → X ud → td k ) is given by where we have neglected all quark masses except the top quark mass. Hereby, θ * denotes the scattering angle andŝ the center-of-mass energy of the partons. As will be discussed later, spin and colour multiplicity factors depending on the spin and colour representation of the initial and final states need to be added correspondingly. The total decay width Γ X ud of the scalar diquark X ud is given by the sum of its partial decay widths, where C X ud is the colour multiplicity of X ud . Similarly, one can obtain the relevant partonic cross section for the resonant dijet signature with the relevant partial decay width, Given the partonic scattering cross section, the experimentally measured hadronic cross section (e.g. at the LHC) can be obtained by employing the relevant parton distribution functions (PDFs) f (x) and summing over all partons [108] where Q is the characteristic scale of the interaction, e.g. the invariant mass in a two-to-two partonic scattering. Furthermore, µ F corresponds to the factorisation scale, which factorises the nonperturbative contributions from the short-distance hard scattering and µ R is the scheme-dependent renormalisation scale. Note that the scales µ F,R are parameters determined for a fixed-order calculation, while in general the total cross section should be independent of these scales at any given order of the perturbative expansion. The Parton Distribution Functions (PDFs) can be taken into account by introducing the parton luminosity factor defined as [108] dL ij dτ = where In the above relations the initial partons are carrying fractions x 1 and x 2 of the hadron momentum and the invariant mass of the two-parton system is defined asŝ ≡ x 1 x 2 s, with √ s being the energy of the colliding hadrons in center-of-mass frame. Furthermore, constraints are also imposed on the rapidities y of the final state partons observed at the collider experiments as jets. It is particularly convenient to express the parton luminosity in terms of the variables τ andȳ as where we use the identity dx 1 dx 2 = ∂(τ,ȳ) ∂x 1 ,x 2 dτ dȳ = dτ dȳ and the rapidity can be expressed as a function of the momentum fractions asȳ = 1/2 ln(x 1 /x 2 ) in the center-of-mass frame. Consequently, the total cross section can be expressed in terms of the parton luminosity factor and the partonic cross section as Even though an exact calculation of the cross section or the decay widths should include all possible Feynman diagrams, for all practical purposes, the experimental searches are principally focused on narrow resonances, expecting resonant peaks on a smoothly falling dijet mass spectrum corresponding to a s-channel decay mode of the resonance. Such a cross section for a resonance decaying via s-channel can be approximated by a Breit-Wigner form where Γ X and m X correspond to the total width and the mass of the resonance, and m √ŝ = √ τ s.
The partial widths Γ(a + b → X) and Γ(X → c + d) correspond to the creation of the resonance from specific initial states and the decay of the resonance to the specific final states, respectively. The different multiplicity factors are taken into account through N , defined as with N S X and N S a,b representing the spin multiplicities of the resonance and the initial state particles, respectively (e.g. for scalar diquarks N S X = 1 and N S a,b = 2 for initial state quarks). C X and C a,b are the relevant colour multiplicities (e.g. for colour sextet diquarks C X = 6 and for initial state quarks C a,b = 3). It is relevant to note here that the cross section above is obtained after integrating over cos θ * , which in practice is constrained by the kinematics of the experimental searches. Therefore it is convenient to express the Breit-Wigner partonic cross section aŝ where BR par corresponds to the branching fraction of the partonic subprocess, and A is the experimental acceptance factor after the cos θ * cut. The Breit-Wigner partonic cross section can further be simplified in the narrow-width approximation (Γ X m X ) leading to the hadronic cross section which, for given possibilities of initial (a, b) and final (c, d) state partons, can be expressed as where the factor (1 + δ ab ) accounts for the possibility of two identical incoming partons (which gets compensated by a factor 1/2 in the partial width of the final state phase space Γ(X → ab) for identical partons). We follow the prescription of Ref. [107] for the estimation of the acceptance parameter: defined as A = A ∆ A η , with A ∆ being the acceptance by requiring |∆η| < 1.3 for the dijet system and A η is the acceptance factor due to also requiring |η| < 2.5 for each jet, separately.
Assuming the decay of the scalar diquark resonances to be isotropic, one has A ∆ = 0.57 for all masses and A η 1 [107]. Note that in case of a single top or anti-top quark in the final state, the top quark can decay before hadronising. This provides the possibility of tagging it to reconstruct the invariant mass of the diquark resonance. To provide some crude benchmark estimates for the t + b final state we use the top tagging efficiency ∼ 0.3 and bottom tagging efficiency ∼ 0.8 for the decay channels with a t or b quark in the final state [109]. Using CTEQ6L1 [110] PDFs and computing the PDF luminosity functions using the ManeParse package [111], we show in Fig. 4.1 (top left) the cross sections of the scattering u + d → X ud → dijet and u + d → X ud → t + b as a function of the diquark mass m X ud for three different benchmark values of the coupling f ud . Given that the current best limits on the t + b final state searches [109] are at best of the order of a dijet final state search, the limits considering the t + b final state are less stringent as compared to that of dijet searches. While we still show the t + b (b + b) final state in our plots, we consider only the dijet constraints in the following. Fig. 4.1 (top right) shows the corresponding cross sections of the scattering processes d+d → X dd → dijet and d+d → X dd → b+b as functions of the diquark mass m X dd for three different benchmark values of the coupling f ud . Again, we show the exclusive two bottom quark final states separately, due to the possibility of tagging the b quarks in light of recent experimental improvements in b-tagging. For reference, we indicate the current experimental limits on dijet searches from the CMS collaboration [107] as well as the SM prediction.
The current search limits already exclude parts of the parameter space for mass ranges as high as 8 TeV, which is expected to be improved significantly with more data and future colliders e.g. HE-LHC with 27 TeV center-of-mass energy and FCC-hh with 100 TeV center-of-mass energy [112,113]. In Fig. 4.1 bottom row we show a comparison of expected cross sections at LHC, HE-LHC and FCC-hh for coupling f ud (f dd ) of order unity. Given the expected several ab −1 luminosity from future LHC upgrades, the collider searches will play a key role in probing the parameter space of baryogenesis complimenting the n-n oscillation searches.

Meson oscillations
FCNC processes such as neutral meson oscillations and flavour changing non-leptonic meson decays provide stringent constraints on the masses and couplings of the scalar diquark states [105]. Since we assume X uu to be very heavy 8 , we mainly focus here on constraints on X dd and X ud . The contributions to neutral meson mixing, M − M with M = B s , B d , K, can be described as where O j are the effective four-quark operators and C SM (NP) j denotes the associated Wilson coefficients for the SM (NP).
For the SM, the relevant operator for the |∆F | = 2 process B q − B q is given by with the associated SM Wilson coefficient [115] whereη B = 0.8393 accounts for QCD-corrections, and S 0 (x t ) = 2.35 is the Inami-Lim function for the SM top quark box diagram, with x t = m 2 t /m 2 W . For ∆S = 2 kaon oscillations K − K, the relevant SM operator is given by with the corresponding SM Wilson coefficient [115] where η tt = 0.5765, η cc = 1.39 (1.29 GeV/m c ) 1.1 , and η ct = 0.47 take into account QCD corrections. The Inami-Lim functions including the charm quark are given by The exchange of the diquarks X dd and X ud gives rise to the following operator [25] where s) for B s , B d , K oscillations, respectively. With X dd being a colour sextet, it contains flavour symmetric couplings to a down quark pair and can therefore mediate ∆F = 2 neutral meson mixing at tree level as well as at one-loop level. However, in the particular case where the coupling matrix f dd in Eq. (38) is diagonal, the one loop contributions vanish, as can be quickly verified by considering the one-loop diagram in Fig. 11 and replacing u by d and X ud by X dd in the loop. The exchange of X dd at tree-level can be related to the Wilson coefficient Furthermore, X ud cannot contribute via a tree-level contribution to the meson mixing. However, in contrast to X dd , it can induce neutral meson mixing at one-loop level even if one starts with a diagonal structure for f ud in Eq. (38) (written in the flavour basis), for instance in the presence of a right-handed analog of the CKM mixing matrix (see e.g. discussion in Refs. [116,117]). The exchange of the X ud in the one-loop box diagrams lead to [25] C NP where we have used the fact that f ud ij is symmetric and we have defined f ud ij = (V R ) ik f ud kj , with V R being the right-handed quark mixing matrix diagonalising the right-handed quark charged current (similar to the CKM matrix for left-handed currents).
Following the prescription of [115,[118][119][120], we define the ratio of the total contribution to the SM for M − M oscillations as where ∆ M is an experimentally determined complex parameter that depends on the respective meson M . It is experimentally constrained to ∆ B s(d) = 1.11 +0.96  [118,119]. Using the above 95% CL experimental limits, we derive the relevant constraints summarised in Tab. 5. In addition to the neutral meson oscillations, the diquark coupling can also give rise to a number of non-leptonic rare meson decays at tree-or loop-level, however the relevant constraints obtained are comparably weaker [25,105].
Generally, the latest experimental data on neutral meson oscillations provide some of the most stringent constraints on the X dd scalar diquark couplings with masses around the TeV scale (see Tab. 5). Therefore, the constraints from neutral meson oscillations will be of great importance in synergy with the direct collider searches for probing low-scale baryogenesis scenarios.

Dinucleon decay
In addition to n-n oscillations, the dinucleon decays can also provide useful constraints for diquark couplings to the SM quarks. For the diquark couplings to the first generation of SM quarks both n-n oscillations and dinucleon decay, e.g. nn → π 0 π 0 , can occur at tree-level 9 . Given the availability of better numerical estimates for the well studied transition matrix elements for n-n oscillations as compared to potentially large hadronic uncertainties for dinucleon decay matrix elements, the former provides more reliable constraints in this case. However, for the scenario where the diquark couples dominantly to the third generation SM quarks, the dinucleon decay modes, e.g. nn → π 0 π 0 (induced at the two-loop level), provide a more stringent constraint as compared to the n-n oscillation (induced at the three-loop level) 10 . In particular, for the scenario of low-scale baryogenesis, where both X dd and X ud masses lie around TeV scale, the dinucleon decays are of particular interest if the diquark couples dominantly to the third generation SM quarks. Even though we will mainly focus on the simplest case of flavour diagonal and universal couplings of diquarks to SM quarks for the study of baryogenesis, below we briefly discuss the case where the dinucleon decay mode nn → π 0 π 0 can occur at two-loop level, when the diquark couples dominantly to the third generation SM quarks. A representative Feynman diagram inducing such a mode is shown in Fig. 12. Given that in our simplified model the diquarks can exclusively couple to right-handed quarks one would require mass insertions for top and bottom quarks in the loops which we shall ignore in the following estimate in the interest of simplification and due to the large uncertainties involved in estimating the transition matrix elements. The relevant dinucleon decay rate is given by [54] where ρ n 0.25 fm −3 is the nucleon density, and m n the neutron mass. The Wilson coefficient C DN induced at the two-loop level is given by Figure 12: A representative diagram giving rise to the dinucleon decay mode nn → π 0 π 0 at two-loop level with the scalar diquarks dominantly coupling to third generation quarks.
where under the simplifying assumption of vanishing external momenta the relevant two-loop integral I can be written as Using partial fraction and rescaling all the masses by the scalar mass m X dd , the integral I can further be expressed as [122,123] where the integral J α,β;µ,ν is defined as with r {α,µ} = (m F {α,µ} /m X dd ) 2 for fermionic and t {β,ν} = (m B {β,ν} /m X dd ) 2 for bosonic masses. Note that the momenta k {1,2} are rescaled with respect to Eq. (105). Now, the integral J α,β;µ,ν can be reduced in terms of the standard integral as where in general D = 4 + dimensions we have which includes an infinite and a finite piece [122]. After collecting the different terms in Eq. (106), the divergent pieces in g(a, b) drop out, leaving the finite piece with the argument of the Spence's function containing the parameters In order to estimate the dinucleon decay rates numerically for the relevant matrix element of the dimension-15 operator O 15 DN in Eq. (103), we take | π 0 π 0 |O 15 DN |nn | ∼ Λ 11 , where Λ corresponds to a scale between Λ QCD and m n . To be conservative we take the most stringent possible constraint corresponding to Λ = m n .
The most stringent current experimental constraint on the nn → π 0 π 0 partial lifetime due to Super-Kamiokande is τ > 4.04 × 10 32 years [124], which improved over the earlier limit from the Frejus experiment of τ > 3.4×10 30 years [125]. The future experiments, such as Hyper-Kamiokande, is expected to improve on the current results by up to an order of magnitude [126]. In Fig. 13 we show the relevant current and future experimental constraints on the parameter space spanned by the coupling and diquark masses. The current constraints are relevant for masses of X ud and X dd around a few TeV. For more robust constraints, a more accurate estimation of the relevant matrix elements is desirable.

Colour preserving vacuum
From a phenomenological point of view, in the absence of a specified UV completion, the effective trilinear coupling of the form µX dd X ud X ud (e.g. induced by the vacuum expectation value of ξ in Eq. (38)) can be constrained by the requirement of colour preserving vacua. As shown in Fig. 14, the effective trilinear coupling of the form µX dd X ud X ud leads to effective quartic interactions for the X ud and X dd fields given by The relevant induced effective corrections are given by [127] In the particular limit in which two masses are similar (m X dd ∼ m X ud ) the following relation holds With the effective couplings being negative (in the absence of a UV embedding), the tree-level Lagrangian must contain terms similar to the effective quartic couplings λ eff , λ eff , λ eff , which are positive and greater in their absolute value to ensure a stable colour preserving vacuum.
Moreover, the requirement that the theory has to be perturbative (λ eff , λ eff , λ eff < 1) imposes the following constraints on the µ parameter: which provide important benchmark constraints fixing the hierarchy between m X ud , m X dd and µ = λv in Eq. (38). Consequently we fix λv ∼ m X dd to a few times m X dd for the subsequent numerical study of the baryogenesis scenarios where the hierarchy m X dd > m X ud is maintained but with both m X ud and m X dd not too far from each other to ensure a stable colour preserving vacuum. For the high-scale scenario on the other hand, a given UV completion must take into account the relevant Higgs multiplet in which the scalar fields are embedded. In the presence of a field at a particularly low scale, one must assure that the relevant interactions involving that field in the low-energy theory preserves the vacuum with a perturbative coupling. For example, when m X ud is around few TeV and m X dd is around the GUT scale, then for a low-energy theory at a temperature T < m X dd one can integrate out X dd , leaving all its effects encoded in the couplings of X ud . In order to make sure that in such a scenario the low-energy theory is in the perturbative regime, one requires λ eff < 1 in Eq. (113). This is ensured by the condition Figure 14: Diagrams giving rise to the effective quartic interaction terms of the fields X ud and X dd .
Besides, the couplings of the heavy state X uu with X dd can effectively ensure that the other radiative corrections λ eff and λ eff in Eqs. (114) and (115) do not exceed the tree-level terms [27]. To conclude, for the high-scale scenario the conditions for a colour preserving vacuum can be naturally satisfied in a UV completion where the diquark field X ud lies around TeV scale and the other diquark fields X dd and X uu are significantly heavy with masses around the GUT scale. On the other hand, for the phenomenologically driven low-scale scenario, where X ud and X dd masses are not fixed by a UV completion, one requires a very mild hierarchy between m X ud , m X dd and λv to satisfy the conditions from colour preserving vacua: λv ∼ O(few) × m X dd ∼ O(few) × m X ud . Therefore, we use the benchmark choices λv = 1.2m X dd and m X dd = 3m X ud for our subsequent numerical analysis of this type of scenario to satisfy the conditions given in Eq. (116).

Comments on other possible constraints
In addition to the constraints discussed above, X dd diquark couplings can also contribute to the neutron electric dipole moment (in the presence of complex phases in diquark couplings) and to electroweak precision observables (e.g. related to Z → d cdc ) as discussed in [105]. In the presence of X dd couplings to both left-and right-handed down-type quark pairs these loop observables can in particular be sizeable due to chiral enhancement and in particular when heavy quarks can be in the loop due to flavour antisymmetric couplings. However, given that in our scenario X dd couples only to the right-handed quarks and flavour symmetrically, one would require additional quark mass insertions to realise such observables leading to both chirality and loop suppressed contributions.

T (GeV) Process
High-scale Low-scale 10 13 Gauge interactions + t Yukawa > 10 13 QCD sphaleron ∼ 10 13 EW sphaleron + b Yukawa ∼ 10 12 τ Yukawa ∼ 10 12 c Yukawa ∼ 10 10 µ, s Yukawa 10 8 u, d, e Yukawa Table 6: Approximate temperatures at which different SM processes come into equilibrium in the early Universe [87]. Green ticks mark the processes that are assumed to be in equilibrium in the high-and low-scale scenarios, and red crosses mark the interactions that are assumed to not be in equilibrium.

Results and discussion
Using the formalism described in Sec. 3, in this section we present our findings regarding the viability of successful baryogenesis for the simplified model described in Sec. 3.2.1, and confront the results with the experimental constraints discussed in the previous section. Hereby, we distinguish two distinct scenarios, which we will refer to as the high-and low-scale scenarios. Assuming the mechanism that generates the BAU is the out-of-equilibrium decay of the X dd diquark, as described in Sec. 3, both the high-and low-scale scenario can lead to a baryon asymmetry close to the observed value, but with very different detection prospects. We furthermore define, for both scenarios, a model-dependent CP-violation parameter , originating from the interference between the tree-level decay of X dd and the 1-loop decay with X dd , another generation of X dd , as introduced in Eq. (42), and simplified here by selecting as a benchmark scenario 11 , f dd = f dd and λ = λ, such that where x = (m X dd /m X dd ) 2 , and r is the branching ratio for X dd → X * ud X * ud . Hereby, to simplify the forthcoming baryogenesis parameter space analysis, we will assume a flavour diagonal and universal structure for (in general 3×3 coupling matrices) f ud(dd) ij (neglecting any flavour effects for baryogenesis) and henceforth will denote them as f ud(dd) for brevity. The maximal value of is = 2r and is the most optimistic scenario, for the other cases we select as a benchmark scenario x = 0.2. For the numerical evaluation of the final baryon asymmetry in two different baryogenesis scenarios, we take into account the appropriate SM processes in equilibrium, as tabulated in Tab. 6. The green ticks show the SM processes which are taken into account at the corresponding temperature while solving the Boltzmann equations. During the computation of the baryon asymmetry in the high-scale baryogenesis scenario, we assume that only the Yukawa couplings of the third generation quarks are in chemical equilibrium, while all three generations are assumed to be in chemical equilibrium in the low-scale baryogenesis scenario. All quark generations are considered to be in thermal equilibrium in both scenarios. Furthermore, we also assume that during high-scale baryogenesis, only the decay of X dd to third generation quarks are relevant, since other two lighter generations do not come into equilibrium with the thermal bath of the Universe until much later than the temperature relevant for high-scale baryogenesis. On the other hand, for the low-scale baryogenesis scenario, we consider the decay of X dd to all three generations of quarks. In the following, we proceed to discuss the two scenarios in more detail and present the numerical results of the baryogenesis mechanism of Sec. 3.
High-scale scenario: In the high-scale baryogenesis scenario, the X dd diquark has a mass O(10 13−14 ) GeV, while the mass of X ud is around O(10 3−4 ) GeV. From a theoretical perspective, the large hierarchy between the masses of X dd and X dd is naturally motivated by gauge coupling unification in the context of a GUT embedding of the model as discussed in App. A. Having a mass near the GUT scale, X dd would remain inaccessible at the direct collider searches; however, conspiring together with a TeV scale X ud it can still lead to an observable rate of n-n oscillations, therefore making the scenario phenomenologically exciting. For instance, taking λ ≈ O(1) × (m X dd )/v and f ud,dd ≈ O(0.1 − 1) together with the above discussed combination of masses, one finds that the n-n oscillation rate is very close to the current experimental bound. On the other hand, since m X ud is O(10 3−4 ) GeV, it can potentially be directly produced at the LHC and future colliders providing complementarity between collider searches and n-n oscillation in constraining parts of the parameter space for a successful high-scale baryogenesis scenario.
For illustration, we present an example of the typical dynamics of the baryon asymmetry generation for the high-scale baryogenesis scenario in Fig. 15. We show the evolution of η ∆B with respect to z, for the benchmark choices f ud = f dd = 0.05, m X ud = 5 TeV and λv = 1.5m X dd . The plots in the top row correspond to the case of zero initial abundance of X dd , the plots in bottom row corresponds to the case of an initial thermal abundance of X dd . The plots in the left and right columns correspond to the choices m X dd = 10 14 GeV and m X dd = 10 13 GeV, respectively. Furthermore, a final baryon asymmetry greater than the one observed in the Universe is achieved for a maximal CP violation = 2r, as well as for the case where is given by Eq. (118) with x = 0.2. The evolution that occurs for m X dd = 10 14 GeV (left column) and m X dd = 10 13 GeV (right column) are similar, with the difference that the asymmetry generation starts earlier and gets more severely washed out in the latter case, since both decay rates and washout processes are less suppressed with respect to the Hubble rate. Whether the initial number density of X dd is in equilibrium or not has little effect on the final baryon asymmetry, as can be seen by comparing the top and bottom plots in Fig. 15.
To illustrate the evolution of the relevant decay and scattering processes leading to the evolution of the baryon asymmetry with z in the high-scale scenario, we show Fig. 16 the decay rates and different scattering rates for the same choice of parameters as in the plots in the left column of Fig. 15. We note that in the high-scale scenario, quark and X ud mediated washout processes S are subdominant and the X dd mediated X s(t) processes provide the most dominant washout for z 10. This can be understood as follows. The quark and X ud mediated scatterings S (0) s(t) contain heavy X dd as one of its external legs and leads to large Boltzmann (or phase space) suppression of these scattering rates for T m X dd . On the other hand, the X dd mediated X s(t) scatterings are much less-severely suppressed due to X dd mass in the propagator for T m X dd . One particular caveat to the above discussion is the situation when f ud > f dd and the coupling f ud is large, O(1). In such a case the washout processes S t , S s , S 0 t and S 0 s processes can still be quite effective for large z and can potentially be comparable to X s(t) . As a benchmark, as long as the coupling f dd does not exceed values of O(0.1) the dominant washout processes X s and X t remain modest leading to the generation of a sizeable final baryon asymmetry. With the washout processes being modest in the high-scale scenario set-up for f dd O(0.1), a baryon asymmetry close to the observed value can be found over a wide parameter range of phenomenological interest, as shown in Fig. 17. This  Figure 16: Equilibrium decay and scattering rates with respect to z = m X dd /T and f ud = f dd = 0.05 for the high-scale scenario with m X dd = 10 14 GeV, m X ud = 5 TeV and λv = 1.5m X dd . The solid orange line shows the decay rate of X dd , while the solid (dashed) blue, green, and purple lines show the scattering rates S s (S t ), X s (X t ), and S 0 s (S 0 t ) respectively. For the notation see Sec. 3.
statement holds for both cases, a vanishing initial X dd number density, and a non-zero initial X dd number density in equilibrium. In Fig. 17, we show the parameter space of the high-scale baryogenesis scenario in the m X ud -f ud plane with f ud = f dd . The generated final baryon asymmetry is shown in yellow with deeper shades corresponding to higher yields. The red lines indicate the observed value 12 . Note that in Fig. 17, we take the diquark couplings f ud and f dd to be equal, ranging from f ud = f dd = 10 −3 to unity. In general f ud and f dd couplings could be treated as independent parameters; however, we note that for a large hierarchy between f ud and f dd , the connection of baryogenesis to the relevant experimental observables gets weakened. This is because on the one hand the LHC direct searches for X ud and n-n oscillation rate largely depend on a high value for the coupling f ud , on the other hand the baryon asymmetry depends to a large degree on f dd . Therefore, we find that fixing the ratio of f ud and f dd to unity in Fig. 17 illustrates the interplay between the final baryon asymmetry and the reach of experiments more clearly. While in Fig. 17 the tri-linear scalar diquark coupling is taken to be λv = 1.5m X dd , we have found that a lower or higher value can also impact the final baryon asymmetry. A lower value for λv would generally result in a larger baryon asymmetry at the cost of lowering the experimental reach of the n-n oscillations, in contrast to a higher value which can run into problems with the constraints from colour preserving vacua (see Sec. 4.4). However, if we relax the constraints from colour preservation of the vacuum, then a higher value of the trilinear coupling λv would lead to an increased decay and washout rate, leading to an asymmetry generation at later times and a more effective washout, resulting in a lower final asymmetry. In all four plots of Fig. 17, the baryon asymmetry is generally lower for very high values of f ud (i.e. when this coupling is O(1)) as the consequence of an increased washout from X s and X t . This effect is more pronounced in the right plots, as m X dd is an order of magnitude lower as compared to the left plots.
For lower values f ud O(0.1), the asymmetry again gets smaller for smaller couplings in the top plots of Fig. 17, while it remains rather constant in the bottom plots. This effect is due to the CP-asymmetry parameter being smaller for small couplings in the top plots 13 . This is in contrast to the plots in the bottom row, where the maximal CP asymmetry = 2r is used. The reach of n-n oscillations and collider experiments, which provide important constraints for the high-scale baryogenesis parameter space, do however, depend on the mass scales m X dd and m X ud . As can be seen in the left column of Fig. 17, for m X dd = 10 14 GeV, the LHC proves to be a better probe than n-n oscillation experiments in parts of the parameter space. In contrast, for a slightly smaller mass, m X dd = 10 13 GeV, (Fig. 17, right column) n-n oscillations almost completely dominate in sensitivity. Furthermore, considering that both the LHC and n-n oscillation rate depend non-trivially on the smallness of m X ud and the generated baryon asymmetry is largely insensitive to small changes in m X ud , a large part of the parameter space in which baryogenesis is successful (e.g. when m X ud O(10) TeV) remains unreachable by any current or near future experiments. Note that the final baryon asymmetry shown in Fig. 17 was calculated under the assumption that X dd and X ud dominantly decay into to third generation quarks, while the n-n oscillation requires X dd and X ud couplings to first generation quarks. In general both, the SM quark Yukawa coupling matrix and the X dd (and X ud ) Yukawa coupling matrix need not to be diagonal in the same basis. Moreover, the inverse decays and washout processes can lead to the generation of an asymmetry in all flavours even if the initial asymmetry is generated in only one of the flavours. If an initial asymmetry is generated in a the third flavour of quarks together with an equal and opposite asymmetry in an unflavoured X ud (comparable to the Higgs asymmetry generated in the leptogenesis [49]), then through inverse decay or washout processes the X ud asymmetry will induce asymmetry in all three quark flavours realising a thermal contact among flavours such that asymmetry in one flavour induces an asymmetry in the other flavours. Therefore, the asymmetry that will "leak" into the first generation will be subject to washout effects from the first generation Yukawa couplings of X dd and X ud inducing n-n oscillations. In a full two-flavour treatment of the high-scale scenario we expect that the first generation interactions relevant for n-n oscillation operators will therefore at least partially washout the asymmetry along the perpendicular direction to the third generation quarks.
To this end, the synergy of LHC searches for a tb vs. dijet final state searches with the n-n oscillation will provide valuable hint towards the underlying flavour dynamics of the high-scale scenario. A comprehensive analysis of these interesting flavour effects is beyond the scope of the present work and will be addressed in the future.
In summary, the high-scale scenario is a case where baryogenesis can be successful, while at the same time, the baryon-number-violating mechanism can be potentially discovered by n-n oscillation experiments. Moreover, the X ud diquark is close to the current reach of the LHC. In case a signal is observed in either of these experiments, the other experiment can be used to set further constraints, perhaps rejecting or confirming this mechanism as the one responsible for baryogenesis.
Low-scale scenario: In the low-scale baryogenesis scenario, we consider the case where the diquark X dd mass is a few times the mass of X ud , with m X ud varying in the range O(10 3−9 ) GeV. The phenomenological choice of a mild hierarchy between m X dd and m X ud is primarily driven by the constraints from colour preserving vacua as discussed in Sec. 4. To illustrate the typical dynamics . The rest of the parameters are fixed to the benchmark values m X ud = 10 8 GeV, m X dd = 3m X ud and λv = 1.2m X dd . The evolution of the baryon asymmetry for the maximum case (orange curve) as well as corresponding to Eq. (118) with x = 0.2 (red curve) are shown in each plot. An illustration of the evolution of the decay rates and scattering rates in the low-scale baryogenesis scenario is shown in Fig. 19, using the same parameters as in the the left plot in Fig. 18. Similar to the high-scale case (cf. Fig. 16), Fig. 19 shows for the low-scale case the overall dominance of the decay rate until some high z when the X s(t) washout takes over. In Fig. 20, we present the parameter space for low-scale baryogenesis scenario in the m X ud -f ud plane for the maximum value of , being = 2r. The final baryon asymmetry is shown in yellow contours, where the low-m X ud end of the plots is stretched out in order to show the complementary sensitivity of the dinucleon decay, neutral kaon oscillation, and collider searches with respect to n-n oscillations. We show two different coupling hierarchies: f ud = f dd (Fig. 20 top) and f ud = 10f dd (Fig. 20 bottom). The tri-linear coupling is chosen as λv = 1.2m X ud , and the mass ratio between X dd and X ud is taken to be m X dd /m X ud = 3. Our findings suggest that in the low-scale baryogenesis : Equilibrium decay and scattering rates with respect to z = m X dd /T and f ud = f dd = 0.05 for the low-scale scenario with m X dd = 3 × 10 8 GeV, m X ud = 10 8 GeV, and λv = 1.2m X dd . The solid orange line shows the decay rate of X dd , while the solid (dashed) blue, green, and purple lines show the scattering rates S s (S t ), X s (X t ), and S 0 s (S 0 t ) respectively. scenario, even for the case of maximal CP violation, successful baryogenesis can only occur when the mass scales m X dd ∼ m X ud are sufficiently heavier (at least O(10 8 ) GeV) than the electroweak scale such that the washout processes remain largely subdominant as compared to the baryon asymmetry generation through X dd decays. Since in this scenario m X dd is not very heavy as compared to m X ud , the quark and X ud mediated washout processes S (0) s(t) do not receive any relative suppression with respect to X dd mediated X s(t) processes, in contrast to the high-scale scenario where the heaviness of X dd leads to Boltzmann suppression of the washout processes S (0) s(t) as compared to X s(t) . For large f ud and f dd couplings and small X dd (X ud ) masses, the washout processes become too strong to generate any sizeable asymmetry. In particular, for a fixed hierarchy of X ud and X dd masses (e.g. our benchmark choice m X dd = 3m X ud ) the final baryon asymmetry gets reduced for smaller masses m X ud , as can be noticed in Fig. 20. This can be understood to be due to X dd being less massive (being related to X ud masses by a fixed hierarchy), thereby falling out-of-equilibrium at a later time when the washout rates are less suppressed with respect to the Hubble rate. Furthermore, a larger coupling f ud should lead to an increased washout via X s and X t , and therefore a lower final baryon asymmetry; however, since the X s and X t washout processes also depend on the coupling f dd , such an effect is noticeably less strong in the bottom plot in Fig. 20, in which f dd is an order of magnitude smaller as compared to the top plot. On the other hand, a higher value of λv would (as in the high-scale scenario) lead to an increase in both, the rate of the decay for X dd as well as in the scattering rate for washout processes S s(t) and X s(t) . While this leads to more asymmetry being generated at an earlier time, it finally would lead to a larger washout for large z once the baryon asymmetry generation freezes out, leading to a smaller final baryon asymmetry.
From Fig. 20, it is very clear that if m X dd ∼ O(few)m X ud lies below O(10 8 ) GeV then the washout processes are too strong to generate any sizeable asymmetry, even if the CP violation is maximal. Therefore, an observation of X dd and X ud at the LHC together with a signal for n-n oscillations would imply that the low-scale baryogenesis scenario is completely ruled out as a mechanism behind the observed baryon asymmetry. On the other hand, for m X dd ∼ O(few)m X ud > O(10 8 ) GeV the low-scale baryogenesis mechanism remains a viable option, however, the small couplings required together with the mass scale make this part of the parameter space inaccessible to all current and near-future experimental efforts. In this region, all washout processes involving scatterings are small because of the Boltzmann suppression due to the heaviness of X dd (and X ud ) and the smallness of the couplings. For instance, the dominant washout processes X s and X t are small due to the Boltzmann suppression coming from the two X ud in the outer legs. Therefore, an initially generated asymmetry due to the CP violating decay of X dd can survive. In the contrary, if the hierarchy between m X dd and m X ud is larger (m X dd /m X ud 4), then the Boltzmann (or phase-space) suppression of the washout processes (e.g. X s and X t ), due to X ud in the external legs, get lifted making the washout very strong, and consequently the baryon asymmetry gets rapidly wiped out 14 .
In case the CP-violating parameter would be given by Eq. (118) rather than being fixed to the maximum value as in Fig. 20, we find that the entire parameter space of the low-scale baryogenesis scenario cannot reach the observed baryon asymmetry of the Universe. The asymmetry is particularly small for lower values of the coupling f dd , as a consequence of less asymmetry being initially generated due to being smaller. Apart from this effect, given by Eq. (118) shows a similar behavior as in Fig. 20 for the low-scale scenario with respect to a decreased asymmetry for lower masses m X ud . This is caused by the washout rates being larger than the Hubble rate for low masses of X dd .
To summarise, if m X dd ∼ O(few)m X ud is less than O(10 8 ) GeV, then the baryon number washout processes are too strong to generate any sizeable asymmetry even if the CP violation is enhanced to its theoretically allowed maximal value. Consequently, the possibility of a low-scale baryogenesis to explain the correct baryon asymmetry can be ruled out in this scenario if both the diquark states X dd and X ud are discovered at the LHC or future colliders together with a positive signal for n-n oscillations. On the other hand, successful baryogenesis can be achieved for some parts of the parameter space of the low-scale baryogenesis scenario, however, then the corresponding masses of both diquarks are generally too high and the couplings are too small to probe such a scenario at the LHC, dinucleon decays or n-n oscillation searches. 1 Figure 20: Low-scale scenario: Yellow shading indicates the final baryon asymmetry η ∆B at the time of electroweak symmetry breaking in the m X ud -f ud plane. The diquark couplings are chosen as follows: Top: f ud = f dd . Bottom: f ud = 10f dd . In both plots, we chose λv = 1.2m X dd , m X dd = 3m X ud , and the maximal CP asymmetry = 2r is used. The red line indicates the observed value of the baryon asymmetry η obs B , while the orange, green, black, and blue lines correspond to experimental constraints coming from dinucleon decay, meson oscillations, the LHC, and n-n oscillations respectively. For the dinucleon and n-n oscillation constraints, solid lines represent current experimental bounds, while dashed lines show projected future sensitivities. For the meson oscillation constraints, dotted and dot dashed lines correspond to tree-level oscillation via a X dd mediator, and box-diagram oscillation involving an X ud diquark, respectively (see Sec. 4 for details). For the LHC limit CMS data with √ s = 13 TeV and 36 fb −1 integrated luminosity was used [107].

Conclusions
Within the SM, baryon number is an accidental global symmetry and only violated at finite temperature via non-perturbative sphaleron interactions. On the other hand, new physics is required in order to explain the observed baryon asymmetry via a mechanism that complies with the three Sakharov conditions by violating B − L, C and CP while leading to a departure from thermal equilibrium. Hence, it is natural to search for new B-violating interactions that might have far-reaching implications on the mechanism behind the baryon asymmetry. While the dim-6 operators with |∆B| = 1 and |∆(B − L)| = 0 are tightly constrained by limits from proton decay, dim-9, |∆B| = |∆(B − L)| = 2 interactions are less stringently constrained. However, in the future, new facilities such as the Deep Underground Neutrino Experiment and the European Spallation Source, will probe those interactions via n-n oscillation experiments to an unprecedented sensitivity and make this possibility timely to be reconsidered. As a discovery of n-n oscillations would directly imply baryon-number violation, it might have a tight link to the mechanism of baryogenesis. Moreover, its oscillation time τ n-n would give us a hint on the scale of new physics Λ NP of the effective interaction.
In this work, we have studied the implications of such a discovery for baryogenesis in a comprehensive manner. First, we performed a model-independent EFT analysis of the new dim-9, |∆B| = 2 operator uudddd. Assuming no inherent CP violation source connected to this new operator, its interaction would lead to washout processes that could wipe out a potential pre-existing baryon asymmetry. Under the assumption of a discovery with a given oscillation time τ n-n , we estimated the corresponding NP scale Λ NP and calculated the corresponding washout strength Γ W /H. In order to relate the low scale, where n-n oscillations are expected to occur, with the high scale, where the NP is active and the washout takes place, we considered the running of the Wilson coefficients from the scale of NP down to the scale µ 0 = 2 GeV, where the nuclear matrix elements are provided by lattice calculations. We demonstrate that the observation of n-n oscillations at future experiments such as DUNE or NNBAR, would imply a strong washout down to around 100 TeV. With the LHC already probing scales up to 1 − 10 TeV, we would predict, under the made assumptions, new physics between 10 − 100 TeV, possibly detectable at a future 100 TeV collider.
In order to accommodate the possibility of different hierarchies of NP within the effective operator uudddd and an additional new source of CP violation, we explored a simplified set-up that realises one of the two possible UV-complete topoligies, featuring two new types of scalar diquarks X ud and X dd with couplings f ud and f dd to SM quarks, respectively. After B − L breaking, both carry a baryon number −2/3 such that the out-of-equilibrium decay of the heavier one, X dd → X ud X ud , violates baryon number and generates (together with a new CP phase) a baryon asymmetry. Due to CPT invariance, also the decay X dd → d c d c contributes indirectly to this asymmetry. In order to describe the final amount of the baryon asymmetry, we presented a comprehensive framework of Boltzmann equations.
Moreover, we discussed in detail all experimental constraints which this set-up could be subject to. We performed a dedicated study of current LHC, but also future HE-LHC and FCC-hh limits on the scalar sextet diquarks. While current LHC limits exclude these new particles for f ud,dd ≈ O(1) up to around 10 TeV, smaller couplings (f ud,dd ≈ O(0.01)) allow for masses heavier than 6 − 7 TeV. Depending on the mass scale of the diquarks, meson oscillations can strongly constrain the respective couplings. As X ud can contribute to meson oscillations only via loop diagrams, limits on the coupling f ud are generally weak. In contrast, due to the direct tree-level contribution of X dd , the coupling f dd O(10 −3 ) is strongly constrained for m X dd ≈ TeV. In particular for third generation couplings and comparable masses of X ud and X dd , the dinucleon decay can set more stringent limits than n-n oscillations, as they contribute at the two-loop instead of the three-loop level. This leads with current limits from Frejus to f ud 33 < 1 for m X ud ≈ m X dd /3 ≈ 2 TeV and is expected to improve to f ud 33 < 0.3 for m X ud ≈ m X dd /3 ≈ 2 TeV or f ud 33 < 1 for m X ud ≈ m X dd /3 ≈ 4 TeV with Hyper-Kamiokande. Moreover, we discussed considerations from colour preserving vacua and commented on further limits from electroweak precision observables and electric dipole moments from which we do not expect any major constraints due to chirality and loop suppression.
Finally, we compared the prospects to generate the observed baryon asymmetry with the allowed parameter space from all considered experimental constraints. In our first scenario, we assumed a heavy m X dd ≈ 10 13 GeV and a lighter m X ud ≈ O(TeV). Depending on the actual size of the CP violation, the observed baryon asymmetry can be produced with respective couplings, e.g. f ud = f dd ≈ 0.3 for = 0.01 or f ud = f dd ≈ 0.4 for maximal . In this scenario, an interesting interplay between the LHC and future n-n experiments can take place: In order to account for the observed baryon asymmetry, a discovery of a diquark at the LHC tb final state would still allow for a successful high-scale baryogenesis scenario. If the diquark is also discovered at the LHC in the lighter generation quark dijet final states, a signal at n-n oscillation experiments would give valuable insight into the flavour dynamics of a possible high-scale baryogenesis scenario, and would motivate a comprehensive study of flavour effects. However, if a corresponding signal at both the LHC and n-n oscillation experiments is discovered and naively applying the maximal washout corresponding to the first generation n-n oscillations to our single flavoured analysis suggests that the high-scale scenario will still remain a viable option.
In our second scenario, we analysed a less hierarchical set-up with m X ud = m X dd /3. In this case, the observed baryon asymmetry can only be generated with m X ud = m X dd /3 ≈ 10 5 TeV, which is out of reach of all experiments. This confirms our naïve EFT estimate of a strong washout Γ W /H even in the presence of a CP violating source at low scales. However, we have shown that the future NNBAR experiment will set stringent limits on the parameter space ranging from f ud = f dd = 1 for m X ud = m X dd /3 ≈ 700 TeV to f ud = f dd = 0.001 for m X ud = m X dd /3 ≈ 10 TeV. This demonstrates the excellent future prospects of this experiment, making it possible to directly observe |∆B| = 2 interactions up to an unprecedented scale.
In summary, the discovery of n-n oscillations can have far-reaching consequences on our understanding of physics beyond the Standard Model. Not only would it point to new B-violating interactions, it could also have significant implications on the mechanism behind the baryon asymmetry. The studied less hierarchical set-up confirms the naïve EFT estimate of a too strong washout in order to generate the observed baryon asymmetry in case of an observation. Moreover, it even demonstrates that a possible CP violation of the new low-scale interactions could not balance out this strong washout. Hence, in order to generate the observed baryon asymmetry, we would expect some other new physics between the LHC exclusion and the scale of strong washout in case the |∆B| = 2 interaction does not involve any new CP-violating interaction. In contrast, the analysis of our high-scale scenario demonstrates that for a large hierarchy of new degrees of freedom that feature additionally a new CP phase, absorbed in the effective |∆B| = 2 interaction uudddd, we are actually able to generate the observed baryon asymmetry in case of some new physics being at a low scale. This would imply that an observation of n-n oscillations might point towards new physics experimentally reachable at the LHC or future colliders. (1, 1, 2, 1) (1, 1, 1) ⊕ (1, 1, 0)   Tables  8, 9 and 10, we present the decomposition of the three representations 45, 54 and 126 of SO(10) via the breaking chain G PS → G LR → G SM , which will be essential for our discussion. In what follows we will briefly highlight how the criterion of gauge coupling unification can naturally motivate one of the choices for diquark mass scales that we subsequently use to study the scenario of high-scale baryogenesis. Before moving onto the discussion of gauge coupling unification let us briefly summarise the possibility of intermediate symmetries which can be realised by choosing appropriate scalar multiplets (which will also contain the scalar diquark fields) relevant for realising a given symmetry breaking chain. The possibility of realising an intermediate gauge symmetry group is dictated by the choice of Higgs multiplets of SO(10) used for symmetry breaking. In particular, a real 45, a real 54, or a real 210 can be used along with a complex 126 to achieve the breaking to the SM gauge group. Using the 210 leads to the intermediate symmetry G PS but with a broken discrete parity symmetry (subsequently, leading to an asymmetric gauge couplings for SU (2) L(R) : g L = g R ), often referred to as D-parity, while 54 can break SO(10) down to G PS symmetry preserving D-parity. On the other hand, the intermediate symmetry corresponds to a left-right symmetric gauge group G LR with a broken D-parity [128].
Effective (B − L)-violating interactions can naturally occur in GUT theories like SO(10), G PS or G LR , where (B − L) is part of the local gauge symmetry and broken by the vacuum expectation value of a scalar field carrying non-vanishing (B − L) charge. In SO(10) and G PS the effective (B − L)-violating couplings can arise, when the (B − L)-symmetry is broken by giving a vacuum expectation value to the complete SM singlet field. For instance, it can be contained in the 126 H multiplet of SO(10) transforming as (1,1,0) under the SM gauge group G SM while transforming as (1,3,10) under G PS , see bold-print in Table 10. Note that this field, denoted by ξ in our Lagrangian in Eq. (38) carries (B − L) = −2. In SO(10) GUT ξ can potentially also generate large Majorana masses for the right-handed neutrinos through the couplings of the form ν c ν c ξ. Note that we assume  where C 2 (G) corresponds to the quadratic Casimir operator for the gauge bosons in their adjoint representation T (R f ) and T (R s ) correspond to the Dynkin indices of the irreducible representation R f,s for a given fermion and scalar field, respectively, e.g.
T (R f,s ) ≡ and d(R f,s ) is the dimension of the representation R f,s under all gauge groups except the i-th gauge group under consideration. An additional factor of 1/2 is multiplied to T (R s ) in the case R s is a real Higgs representation. Let us now consider a SO(10) GUT scenario where X ud is the only light field (m X ud ∼ O(TeV) scale) affecting the running of gauge coupling constants at one-loop level (in addition to the SM field content), while the remaining diquarks (X dd and X uu ) and other new fields lie around the unification scale. This scenario is particularly interesting because it allows for a high-scale baryogenesis scenario to be realised through the decay X dd → X * ud X * ud as discussed in Sec. 3. Additionally, a relatively light X ud can lead to a sizeable contribution to the n-n oscillation rate and can be directly searched for at the LHC and future colliders. To simplify the discussion, let us first assume that (if present) the intermediate gauge symmetry G PS also lies very close to the SO(10) unification scale. In such a scenario where X ud is the only light field with a mass around O(TeV), the gauge couplings do not unify. However, if one also includes two copies of ∆ : (1, 3, 0) around a O(TeV) scale a successful SO(10) unification can be achieved as first pointed out in [27]. Note that the field ∆ : (1, 3, 0) can be obtained from a 45 or 54 representation of SO(10), see Tab. 8 and 9. In Fig. 21 (left panel), we show the evolution of the gauge couplings for m X ud ∼ m ∆ = 3 TeV. With increasing m X ud ∼ m ∆ (beyond a few TeV) a unification triangle starts to develop (not visible in the figure due to the scaling). However, subject to the threshold corrections, unification can still be considered as a marginally successful. This leads to m X ud 10 TeV as a necessary condition for a successful SO(10) gauge coupling unification.
Another alternative scenario that provides successful unification is X ud as the only light field with mass O(TeV) and an additional field Σ : 6, 3, 1 3 at an intermediate scale m Σ ∼ 10 4 TeV. In Fig. 21 (right panel), we show the evolution of the gauge couplings for m X ud = 3 TeV and m Σ = 10 4 TeV.

B Chemical potential relations
In the early Universe, the chemical potential of a particle species a can be related to the ratio of number density η a to equilibrium number density η eq a , such that [42,43] where µ a is the chemical potential of the species a and T is the temperature of the Universe. Furthermore, in writing Eq. (125), we have assumed that the number density of the particle species is close to equilibrium such that |η a − ηā| η eq a , whereā is the antiparticle of a, and where the chemical potential of a particle species a is related to the chemical potential of its antiparticle via the relation µ a = −µā. If an interaction is in equilibrium, a relation can be found between the chemical potentials of the particle species involved. For example, the reaction W − ↔ū L + d L being in equilibrium leads to the relation µ W +µ u L = µ d L . At temperatures above the scale of electroweak symmetry breaking T Λ EW , the total electric charge Q and the z-component of isospin T 3 must both be zero. The latter constraint leads to µ W = 0 [42], such that (suppressing flavour indices) Similarly, from the reactions W − + ν L ↔ e L and W − ↔ h − + h 0 , we find the relations From the Yukawa interactions h 0 +ē c ↔ e L , h 0 + u L ↔ū c , and h 0 +d c ↔ d L , we obtain and the electroweak sphaleron interaction in equilibrium leads to the relation where i denote the flavour indices. Assuming flavour universality, henceforth we replace the sum over flavours by the number of SM quark generations in thermal equilibrium N . Lastly, since we are interested in the evolution of the number density of the diquark X * ud in Sec. 3, we will consider the reaction X * ud ↔ū c +d c to be in equilibrium, thereby obtaining the relation µ X * ud = µūc + µdc.
In terms of chemical potentials, the total electric charge Q of the Universe can then be expressed as Q = 2N (µ Q + µūc) − N (µ Q + µdc) − N (µ L + µēc) + 2µ H + 2 where C X * ud is the colour multiplicity of X * ud , N is the number of fermion generations in thermal equilibrium 15 . Since the total electric charge Q and the z-component of isospin T 3 of the Universe must be zero above the scale of electroweak symmetry breaking, the individual chemical potential relations can be combined with equation (131) to relate the chemical potential of the quark singlets u c andd c to the diquark X * ud . Using the notation for any particle species a, we obtain the relations To define the net baryon number density of the Universe, we include the particle species X ud ,d c andū c as well as their corresponding antiparticles, and assume that the other SM fields do not contribute to a baryon number asymmetry generation or washout 16 . Such a prescription leads to the net baryon number difference of the Universe, cf. Eq. (46), given by Under similar assumptions, the equilibrium baryon number density then can be defined in terms of the number densities of the relevant species as n eq b ≡ 2 3 C X * ud n eq X * ud + N n eq d c + N n eq u c + N n eq u L + N n eq To derive the two above equations, we have used with ζ(s) denoting the Riemann zeta function. Combining equations (134) and (135), we finally obtain a relation between baryon-to-photon density and the chemical potential of X * For the relevant cases in our analysis, we have C X * ud = 6, and depending on the number of generation in thermal equilibrium (N ), C B takes the values C B = π 2 3 ζ(3) × with M i denoting the matrix element corresponding to the process and the symmetry factor in front with δ ab takes care of identical initial (final states), and the full reaction rate for (inverse) decays is the sum of all allowed channels, The reaction rate γ 1,2↔3,4 due to the s-and t-channel scattering 1, 2 ↔ 3, 4 is given by with s min = max((m 1 + m 2 ) 2 , (m 3 + m 4 ) 2 ), andσ (s) is given by [136] σ (s) = 1 8πs where the limits on t can be obtained as where t and s are the usual Mandelstam variables, and Γ i is the decay width of particle i. Similarly, for the diagrams meditated by X ( * ) dd the matrix elements are given by