The role of CP violating scatterings in baryogenesis - case study of the neutron portal

Many baryogenesis scenarios invoke the charge parity (CP) violating out-of-equilibrium decay of a heavy particle in order to explain the baryon asymmetry. Such scenarios will in general also allow CP violating scatterings. We study the effect of these CP violating scatterings on the final asymmetry in a neutron portal scenario. We solve the Boltzmann equations governing the evolution of the baryon number numerically and show that the CP violating scatterings play a dominant role in a significant portion of the parameter space.


Introduction
A baryon symmetric universe would result in a baryon density around ten orders of magnitude smaller than observed [1]. Deductions of the baryon asymmetry from the cosmic microwave background give a baryon-to-entropy density ratio of [2,3], Y B = (0.86 ± 0.01) × 10 −10 , (1.1) which is in agreement with the baryon asymmetry required to match observed primordial element abundances: 0.4 × 10 −10 Y B 0.9 × 10 −10 [4]. The required baryon asymmetry can be created dynamically by scenarios satisfying the Sakharov conditions [5]. One way of producing an asymmetry is through the out-of-equilibrium and CP violating decays of a heavy particle such as in leptogenesis or grand unified theory baryogenesis [6][7][8][9].
Scatterings of particles in 2 ↔ 2 type processes can also violate CP and potentially produce asymmetries. In such scenarios, the out-of-equilibrium condition can be satisfied in two ways. The first option is to have a visible and a hidden sector at different temperatures [10][11][12]; the second is to have one or more of the particles freeze out, such as the annihilating DM particles in WIMPy baryogenesis scenarios [13][14][15][16][17]. The general principles behind such a freeze out scenario have been investigated in ref. [18]. CP violation in 2 ↔ 2 processes could also provide a mechanism for producing asymmetric dark matter [19][20][21][22] in both freezein [10][11][12] and freeze-out type scenarios [18,23].
The presence of CP violation in scatterings is, however, not confined to models specifically constructed to exploit this feature: it can also appear in baryogenesis-via-decay type scenarios. For example, scatterings of the heavy Majorana neutrinos in leptogenesis also violate CP [24][25][26][27][28]. Numerical calculations show including CP violating scatterings in leptogenesis has a large effect on the baryon asymmetry at high temperature but only a negligible effect at low temperature [24]. The purpose of this paper is to study the effects of CP violating scatterings on the generation of the baryon asymmetry more generally. We choose to study a neutron portal model in which the parameters are not restricted by having to explain low energy neutrino data. Baryogenesis via CP violating decays in the neutron portal have been studied previously [29]. CP violating scatterings in a neutron portal model were discussed in ref. [23]; however, the unitarity constraint was not properly taken into account. We will extend the analysis to include a full numerical treatment of the Boltzmann equations including CP violating scatterings. Such effects are not constrained to this particular model and will play a role in baryogenesis scenarios more generally.
The layout of the paper is as follows. In section 2 we introduce our notation. In section 3 we introduce the neutron portal Lagrangian and discuss the relevant decays and scatterings. In section 4 we write down the Boltzmann equations for the evolution of the baryon asymmetry and show example numerical solutions. Constraints on the model are discussed in section 5. Possible UV completions are discussed in section 6.

Notation
In an expanding universe with no collisions, a particle number density, n X , drops as R −3 where R is the scale factor. Taking into account collisions the density evolves as: where H =Ṙ/R is the Hubble expansion rate and C(X) is the collision term, i.e. the rate of change of n X due to interactions with other particles. Let us imagine now such an interaction from state α to state β and its reverse process which changes X number by one unit. The collision term for this transition can be written as [30]: where f ψ = Exp[(µ ψ − E ψ )/T ] is the phase space density of species ψ with chemical potential µ ψ at energy E ψ and temperature T , is the normalised volume element of the three momentum and g ψ counts the degrees of freedom of ψ. As an approximation we have used Maxwell-Boltzmann statistics and so have dropped the stimulated emission factors and Pauli blocking terms appropriate for the full quantum statistics. We assume kinetic equilibrium, i.e. a common temperature T , throughout. We denote the equilibrium reaction rate density for a process α → β as: ≡n eq α1 ...n eq αn vσ(α → β) , (2.4) where vσ(α → β) is the thermally averaged cross section and f eq αi (n eq αi ) denotes the phase space (particle number) density in the absence of a chemical potential. Using Maxwell-Boltzmann statistics: where M αi is the mass of the species and K 2 (x) is the modified Bessel function of the second kind of order two. In the Maxwell-Boltzmann approximation, the non-equilibrium rate is found using the appropriate re-weighting: Unitarity of the S-matrix and invariance under charge parity time (CPT) gives a condition for the equilibrium rate densities [6,7]: where the sum runs over all possible final states and α denotes the CP conjugate of α. This ensures no particle number asymmetry can be generated without a departure from equilibrium [31][32][33].

Lagrangian
Consider the interaction Lagrangian [29], where the κ aijk are couplings with units of (mass) −2 , u Ri (d Ri ) is the right chiral up (down) type quark field of flavour i, Ψ c denotes the conjugate of field Ψ and the X α are Majorana fermions. We impose a global lepton number symmetry forbidding coupling of the X α to the standard model (SM) Higgs and lepton doublets. Spinor indices are contracted between the first two and last two fermion fields in each of the terms. We have suppressed colour indices. Due to the antisymmetric colour index and the identity Ψχ c = χΨ c the couplings κ 1ijk and κ 2ijk are necessarily antisymmetric in down type quark flavour.
Other operators in addition to the ones of eq. (3.1) are allowed by the symmetries of the theory e.g. d Ri X 1L X 2L d Rj . For simplicity we take such operators to be sufficiently suppressed compared to those of eq. (3.1) as to be negligible. Even if present these would not qualitatively alter the outcomes below.
The decays of the heavier X α -from now taken to be X 2 -can be CP violating and if they occur out-of-equilibrium will lead to a baryon asymmetry. However, as we will see below, scatterings of the form u i X αL → d j d k also violate CP: we will study the effects of this CP violation on the final asymmetry.
The above Lagrangian contains 28 physical phases if one considers couplings to all possible flavour combinations. As a demonstration and in order to simplify the analysis we will restrict ourselves to considering only couplings to the first generation up quark and second and third generation down quarks, s and b. From now on we will suppress flavour indices on the above couplings. After rephasing we are left with two physical phases. For simplicity we set: while we set the other phase (of κ 3 ) to zero. This choice for the phases leads to the largest CP violation: if the final asymmetry exceeds the observed one the asymmetry can always be reduced by decreasing the κ 1 phase. The other couplings are taken to be real.

Decays
The tree level decay rate Γ(X α → usb) + Γ(X α → usb) is: where we have ignored the masses of the final state particles. X 2 has an additional decay channel with rate Γ(X 2 → X 1 uu): where we have ignored the masses of the final state particles (the full integral expression for massive X 1 can be found in appendix B). CP violation in the decays of X 2 arises from the interference between the tree and one loop diagrams depicted in figure 1. We parameterise the CP violation in the X 2 decay in the following way: The CP violation can be calculated in terms of the underlying parameters of the theory by using the Cutkosky rules to extract the imaginary component of the interference term [34]. The results are shown in figure 2 and further details can be found in appendix B. From dimensional arguments one expects, in the massless limit for X 1 , where the second relation follows for couplings of a similar magnitude O(κ a ) ∼ κ and an order one phase. This is confirmed by our detailed calculation, as can be seen in figure 2. One must also take into account CP violation in the scattering process usb → usb mediated by an X 2 with the real intermediate part of the scattering subtracted [7]. We define this rate as: where the CP conjugate can be found by taking OS → − OS . The unitary condition (2.7) -taking usb as the initial state and summing over all possible final states -enforces the relation: This CP violating rate balances the CP violation in the decays when X 2 is in thermal equilibrium [7] and is included in our Boltzmann equations below. The CP symmetric washout rates W (usb → usb), mediated by off shell X 1 and X 2 , are O(κ 4 ) and are negligible for the parameters we are interested in below: these have been omitted from the Boltzmann equations.

Scatterings
Next we turn our attention to the scatterings. CP violation arises due to the interference between the tree and one-loop diagrams such as those depicted in figure 3. We parameterise the relevant CP violating equilibrium collision terms as: where the CP conjugate interaction rate can be found by making the substitution i → − i . These interactions will contribute to the baryon asymmetry. The unitarity conditions for these interactions give: 1 CP conserving (to the order we are working) collision terms are also present; we denote them in the following way: We calculate the relevant cross sections and find the thermal averaged cross sections numerically by making use of the single integral formula [35]: whereŝ is the centre-of-mass energy squared, p ij is the initial centre-of-mass momentum, K 1 (x) is the modified Bessel function of the second kind of order one and Λ is the effective theory cut-off. The expressions for the cross sections can be found in appendix B. The CP violation in the scatterings at the cross section level i.e. the difference between the cross section and its CP conjugate, (σ − σ), can be calculated in terms of the underlying  parameters using the Cutkosky rules [34]. To find the CP violation in the equilibrium rates one must then take a thermal average: The results for such a calculation are depicted in figure 4. At high temperature T M X2 , one expects on dimensional grounds the CP symmetric thermally averaged cross section to scale as, where κ is the general scale of the relevant couplings. Similary the CP violation is expected to scale as: Unlike the CP asymmetry in the decays, the CP asymmetries in the scatterings are temperature dependent. Note 1 is necessarily suppressed on kinematic grounds at low temperatures (W 1 does not depend on the heaviest particle mass, M X2 , while the CP violation is sensitive to M X2 from the on-shell particles in the loop). The other CP violating parameters take values close to  The simplifying approximation of Maxwell-Boltzmann statistics allows one to factor out the chemical potential from the overall reaction rate density. Let us define: The non-equilibrium reaction rate for the process ij → kl is then given by: The Boltzmann equations can then be conveniently expressed using this notation. The Boltzmann equation for n X1 is given by The Boltzmann equation for n X2 is given by: The Boltzmann equation for the baryon-minus-lepton number, n B−L , is given by: The terms proportional to 3 and D are the source terms coming from scatterings and decays respectively. These lead to the generation of a baryon asymmetry once there is a departure from thermal equilibrium. Note we have used the unitarity conditions to express the source terms in a form in which one sees they are explicitly zero when all r Ψi = 1 (thermal equilibrium) [18]. The other terms are the washout terms: these drive the solutions back to the equilibrium values. We use the standard change of variable to obtain the Boltzmann equations in a radiation dominated universe in terms of temperature using: where g eff counts the effective radiation degrees of freedom and M P l ≈ 1.22 × 10 19 GeV is the Planck mass [30]. The total entropy is conserved in the absence of first order phase transitions in a radiation dominated universe and hence it is useful to normalise number densities to the entropy density Y ψ ≡ n ψ /s. The entropy density is, where h eff counts the entropic degrees of freedom [30]. The resulting differential equations are of the form, where C(ψ) is the collision term for ψ, e.g. for ψ = X 1 it corresponds to the right hand side of eq. (4.2) excluding the expansion term proportional to H.

Chemical potentials
The chemical potentials of the quarks and leptons depend on the baryon and lepton asymmetries. Here we are interested in determining the chemical potentials of the right chiral up and down type quarks. We consider two regimes: above and below the electroweak phase transition. We approximate the transition to occur instantaneously at a temperature T EW at which point we take the SU(2) sphalerons to have switched off and the SM quarks and leptons to have gained masses through the Higgs mechanism.
Field Relativistic degrees of freedom chemical potential Table 1. SM fermionic fields and relativistic degrees of freedom for temperatures below the electroweak phase transition. We set N ν = 3 and use a simple step function in determining N u , N d and N e .

Above T EW
In determining the chemical potentials in terms of n B−L in this regime one takes into account: (i) the SU(3) sphalerons, (ii) the SU(2) sphalerons, (iii) the SM Yukawa interactions, (iv) conservation of the weak hypercharge. The analysis is the same as for the usual leptogenesis [36]. The chemical potential of the right chiral up type quark is given by: The right chiral down type quark chemical potential is given by: The baryon asymmetry is related to n B−L : The chemical potentials of the up and down type quarks are therefore determined in terms of n B−L in this temperature regime. The final n B in this temperature regime is used as an initial condition for the Boltzmann equations at lower temperatures.

Below T EW
Below the electroweak phase transition the electroweak sphalerons have switched off and the weak hypercharge is no longer conserved. However, lepton number is now conserved and the overall electric charge must still vanish. The fermionic fields present are listed in table 1. The baryon number asymmetry is approximately given by [36]: where B = 2N u µ u + 2N d µ d and N u (N d ) is the number of relativistic up (down) type quark flavours. The lepton number asymmetry is approximately given by [36]: where L = 2N e µ e + 3µ ν and N e is the number of relativistic charged lepton generations. The chemical potentials are determined using the following constraints: • Net electric charge vanishes. Once the temperature drops below the mass of a field its number density is very rapidly suppressed, if the field is in chemical equilibrium, as is the case for the SM fermions. We approximate this by only counting fields with mass below the temperature in contributing to the net charge of the plasma. This yields a constraint: (4.13) • Rapid interactions involving W boson exchange give the following relation: (4.14) • Y L ≡ n L /s is constant after T EW . This gives the constraint: where L is computed using the lepton asymmetry Y L at T EW .
Combining these three constraints gives the up quark chemical potential: The down type chemical potential can neatly be expressed as: In this temperature regime the chemical potentials of the quarks can therefore be determined in terms of n B and the net lepton number n L at the sphaleron switch-off temperature.

Numerical solutions
We solve the Boltzmann equations numerically using Mathematica [37]. This involves evolving the three coupled ordinary differential equations from a suitably high temperature to a temperature at which the B violation is no longer effective. The initial temperature is chosen to be high enough so the 2 ↔ 2 interaction rates are well above the expansion rate. The initial conditions correspond to equilibrium distributions for X 1 and X 2 and zero n B−L . 4 At the electroweak phase transition the sphalerons switch off and the lepton number remains constant. We then proceed to track n B rather than n B−L using the value of n B obtained at T EW as an initial condition. Finally a dilution factor is applied to the final Y B if the X 2 and X 1 particles are sufficiently long lived. This dilution factor comes about because the X 2 and X 1 can come to dominate the energy density of the early universe. The universe switches from a radiation to a matter dominated epoch until the X 2 and X 1 decay. The decays of these particles then produces a large amount of entropy -slowing down the rate at which the universe is cooling -and lead to a departure from the assumption of constant entropy in the early universe. The dilution factor can be approximated as [38]: where Y α | T fo is the density to entropy ratio of the decaying particle at freeze out and Γ α is its decay rate. Being more abundant at freezeout and having a longer lifetime X 1 tends to dominate the dilution factor. An example solution to the Boltzmann equations is shown in figure 5. The cosmological history of the universe proceeds in the following way. At high temperature the 2 ↔ 2 interactions are rapid, keeping the n X1 and n X2 close to their equilibrium values. Due to the expansion of the universe the particles are never exactly in equilibrium. As the temperature decreases the interaction rates drop and r X1 , r X2 and Y B continue to increase. The size of the source term will depend not only on the CP violation but also on how far away the temperature is from M X2 and M X1 : in the massless limit r X1 = r X2 = 1 even in the absence of interactions (assuming a common T ). Eventually, at Γ α ∼ H, the decays take over. Excess X 2 and X 1 decay away and the X 2 decays also contribute to the baryon asymmetry. The freezeout temperature is determined by the coupling size; numerically we find the maximum asymmetry for freezeout at T ∼ M α , i.e. before the number densities of X 1 and X 2 become Boltzmann suppressed.
The CP violation in the decays and the scatterings increases ∼ κ. However, as κ increases the freeze out of the 2 ↔ 2 interactions occurs closer and closer to the decay temperature, i.e. at Γ α ∼ H, and the particles are kept closer to thermal equilibrium. The  resulting effect on Y B is shown in figures 6 -8 in which the final asymmetry first increases as the couplings increase but eventually becomes suppressed.
Also shown in figures 6 -8 is the domination of the CP violating scatterings in determing the final Y B for the majority of the parameter space. This is because the CP violation due to the scatterings -which scales as κT 2 -can be relatively high at the freezeout temperature. This balances the relatively small departure from equilibrium just prior to freezeout, when the 2 ↔ 2 interactions generate the majority of the asymmetry, compared with the departure from equilibrium at Γ 2 ∼ H, when the decays contribute to the asymmetry. The CP violation in the decays is typically much smaller than the CP violation in the scatterings at freeze out -at least for couplings small enough for a significant departure from equilibrium to take place. As M X1 is decreased the CP violation in decays becomes more important. The reasons are twofold: smaller M X1 results in a greater | D | and in larger washout effects after X 2 freezes out for kinematic reasons. Note that the scaling of the CP violation with temperature is different in leptogenesis. The dimensionless couplings in leptogenesis mean the CP violation in scatterings will not grow as strongly with temperature, as it does here, but remain mostly constant [24][25][26][27][28]. This explains why the CP violating scatterings can play a crucial role in the neutron portal but only a negligible role in resonant leptogenesis in determining the final Y B [24].
The precision of these calculations can be improved by taking into account departures from kinetic equilibrium, quantum statistics and thermal masses [39][40][41][42][43]. The corrections are expected to be at most O(1), as in the case for leptogenesis. We leave the inclusion of such effects to further work.

Constraints
Constraints on the neutron portal have previously been discussed in ref. [29]. The decay of relic particles after t ∼ 1 s can disrupt big bang nucleosynthesis (BBN) [44][45][46]. Considering the X 1 lifetime one finds a constraint [29]: Similarly from the X 2 lifetime one requires [29]: where we have ignored the final state M X1 mass. The κ 3 bound becomes more stringent as M X1 is increased. Operators of the form κX αL d R u c R d R or κX αL d R Q c L Q L , where X α is Majorana, are constrained by limits on neutron-antineutron oscillations (see figure 9). The oscillation period is estimated as [29], where we have approximated the nuclear matrix element with Λ QCD . Comparison to the experimental limit from bound neutrons τ n−n ≥ 2.4 × 10 8 s [47] (or from free neutrons τ n−n ≥ 8.6×10 7 s [48] -which has smaller theoretical uncertainty) shows broad compatibility in the parameter range of interest. The operators we considered, however, couple only offdiagonally in down quark flavour. The oscillation period is therefore further suppressed and certainly does not pose any problems for the parameter choices we have been interested in above. Similar conclusions hold for the loop induced meson mixing operators [29]. Flavour off-diagonal operators such as κX αL u R s c R b R + κX αL u R s c R d R will lead to meson decays such as B + → π + K 0 (see figure 5). The contribution to the branching ratio is estimated as: which is far below the experimental observation Br(B + → π + K 0 ) = (2.3 ± 0.07) × 10 −5 [49][50][51]. The most important constraints on this scenario therefore come from BBN and neutron-antineutron oscillations.

UV completions
One or more complex scalars transforming as S ρ ∼ (3, 1, 4/3) under the SM gauge group may be added in order to UV complete this model. The high energy interaction Lagrangian takes the form: where y ραi and λ ρij are dimensionless couplings and the latter is antisymmetric in quark flavour. For couplings y ραi ∼ λ ρij ∼ O(1) the effective operators of eq. (3.1) are generated with κ ∼ 1/M 2 Sρ , where M Sρ is the mass of the scalar S ρ , once the S ρ are integrated out.
If the couplings are only to the u, s and b quarks two copies of S ρ must be added in order to obtain complex phases. However, in a realistic model -allowing couplings to all possible quarks -one finds two complex phases with only one S ρ . The corresponding 28 complex phases of the low energy effective field theory are therefore not independent. This can lead to cancellations of phases in some of the interference terms and may reduce the overall CP violation. We leave investigation of such effects to further work.
If we do not impose a global lepton number symmetry the X α will mix with the SM neutrino degrees-of-freedom through Yukawa couplings of the form ab l La Φ * b X c αL , where ab is the Levi-Civita symbol, Φ is the SM Higgs doublet and l L the SM lepton doublet. If the X α are to be responsible for the observed neutrino masses, the mixing between the mostly-active light neutrinos and mostly-sterile heavy states must be of the order θ ∼ m ν /M Xα , where m ν is the light neutrino mass. Consequently S ρ mediated proton decay p → K + ν will occur with a partial lifetime estimated as: (6. 2) The partial proton lifetime will conflict with the experimental limit τ p→K + ν 2.3 × 10 33 years [52] for the parameter range of interest if the X α are responsible for the neutrino mass through the usual seesaw mechanism.
Furthermore Yukawa couplings of the form ab l La Φ * b X c αL would lead to leptogenesis style interactions and subsequent complications in the analysis of the cosmological history. The non-zero neutrino masses can be explained, even with the global lepton number symmetry we have imposed, by introducing right chiral singlets N R with Yukawa couplings ab l La Φ * b N R : these give the SM neutrinos Dirac masses once the SM Higgs gains a vacuum expectation value. Majorana masses for the N R are, of course, forbidden.
In this section we have briefly discussed UV completions for the the neutron portal effective operators of section 3. Other neutron portal type operators are possible and we note that similar conclusions for their UV completions, in particular regarding proton decay if no global lepton number conservation is imposed, also hold.

Conclusion
We have studied the effects of CP violation in 2 ↔ 2 interactions for baryogenesis. As a case study we have taken a neutron portal model in which two Majorana fermions X 1 and X 2 are coupled to the neutron operator. In order to find the final baryon asymmetry, we calculated all the relevant decay rates and 2 ↔ 2 scattering rates and solved the Boltzmann evolution equations numerically. From dimensional grounds the CP violation in the decays scales as κM 2 X2 where κ is the order-of-magnitude of the relevant couplings. The CP violation in the scatterings scales as κM 2 X2 at low temperatures but grows as κT 2 for T M X . Consequently CP violating 2 ↔ 2 scatterings can play a crucial role at high temperature even when the departure-from-equilibrium is relatively small. Indeed for many areas of the parameter space the CP violating scatterings play a dominant role in determining the final baryon asymmetry. This is to be contrasted with leptogenesis -in which the CP violating scatterings do not scale as strongly with temperature -and have only a negligible effect on the final asymmetry [24][25][26][27][28]. We also discussed constraints on the model from experiments and BBN and discussed possible UV completions. The techniques learned can be applied to other baryogenesis scenarios in order to take into account the -possibly dominantcorrections due to CP violating scatterings.
Here we discuss the unitarity constraint in the case of couplings to all possible quark flavours and show that at least two X α are required to obtain a non-zero source term from the CP violating scatterings. 5 We denote the equilibrium reaction rates as: where the CP conjugate can be found by taking αijk → − αijk , ζ αiβj → −ζ αiβj or a ijkl → −a ijkl . Note we have also included quark flavour changing interactions in order to make the argument more general. The unitarity constraint then yields: Now consider the case of only one Majorana X α = X 1 . (The argument may be easily modified for Dirac X 1 and the conclusions remain unchanged.) Taking r d = r s = r b and r u = r c = r t the Boltzmann equation for the baryon asymmetry (excluding the decay terms) is: Now consider only the source term: A departure from equilibrium due to the expansion of the universe gives r X1 = 1 and r u = r d = r u = r d = 1 and the source term may be written: where the first equality follows from eq. (A.4) and the second as ζ αiαj = −ζ αjαi . Hence the generation of an asymmetry from scatterings with only a single X α is not possible.

B.1 Cross sections
Here we collect the results of our calculations for the cross sections used in the above analysis. The sums over initial and final colours have been performed. The cross sections have been calculated in the centre-of-mass frame. E a and E b (E c and E d ) denote initial (final) state energies of the particles in the order listed. The initial momentum is denoted p i and the final momentum p f . The centre-of-mass energy is √ŝ .

B.2 Cross sections -CP violation
CP violation for u R + X 1 → u R + X 2 : This CP violation arises from the interference of two tree level diagrams (corresponding to different continuous fermion lines) with a loop level diagram. where: For this decay we take into account the mass of X 1 . In integral form the width is given by: where µ 12 ≡ (M X1 /M X2 ) 2 and the limits of integration are: This CP violation arises from the interference of a tree level diagram with two different loop level diagrams. The tree level amplitude for X 2 → usb with momenta p i , p a , p b , p c respectively is given by: Here L (R) denotes the left (right) projection operator. The first loop diagram has amplitude: Tr R( / p t + / k)/ k The total CP violation can be expressed as: where γ 1 (γ 2 ) comes from the interference between M * T and M L1 (M L2 ).