Leptogenesis and muon $\boldsymbol{(g-2)}$ in a scotogenic model

We present a detailed study of a scotogenic model accommodating dark matter, neutrino masses and the anomalous magnetic moment of the muon while being consistent with the existing constraints on flavour violating decays of the leptons. Moreover, this model offers the possibility to explain the baryon asymmetry of the Universe via leptogenesis. We determine the viable regions of the model's parameter space in view of dark matter and flavour constraints using a Markov Chain Monte Carlo setup combined with a particular procedure to accommodate neutrino masses and the anomalous magnetic moment of the muon at the same time. We also discuss briefly the resulting collider phenomenology.


Introduction
The Standard Model (SM) gives an accurate description for most of the data up to the TeV scale. Despite its successes, it should nevertheless be considered as an effective theory, which has to be embedded in a more fundamental framework. One reason is the flavour hierarchies in the fermion sector for which we do not know the underlying principle governing its structures. Moreover, there are several experimental observations which require an extension of the SM. This includes neutrinos, which are massless in the SM, but need to be massive in view of neutrino oscillations experiments [1]. Strong arguments from cosmology underline the call for new physics beyond the Standard Model (BSM), such as the presence of dark matter (DM) [2], as well as the baryon asymmetry observed in the Universe [2]. The SM is also challenged by precision measurements of the anomalous magnetic moment of the muon [3,4].
Most of these deviations are located in the lepton sector. Moreover, despite the fact that it ultimately concerns the hadronic sector, the baryon asymmetry can be explained through leptogenesis [5][6][7], a mechanism stemming from the leptonic sector and translating the generated lepton asymmetry to the hadronic sector through the sphaleron processes. Finally, it is worth noting that generating neutrino masses generally leads to the opening of lepton flavour violating processes, involving, for example, transitions from electronic to muonic states, which are strongly constrained by very precise experimental data. We take this as a motivation to consider models featuring new contributions to the lepton sector, while the hadronic sector may be less relevant to explain the above shortcomings of the SM.
One potential class of such frameworks are the so-called scotogenic models, originally aiming at the simultaneous explanation of neutrino masses and cold dark matter. The two are linked in the sense that neutrino masses are generated radiatively through particles and couplings from the dark sector. After the first works on minimal scotogenic realisation [8][9][10][11][12], more complex models have emerged in recent years, studied mainly at the level of dark matter phenomenology and lepton flavour violating observables, see for example [9,[13][14][15][16][17]. A general classification of viable scotogenic frameworks can be found in Ref. [18]. Recently, two of us have studied a particular framework, the so-called 'T1-2-A' model, where the SM is extended by a scalar doublet, a scalar singlet, a fermionic Dirac doublet, and a fermionic singlet [19,20]. This setup features a very predictive dark matter phenomenology, especially for fermionic dark matter [20] and, in principle, it can explain the recent measurements of the anomalous magnetic moment of the muon. However, the corresponding region in parameter space is excluded by the the constraints on flavour violating decays of the leptons. Moreover, the 'T1-2-A' setup fails to accommodate leptogenesis as an explanation for the observed baryon asymmetry.
In the present work, we extend the 'T1-2-A' setup by adding an additional fermionic singlet. The additional degrees of freedom allow for the successful generation of three nonzero neutrino masses, while the couplings can be chosen such that the deviation related to the anomalous magnetic moment of the muon can be accommodated while being consistent with the bounds on flavour violating lepton decays. Moreover, this set-up allows for an explanation of the observed baryon asymmetry via leptogenesis, as we will demonstrate below. Our paper is organised as follows: in Sec. 2, we start by introducing the scotogenic model under consideration. Sec. 3 is then devoted to a discussion of the anomalous magnetic moment within our model and the required coupling hierarchies. In Sec. 4, we discuss the applied constraints and the observables of our interest. The results from our Markov Chain Monte Carlo (MCMC) analysis are presented in Sec. 5, where we analyse the parameter space, charged lepton flavour violating decays, dark matter observables and discuss colliderrelated aspects. In Sec. 6, we present our findings concerning leptogenesis as a means to generate the baryon asymmetry. Conclusions are drawn in Sec. 7.

Model
We consider a scotogenic framework extending the SM by two Weyl fermion SU (2) L doublets, Ψ 1 and Ψ 2 , two Majorana fermion singlets, F 1 and F 2 , a scalar SU (2) L doublet, η, and a real scalar singlet, S. In addition, we assume a Z 2 -symmetry under which the SM fields are even and the additional ones are odd. This ensures neutrino mass generation at the one-loop level together with the existence of a stable dark matter candidate. We note for completeness that the additional fields are singlets with respect to SU (3) C .
The new field content including their respective representations under SU (2) L ×U (1) Y is summarised in Tab. 1. In the following subsections, we briefly summarise the different sectors, present the corresponding Lagrangian, and set the notation.

The scalar sector
The scalar sector of the model consists of the SM Higgs doublet H, an additional real singlet S, and a SU (2) L doublet η. Their charges are given in Tab. 1. Upon electroweak symmetry breaking (EWSB), which involves the Higgs doublet only, the doublets can be expanded into components according to Here, h 0 is the SM Higgs boson, G 0 and G + are the would-be Goldstone bosons, and v = √ 2 H ≈ 246 GeV denotes the vacuum expectation value (VEV). Moreover, η 0 and A 0 are CP -even and CP -odd neutral scalars, and η + is a charged scalar. Neither S nor η may obtain a VEV due to the assumed Z 2 -symmetry.
The scalar potential of the model is given by (2. 2) The first two terms are the SM part related to the Higgs doublet H. We assume here for simplicity that λ η and α are real. After EWSB, the usual minimisation relation in the Higgs sector, allows to eliminate the mass parameter M 2 H in favour of the Higgs self-coupling λ H . Imposing m h 0 ≈ 125 GeV leads to a tree-level value of λ H ≈ 0.13.
The mass matrix of the neutral scalars in the basis {S, η 0 , A 0 } reads as after EWSB. Here, we have defined λ L,A = λ η + λ η ± λ η . We order the mass eigenstates as follows The corresponding squared masses at tree-level read as where m φ 0 1 < m φ 0 2 . Finally, the tree-level mass of the charged scalars is given by

The fermion sector
The Lagrangian for the additional fermions presented in Tab. 1 reads with i, j = 1, 2 and k = 1, 2, 3. L k and e c k denote the left-handed and right-handed leptons, respectively. Moreover, we have introduced the notationφ = iσ 2 φ * for φ = H, η. Without loss of generality we work in a basis where M F 12 = 0. Moreover, we impose |M 1 | ≤ |M 2 |, where we have simplified the notation by setting M i = M F ii for i = 1, 2. Finally, we adopt the phase-convention Ψ 1 = (Ψ 0 1 , Ψ − 1 ) and Ψ 2 = (Ψ + 2 , −Ψ 0 2 ) for the SU (2) L doublets. After EWSB, we have a charged heavy Dirac state Ψ + with mass M Ψ and four neutral Majorana fermions. Their mass matrix is given in the basis This matrix is diagonalised by a unitary matrix U χ according to

Neutrino masses
The main difference with respect to the T1-2-A model discussed in Refs. [19,20] is the extra copy of the singlet fermion. Although the mechanism is very similar, in this case, due to the extra degree of freedom, the neutrino mass matrix has rank three instead of two, and consequently, all three active neutrinos will acquire a non-zero mass. After EWSB, rotating to the mass eigenbasis, a Majorana mass term is generated at the one-loop level via the diagram where the neutrino mass matrix can be expressed as This is a well-known structure common to most of the scotogenic models and similar to the type-I seesaw, where the matrix G contains the couplings defined in Eq. (2.8) ordered as and M L is a 3 × 3 symmetric matrix which encodes the information of the loop function, and the mixing in the neutral scalar and fermion sectors, defined in Eqs. (2.5) and (2.10), respectively. For completeness, we explicitly write the expressions for the components of M L , where k = 1, 2, 3, 4 and n = 1, 2, 3. Moreover, the loop integrals are encompassed in the functions (2.20) We make use of the Casas-Ibarra parametrisation [21,22] to express the couplings in Eq. (2.13) in terms of neutrino oscillation data [1,23] according to where D L is the diagonal matrix defined by and D ν is the diagonal matrix containing the neutrino mass eigenvalues. Finally, U PMNS is the usual unitary matrix relating neutrino flavours to their mass eigenstates, assuming that the charged leptons are already in their mass eigenbasis. Moreover, as even a precise knowledge of all the parameters and observables in M L and M ν does not univocally define G, the extra degrees of freedom are encoded in the orthogonal 3 × 3 matrix R. This matrix can be parameterised as depending on three complex angles θ i with s i = sin θ i and c i = 1 − s 2 i . Note the importance of these degrees of freedom, since they modify the flavour structure of the Yukawa matrix. The latter is of great relevance when considering charged lepton flavour violation and the anomalous magnetic moment, as we will show in the next sections.

The anomalous magnetic moment of the muon
As mentioned in the introduction, a deviation persists between the SM prediction and the experimental value of the anomalous magnetic moment of the muon, defined as a µ = (g − 2) µ /2. The discrepancy amounts to a significance of 4.2σ, and leads to the following range for the new physics contribution to a µ [4], 1 In general, every scotogenic-like model will contribute to the anomalous magnetic moment of leptons at one-loop level. These contributions can be encoded in the effective electromagnetic (EM) dipole moment operator c ij R¯ i σ µν P R j F µν , coming from the operator O eB ≡ (Lσ µν e R )HB µν before EWSB [27]. The diagonal part of the Wilson coefficient c R is related to (g − 2) and the electric dipole moment (EDM), and the off-diagonal part is associated with charged lepton flavour violating (cLFV) processes [28]. For more details see Appendix A.
The contribution to (g − 2) µ is generally suppressed by the muon mass. Moreover, as the EM dipole operator connects the left-and right-handed parts of the leptons, while neutrino mass models contain usually only couplings of BSM fields to the left-handed components, such an operator will always be chirally suppressed. Consequently, new physics explanations of (g − 2) µ are pushed towards low mass scales and large non-perturbative couplings. A possible way out is to add new fields outside the neutrino mass mechanism, that couple to µ R , in order to enhance the contribution to (g − 2) µ and be able to fit the anomaly within a phenomenologically reasonable parameter space [29]. Note that this situation is realised in the T1-2-A model [18], and consequently its extension under consideration here. In both models, there is a coupling g R of the lepton singlets to η and Ψ 1 , see Eq. (2.8). The latter two also participate in the generation of the neutrino mass matrix. Note that in this way no extra BSM field is needed on top of those involved in the neutrino mass mechanism to have a chirally enhanced contribution to (g − 2). The new leading contributions to the anomalous magnetic moment are shown in Fig. 1.
We note for completeness, that in the original T1-2-A framework the coupling matrix G, see Eq. (2.13), is a 2 × 3 matrix, where the relative sizes between the various entries are fixed by the neutrino mixing angles up to one complex angle. An explanation of the muon (g −2) in this model implies large couplings which in turn lead to too large flavour violating decays of the leptons. In our extension of this model, we have more freedom allowing us to circumvent this problem, as we will show in Sec. 5.2.
Both diagrams depicted in Fig. 1 also generate a sizeable contribution to strongly constrained LFV processes in the charged sector, in particular µ → eγ with an upper limit to its branching ratio of 4.2 × 10 −13 from the MEG collaboration [30]. Although they seem unavoidable, given that the off-diagonal part of the Yukawa matrix G is connected to the neutrino mixing (see Eq. (2.21)) there are several strategies to get a sizeable contribution to (g −2) µ while keeping charged LFV under control. For example, cancellations can be found Figure 1: Dominant one-loop contributions to (g − 2) and charged LFV processes before EWSB. Arrows indicate the flow of quantum numbers. Couplings are given for clarity, see their explicit definitions in Sec. 2. A photon should be attached to the respective charged components.
among the several independent contributions to the EM dipole operator. However, such a scenario is not very appealing, as one should reproduce a difference of more than five orders of magnitude between the diagonal and off-diagonal components of the Wilson coefficient c R [28]. Another possibility is to assume certain flavour structures for the Yukawa couplings, which suppress the off-diagonal components in favour of the diagonal [31]. Following the latter approach, we focus for simplicity on a region of the parameter space where the first diagram in Fig. 1 dominates over the second, as the flavour structure of the diagram is simpler having just two three-component Yukawa vectors involved. We extend the usual Casas-Ibarra parameterisation by the following elements so that these constraints can be easily fulfilled. To do so, we consider y 1,2 to be small and push the trilinear coupling α to larger values, i.e. we suppress the mixing in the neutral fermion sector, while enhancing the one in the neutral scalar sector. Note that, while g R is mainly free, g F and g Ψ are constrained by the fit to neutrino oscillation data, see Eq. (2.21). This means that changing y 1,2 and α not only directly modifies the dominant contributions depicted in Fig. 1, but also indirectly suppresses g F and enhances g Ψ through the neutrino fit. We are looking for a Yukawa matrix G featuring a coupling hierarchy as shown in Fig. 2. Making use of the freedom on the components of g R as well as on the remaining degrees of freedom in g ψ , stemming from the rotation matrix R appearing in Eq. (2.21), we fit the value of a BSM µ while keeping the contributions to the lepton flavour violating decays µ → eγ and τ → µγ under control.
In practice, for each point of our numerical scan, in the region of the parameter space where y 1,2 are small, we use the angles of the matrix R given in (2.23) to suppress the dominant contribution to cLFV processes while enhancing the diagonal contribution associated to (g − 2) µ 2 . Ultimately, we fit the experimental value of the muon (g − 2) within its limits by solving for g 2 R . With this method, we obtain for each point the correct anomalous 2 Actually, solving for two angles is sufficient, such that one angle is left as a free parameter and scanned over for generality. Figure 2: Hierarchy of the Yukawa matrix G of Eq. (2.13) that generates the neutrino masses. This hierarchy is realised for y 1,2 small and solving for the angles defining the rotation matrix R, such that (g − 2) µ is maximised while charged lepton flavour violating decays are kept under control. See text for more details.
magnetic moment, while fulfilling the current limits for charged LFV decays. We note that this is not the most general approach, as we are selecting a specific region of the parameter space. However, given the complexity of the system, we were not able to find a more general approach that could deliver results within a reasonable computing time. Moreover, the parameter space discussed previously is also preferred for low-scale leptogenesis, as we will show later in the paper.

Constraints and observables
In the spirit of the analysis presented in Ref. [20], we use an MCMC scan [32] based on the Metropolis-Hastings algorithm [33,34] to efficiently scrutinise the parameter space of the model in view of the numerous constraints presented above. This technique, especially powerful for high-dimensional spaces, explores the parameter space iteratively, restricted by a set of constraints through the computation of the likelihood. We refer the reader to Ref. [20] for further details about the implementation of the MCMC.
In addition to implicitly satisfying the constraints from neutrino masses and the anomalous magnetic moment of the muon (see above), we explicitly impose constraints coming from various sectors, comprising dark matter, lepton flavour violating processes, and the mass of the Higgs boson. All constraints are listed in Tab. 2 together with their associated experimental limits, as well as applied uncertainties applied in our study. Note that for the Higgs mass m H and the dark matter relic density Ω CDM h 2 , the theory uncertainties 3 [35][36][37][38][39] are larger than the experimental ones, and, consequently, we apply the theory uncertainties. We also ensure that the lightest Z 2 -odd particle, is electrically neutral in order to have a viable, stable DM candidate and to avoid stable charged relics, essentially excluded in the mass range of 1, 10 5 GeV [40][41][42][43].
In total, our MCMC scan runs over 20 free parameters: Eight couplings in the scalar potential, six Yukawa couplings, five masses, the lightest neutrino mass, and the unconstrained angle of the rotation matrix R, which is assumed to be real. The ranges of the scalar and fermion mass parameters are chosen such, that they could be in principle in the reach of high luminosity LHC. The exception is those for the singlet fermions for which we Table 2: Constraints considered in the MCMC analysis: Higgs mass and cLFV observables [40] and DM relic density [2]. The limits from XENON1T [44] to the direct detection crosssection are also taken into account. Note that the errors given for m H and Ω CDM h 2 are not the experimental uncertainties but estimates of the theoretical ones, see text for details. We implement 1σ intervals using a Gaussian function and 90% C.L. for the limits via a single-sided Gaussian, allowing for a 10% uncertainty.
allow a larger range. The reason is that this model can also explain the baryon asymmetry of the Universe via the leptogenesis mechanism as discussed in Sec. 6. The sign of the quartic couplings λ H , λ 4S and λ 4η is fixed from the requirement that the scalar potential is bounded from below. We vary all parameters on a logarithmic scale and assign possible signs on a random basis.
The scan is performed over the parameter ranges specified in Tab. 2 with 75 chains of 200 points each. The first 35 points of each chain have been deleted in order to keep only the points for which the chains were already well initialised, i.e. presenting a phenomenologically viable likelihood value. For the scan we implemented the model in SARAH [45] and generate code for SPheno [46], FlavorKit [47] and micrOMEGAS [48]. The former two compute the mass spectrum and low energy observables, while the latter evaluates the DM relic density and the direct detection (DD) cross-sections.
Assuming a Gaussian likelihood of uncorrelated observables, the likelihood associated with a given parameter point n is computed as where the product runs over the imposed constraints and individual likelihood value L n i associated to each constraint. In the case of a two-sided limit, i.e. for the Higgs mass m H and the DM relic density Ω CDM h 2 , the likelihood is computed following Parameter Interval where O n i is the calculated value of the considered observable for the parameter point n, O exp i is the associated experimental value given in Tab. 2, and σ i is the associated uncertainty. The likelihood computation for upper limits is implemented as a step function, which is smeared as a single-sided Gaussian with a width of 10% of the value corresponding to the experimental upper limit. In this case, we have L n i = 1 if the predicted value O n i is below the upper limit O exp i . In the opposite case, L n i is computed according to Eq.

Results
In this section, we present the main outcome of our MCMC analysis. We shall first show the resulting parameter space for the couplings and then discuss certain observables of interest. In addition, we will discuss possibilities to test part of the available parameter space at the LHC.

Couplings
We are interested mainly in the Yukawa couplings that connect the SM particles with the new fields, i.e. g F 1 , g F 2 , g Ψ , and g R , which are all three-component vectors, see Eq. (2.13). These are relevant for neutrino masses, the anomalous magnetic moment and flavour violating decays of the leptons. Fig. 3 shows the correlations among the different components for each coupling vector, while in Fig. 4 we show the correlation of selected components with the trilinear coupling α.
We clearly see that all components of g F 1,2 behave in a similar way, with an approximate upper limit of |g i F 1,2 | 10 −3 for i = 1, 2, 3. As already explained in Sec. 3, this upper limit is due to our fit of neutrino masses together with the anomalous magnetic moment (g − 2) µ and the constraints coming from µ → e transitions. At the same time, the overall scaling behaviour of all the components of g F 1,2 is caused by the trilinear coupling α, as shown in the left plot of Fig. 4 for g F 1 . An analogous behaviour is found for g F 2 (not shown). Larger values of α imply a larger scalar mixing, which then suppresses the scale of g F 1,2 via the neutrino mass fit.
For g Ψ , as already described in Sec. 3 and depicted in Fig. 2, a specific hierarchy among its components is realised to fit the muon anomalous magnetic moment and be below the  : Correlation of selected Yukawa couplings with the trilinear coupling α. The couplings g Ψ and g F 1 are connected to the trilinear couplings α through the fit of the neutrino masses, while the connection of g 2 Ψ and g 2 R with α stems from the fit of the anomalous magnetic moment (g − 2) µ . Figure 5: Results for the relevant cLFV decays with the current limits from MEG collaboration [30] and Belle [52,53] (full lines) and expected sensitivities (dashed lines) from MEGII [49], Mu3e [50] and Belle II [54]. The other decays not shown here lay below the expected future bounds.
limits of charged LFV searches. This hierarchy is linked to that of g R , as both contribute equally to these processes, see Appendix A. While both g 2 Ψ and g 2 R have to be large to fit (g − 2) µ , g 1 Ψ and g 1 R must remain small to not exceed the current limit of BR(µ → eγ). On the same grounds, g 3 Ψ and g 3 R are similarly constrained by the upper limit on BR(τ → µγ). It is worth noting that the fit of (g − 2) µ links the components of g R and g Ψ with the trilinear coupling α, as can be seen in Fig. 4 (right). As discussed in Sec. 3, the dominant contribution to (g − 2) µ and charged LFV decays comes from the left diagram in Fig. 1, proportional to α. For example, smaller values of α imply larger values of g 2 Ψ and g 2 R in order to fit the anomalous magnetic moment (g − 2) µ , as can be seen in the upper corner on the right panel of Fig. 4.
The perturbativity requirement for both the Yukawa couplings and α sets then a lower and upper limit on the trilinear coupling of roughly 30 GeV α 4 m φ 0 1 . Note, that the upper bound is actually given for α/M φ where M φ is the average masses of the scalars involved in this coupling.

Charged lepton flavour violating decays
Charged lepton flavour violating decays rank among the most stringent constraints for neutrino mass models, as fitting the neutrino mixing angles, in general, requires nondiagonal Yukawa matrices that connect also to the charged leptons and allow for transitions between different lepton flavours. While the limits to the branching ratios of these processes are already remarkable, especially for the limit on the decay µ → eγ from the MEG collaboration [30], there is a renovate interest with new experiments expected to take place in the near future, such as MEGII [49], Mu3e [50], or COMET [51], with an expected improvement on the sensitivity of even four orders of magnitude for certain processes like µ → 3e, or Belle and Belle II for the tau decays [52][53][54].  Figure 6: Histograms of the mass and nature of the dark matter candidate. The separation into fermionic and scalar dark matter candidates clearly exhibits a preference for fermionic dark matter with a mass around 1100 GeV.
Although the charged LFV decays are considered as constraints in our analysis, see Tab. 2, it is worth checking how these processes behave in the present model and exploring to which extent future experiments may restrict the parameter space. In Fig. 5 we show the branching ratios of the most relevant charged LFV decay channels for the muon and the tau, together with their current limits and future expected sensitivities. Note that the muon decays are completely dominated by the dipole contribution, i.e. the diagrams depicted in Fig. 1, while there is an important contribution from box diagrams to the tau decays τ → 3µ. This is due to the large values for g 2 Ψ and g 2 R needed to the fit of (g − 2) µ , which makes the box diagram proportional to g 3 R g 2 * Ψ g 2 * R g 2 Ψ dominate over the dipole contribution with the off-shell photon.
For conciseness, we do not show the results for the charged LFV decays of the tau to electrons. For these, we find that, in the best case, the branching ratio of τ → µe + e − is of the order of 10 −9 , just on the border of the future expected sensitivity. For the remaining processes τ → eγ and τ → 3e, we obtain branching ratios below 10 −17 , not observable in the foreseeable future.

Dark matter observables
Let us recall that the model under consideration includes three possible candidates for cold dark matter (CDM), the lightest Z 2 -odd neutral fermion χ 0 1 , the lighter scalar φ 0 1 , and the pseudo-scalar A 0 , depending on the mass hierarchies in a given parameter configuration. The dark matter relic density is taken as a constraint in our MCMC analysis, with a theory uncertainty of 10%, such that Ω CDM h 2 = 0.120 ± 0.012 [2].
Starting with the overall situation, Fig. 6 shows the obtained distribution for the DM mass, separating fermionic (χ 0 1 ) and scalar (φ 0 1 ) DM. The shown results exhibit similar behaviour to that found in Ref. [20] for the simpler T1-2-A scotogenic model. Fermionic DM dominates the model parameter space, with a preferred mass of around 1100 GeV. Scalar DM accounts for about 38% of the viable parameter points, with a preferred masses of roughly 600 to 1000 GeV. As in the T1-2-A model [20], fermionic DM is essentially doublet-dominated. This can be seen from Fig. 7 showing the split into doublet and singlet-dominated DM candidates together with the singlet content of the fermionic DM. This feature can be traced to necessary co-annihilations which occur naturally in the doublet case between χ 0 1 and χ ± and χ 0 2 due to the very small mass splitting between these states. The (co-)annihilation processes are dominated by gauge interactions similar to the case of pure higgsinos in supersymmetric models. As already mentioned, sizeable Yukawa couplings to the muons are necessary to explain the potential deviation of its anomalous magnetic moment. This gives additional annihilation channels into muons via scalars, stemming from the doublet η, in the t-channel.
In the case of a fermionic singlet-dominated state χ 0 1 , it is much harder to satisfy the relic density requirement from Planck data [10]. While singlet fermions F i are produced thermally, they can only annihilate via the Yukawa g F i or through the mixing with Ψ 1,2 , both small because of the charged LFV constraints and our fit of (g − 2) µ .
In the case of scalar DM, the doublet-like states also dominate the phenomenologically viable parameter regions with preferred masses around 700 GeV. For the same reason as in the case of fermionic DM discussed above, doublet-like DM is preferred, as it allows meeting the relic density constraint through efficient co-annihilations between the different doublet states mediated through gauge interactions. Note that the associated mass peak is wider here, as the mass splitting in the doublet can be larger in the scalar case as compared to the fermionic one.
Finally, we note that we do not find any pseudoscalar DM in this model in contrast  to Ref. [20]. The reason can be easily understood by inspecting the mass matrix given in Eq. (2.4). The mass splitting between the scalar doublet component and the pseudoscalar is given by 1 2 λ η v 2 and the mixing between the doublet scalar and the singlet is given by |αv|. For a pseudoscalar DM candidate one needs the doublet to be lighter than the singlet. As can be seen in Fig. 9 the mixing between the scalar components is always larger than the mass splitting between scalar and pseudoscalar implying that the scalar will be lighter than the pseudoscalar. Note that this feature is less pronounced in the T1-2-A model studied in Ref. [20], as in their study the trilinear coupling α is rather restricted and, in Figure 10: Spin-independent direct detection cross-section versus the mass of the DM in the scalar case, differentiating between singlet (blue) and double-like (orange). The current limit from XENON1T [44] and the future limits from XENONnT [55] and DARWIN [56] are given, as well as the corresponding line for the neutrino floor [57]. The fermionic DM case is not shown as the DD cross-section lays below the neutrino floor, around 10 −60 cm 2 . See text for details. addition, the constraint from (g − 2) µ has not been taken into account.
In Fig. 10 we show the results for the spin-independent direct detection cross-section for the scalar DM case. As already said in Sec. 4, the XENON1T [44] limit is taken as a constraint, such that points not satisfying the current limits are excluded. Most of the remaining viable points can be tested by future experiments like DARWIN [56]. On the other hand, we do not find any constraints in the case of fermionic dark matter. The reason is that our scan requires the modulus of the relevant Yukawa couplings, |y ij |, to be smaller than 10 −4 , suppressing the dominant contribution and pushing the cross-sections well below the neutrino floor.
We note that, in both cases, the direct detection cross-section is mainly dominated by Higgs exchange, since we actually have an inelastic dark matter candidate. Inelastic dark matter refers to DM candidates with a mass splitting between the CP -even and CP -odd components of a neutral state. As the Z-boson couples always between the CP -even and CP -odd components, for the part of the parameter space where the mass splitting between these two states is larger than the kinetic energy of the DM, the contribution from Z channel to the DD cross-section is kinematically forbidden. Since the coupling of the DM to the Z-boson has typically gauge strength, if it is active, then it will be excluded by direct detection. We note here that this contribution was added by hand as micrOMEGAS does not include inelastic channels. Nevertheless, this excluded very few points, as in practice given Figure 11: Mass of the charged fermion χ ± as function of the mass of the neutral fermion in the scenarios with a fermionic DM candidate.
the typical DM average relative velocity, the mass splitting needs to be only larger than O(100) keV to kinematically close the Z-channel [58].

Collider aspects
We have seen in Sec. 5.3 that the preferred range for fermionic DM is between 700 GeV and 1.4 TeV with most of the points having a mass between 1 and 1.2 TeV. Moreover, they are essentially always SU (2) L doublets with the same quantum numbers as higgsinos in supersymmetric models. From Ref.
[59] we can thus infer that the cross-section σ(pp → χ + χ 0 ) at the LHC with √ s = 14 TeV varies between 0.43 fb (1 TeV) and 0.14 fb (1.2 TeV) assuming that both states have the same mass. The signatures depend on the mass difference between the charged and the neutral state which is displayed in Fig. 11. In case the DM candidate stems from SU (2) L , one finds only a small difference and the dominant decay mode is χ + → π + χ 0 . We can infer from the corresponding supersymmetric scenarios that the LHC will not be able to discover the corresponding states, see for example [60] and references therein. In case the DM candidate is singlet-like we find mass splittings between 10 and 150 GeV. The main decay modes are via off-shell neutral scalars The interesting point is that the requirement of explaining the (g − 2) of the muon implies that nearly all cases one has a muon in the final state. While this is a potentially interesting final state, one should keep in mind that the case of singlet fermionic DM seems to be rarely realised in this model, see Sec. 5.3.
In scenarios with a scalar DM candidate, the situation looks somewhat more promising. We see from the left of Fig. 12 that in a large portion of the corresponding parameter space the charged fermion has a significantly larger mass than the DM candidate. Note, that the SU (2) L have a similar mass as the charged fermion. Both will decay into SM-leptons and a Z 2 -odd scalar, where the j corresponds mainly to the mass eigenstates which are dominantly SU (2) L fermions. As above, we expect the decays into muons to be dominant. Thus, the signal will dominantly consist of muons in combination with missing transverse energy. We note for completeness, that one has of course also direct production of scalar doublets. However, already with a mass of about 700 GeV the cross-section is about 0.1 fb as can be inferred from the production of left-sleptons [59] which have the same quantum numbers. We see from the left of Fig. 12 that m η ± > ∼ 700 GeV in scenarios with a sizable mass splitting in the scalar sector. Thus, the direct production will hardly contribute to an LHC signal for this model.

Leptogenesis
This model features heavy Majorana fermions, lepton number violation and complex couplings which are all ingredients for leptogenesis. In this section, we investigate to which extent one could also explain the observed baryon asymmetry of the Universe via the leptogenesis mechanism 4 , in the region of parameter space discussed in the previous section. We present here the main results and collect in Appendix B further details.
We are in a region of parameter space where the couplings y ij , which determine the mixing between the SU (2) L doublet and singlet fermions in Eq. (2.9), are small. Thus, in practice, only the singlet-like fermions will contribute to a possible lepton asymmetry.
These states decay dominantly according to F i → ηL , η †L and F i → HΨ , H †Ψ . (6.1) The decays will occur at a mass scale which is significantly above the scale of electroweak symmetry breaking and, thus, it is more convenient to work in the gauge basis. At tree-level the former is governed by the couplings g i F and the latter by the couplings y ij , where we have neglected the masses of the SM leptons and Higgs boson. The asymmetry is generated at the one-loop level [5] by the diagrams displayed in Fig. 13. Similar to the type-I seesaw model [6,61], there are the typical wave-function and vertex diagrams depicted in the upper row involving the other singlet fermion in the loop. Moreover, there are additional possible vertex diagrams, which are depicted in the lower row. These involve additional couplings like g k R and α which turn out to be important. The CP asymmetry parameters i are given by Explicit expressions for the various contributions to i in this model may be found in Appendix B.
The MCMC yields that in general max(|g k F i |) max(|y ij |) in most of the parameter space and in the remaining part, they are of equal size. Thus, we will base our estimations in the text below on g k F i , but we stress that all parameters were properly taken into account in the numerics.
An important question is to which extent a generated asymmetry gets washed out in the thermal history of the Universe. To answer this, one defines the decay parameters with H being the Hubble parameter. The weak washout regime is realised for K i 1, the strong washout regime for K i 3, and an intermediate regime in between these values [6,7]. In the parameter space discussed in Sec. 5, we are always in the strong washout regime, as can be seen from the second part of Eq. (6.5). As an example we display K 1 in Fig. 14.
The very large values of the decay parameter allow us to neglect washout through scattering processes, as the inverse decays are the dominant sources of the washout. Thus, we may treat leptogenesis as competition between decays and inverse decays [6], with the  B − L asymmetry being generated when the inverse decays freeze out and the surviving number density N F i of the F i decays.
We solve numerically the corresponding Boltzmann equations given in Appendix B with the following initial conditions at high temperatures, and track these number densities down to lower temperatures. As mentioned earlier, we do not assume a large hierarchy between the masses of F 1 and F 2 . Consequently, we take into account both contributions. Note that the inverse decays of the singlet fermions must freeze-out at a temperature before the sphalerons fall out of equilibrium (T ∼ 100 GeV), otherwise, the B − L asymmetry generated in its decay is not converted into baryon asymmetry. As we are in the strong washout regime, we can estimate the freeze-out temperature for the inverse decays as [62], This gives an approximate lower bound on the mass of the lightest decaying singlet M i 2 TeV for which the sphalerons remain active. At lower temperatures, T min{M 1 , M 2 }, we obtain the final B − L asymmetry, N fin.
B−L , which is converted to the baryon-to-photon ratio, η B , via the sphaleron process [63,64] using C sphal. = 8/23 is the sphaleron conversion factor [65], g 0 * = 43/11 the present value of the number of relativistic degrees of freedom (DOF) and g * the relativistic DOF of the full model at high temperatures. Figure 14 depicts the final baryon-to-photon ratio obtained from solving the Boltzmann equations, using the sets of parameters mentioned before, against the mass of the singlet fermion driving leptogenesis. We observe that in contrast to the typical case of strong washout in the type-I seesaw model and the 'classic' scotogenic model [66], the final value of η B has a tendency to decrease with increasing M 1 . Besides, we also find that the CP asymmetries generated in the decays of lighter F 1 are much larger. We note large contributions to i come from the loop diagrams in the lower row of Fig. 13.
In the minimal scotogenic model, it is possible to express the decay parameter (6.5) as a function of the lightest neutrino mass and the λ η parameter [66]. In our model, the link is not so direct due to the additional couplings and particle states present. Moreover, the requirement to explain the potential deviation of the (g − 2) of the muon while being consistent with the bounds on the LFV lepton decays requires, for example, larger values of the trilinear coupling α. As can be seen from the left plot in Fig. 14, the decay parameter decreases with increasing |α|. Note that this coupling does not appear directly in the calculation of the decay parameter (6.5), but reduces the value of the coupling g F through the neutrino fit, which in turn decreases the tree-level decay width and, hence, lowers the value of K i .
We also see from the right plot in Fig. 14 that only a few points are able to explain the observed baryon asymmetry. Investigating the details of the required parameter combinations to obtain the correct baryon asymmetry would be highly interesting, but is beyond the scope of this paper and left for a future publication. However, one feature that we have observed is that nearly all points present a fermionic doublet as dark matter candidate. Only one out of 25 points in parameter space contains a singlet scalar dark matter candidate.

Conclusion
We have investigated a scotogenic model with a very rich phenomenology. We have presented a complete analysis of the associated parameter space, taking into account constraints from the Higgs sector, the neutrino sector, lepton flavour violating processes, the muon anomalous magnetic moment and dark matter observables.
Neutrino data governs the couplings of the new particles to the left-handed leptons and the requirement of explaining the observed deviation of the anomalous magnetic moment of the muon (g − 2) requires that couplings to muons are sizeable. This in turn implies that the decays µ → eγ, µ → 3e, τ → µγ and τ → 3µ are in the reach of upcoming experiments in a sizeable part of the parameter space.
We have found that the dark matter relic density constraint leads to a preference of fermionic dark matter candidates, in most cases the neutral component of an SU (2) doublet in the mass range 1 to 1.2 TeV. Scenarios featuring a scalar dark matter candidate can be tested by future direct detection experiments like XENONnT or DARWIN, whereas the corresponding cross-sections for fermionic dark matter are well below the so-called neutrino-floor.
We have also briefly discussed the LHC phenomenology in the relevant parameter space. The requirement of explaining the observed deviation in (g −2) leads to a preference for decays into final states containing muons. In particular, in case of a scalar dark matter candidate, we expect muons plus missing transverse energy as dominating signal at the LHC. This signal is also expected in case of supersymmetric models due to the decays of the so-called smuons. However, in our case final states with other leptons or jets in combination with missing transverse energy will be (strongly) suppressed. In case of fermionic dark matter, the mass differences are so small that the charged fermion will decay into a pion and the neutral fermion. The discovery of this final state will be very challenging at the LHC as the required mass imply a relatively low cross-section and, thus, it is likely that the LHC will not be able to observe these particles even in the high luminosity phase.
Finally, we have seen that the available parameter space gets severely constrained if one requires in addition an explanation of the observed baryon asymmetry of the Universe via leptogenesis. We have found that nearly all viable points in the parameter space feature a fermionic dark matter candidate. A detailed analysis of the features of the relevant part in the parameter space will be presented in a future work.

Acknowledgments
The work of M. Sarazin  A New contributions to the electromagnetic dipole moment operator In this appendix, we collect the additional contributions to the Wilson coefficients of the dipole operator using the notation of Ref. [28]. The anomalous magnetic moment, as well as other observables like the electric dipole moment (EDM) and the charged lepton flavour violating (cLFV) decays, is directly connected to the electromagnetic dipole moment operator, i.e.
As already explained in Sec. 3, the diagonal of the Wilson coefficient c R is linked to the (g − 2) and the EDM by The real part of the Wilson coefficient will contribute to the anomalous magnetic moment (g −2), whereas the imaginary part is a contribution to the electronic part, the EDM. With this, the anomalous magnetic moment is defined as, On the other hand, cLFV decays can be computed directly from Eq. (A.1), for which their branching ratio is given by, valid for m i >> m j and where Γ i is the total decay width of i . In a general framework, we can consider a coupling between new fermions Ψ, new scalars Φ, with the Standard Model lepton i . In this scenario, we have the Lagrangian, The associated Wilson coefficients are given by, where Q is the electric charge of the fermion in the loop. The functions f, g,f ,g are defined as, In order to compute c R , we need to determine the vertices Γ i L and Γ i R . In our model, as described in Sec. 2, there are two contributions depicted in Fig. 1. The Wilson coefficient c R is then the sum of these two contributions.
For the first diagram (left) in Fig. 1, after EWSB the fields in Eq. (A.5) correspond to Ψ ≡ χ + and Φ ≡ φ 0 k , with couplings We note here that for simplicity we assumed a sum over the new fields Ψ and Φ in Eq. (A.6).
Here, a sum over the index k = 1, 2, 3 must be performed when computing c R , taking into account that several scalars participate in the diagram, i.e. one should replace M Φ by m φ 0 k . Also, given our definition of the neutral scalar basis and their mixing as defined in Eq. (2.5), where for conciseness we included the pseudo-scalar A 0 , here the third φ 0 eigenstate should be considered as φ 0 3 ≡ A 0 . For the second diagram (right) in Fig. 1, after EWSB Ψ ≡ χ 0 k and Φ ≡ η − , with couplings Again, a sum over the index k = 1, 2, 3, 4 must be performed in Eq.

B Details on the calculation of the baryon asymmetry
We collect here details of the calculation of the baryon asymmetry of the Universe which has been presented in Sec. 6. We largely follow the convention used in [6].

B.1 Boltzmann equations for leptogenesis
As a starting point, we write down the Boltzmann equations for thermal leptogenesis, which track the number densities per co-moving volume N F i and N B−L of the singlet fermions and the B − L asymmetry. Here, we have defined the dimensionless variables z i = M i /T where T is the temperature. As the variables are related by z 2 = z 1 (M 2 /M 1 ), we solve the above equations in terms of z 1 . Here we take the contributions of both singlet fermions into account as we do not necessarily have a large mass hierarchy between the two states, see e.g. Refs. [7,67]. The quantities K 1,2 (z i ) refer to the modified Bessel functions of the second kind and N eq.
where M Pl = 1.22 × 10 19 GeV is the Planck mass. The quantity g * refers to the effective number of relativistic degrees of freedom. With the additional particle content, we have g * = 122.25 compared to the SM value of g SM * = 106.75 [68]. We solve the set of coupled Boltzmann equations (B.1) and (B.2) numerically with the following initial conditions at T = 10 5 M 2 implying z 1 1 The sphaleron conversion factor enters the calculation of the baryon asymmetry which is given according to Ref. [65] as where n D is the number of scalar SU (2) L doublets in the model.

B.2 Computation of i
Here we collect the expressions for the various diagrams that contribute to the CP asymmetry parameter i : Evaluation of the Dirac traces and the required Passarino Veltman reduction of the loop integrals was carried out using FeynCalc [69][70][71]. The imaginary part of the B 0 function is given by where λ(x, y, z) = x 2 + y 2 + z 2 − 2xy − 2yz − 2zx is the Källén function and Θ(x) is the Heaviside-step function, which enforces the fact that the imaginary part exists only when particles in the loop can go on-shell. The imaginary part of the C 0 function was calculated numerically using PackageX [72]. We indicate all possible sums over the particles in the loop and have summed over the final leptonic states, i.e. we do not consider flavour effects. Finally, we have defined Im (g j F i ) * (y SM jk ) * (g k R ) * y 2i 1 + M 2