1 Introduction

In the search of Physics beyond the Standard Model (SM), disagreement between measurements and theoretical expectations, that is “anomalies”, can play the role of beacons to guide our explorations. One longstanding anomaly concerns the anomalous magnetic moment of the muon \(a_\mu =\frac{g_\mu -2}{2}\). The Muon g-2 experiment at Brookhaven [1] and its successor at Fermilab [2, 3] have produced the following result

$$\begin{aligned} \delta a_\mu ^{\mathrm{Exp}}=a_\mu ^{\mathrm{Exp}}-a_\mu ^{\mathrm{SM}}=(2.5\pm 0.6)\times 10^{-9}. \end{aligned}$$
(1)

where \(a_\mu ^{\mathrm{Exp}}\) is the experimental observation and \(a_\mu ^{\mathrm{SM}}\) the SM theoretical expectation [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24]. Although there are unsettled discrepancies concerning Hadronic Vacuum Polarization (HVP) contributions to \(a_\mu ^{\mathrm{SM}}\) [25,26,27], we interpret \(\delta a_\mu ^{\mathrm{Exp}}\) in Eq. (1) as a signal of New Physics (NP).Footnote 1

Besides the muon, recent results concerning the anomalous magnetic moment of the electron might also be interpreted as NP hints [32]. On the one hand, perturbative calculations of \(a_e=\frac{g_e-2}{2}\), which have reached impressive levels [5, 33,34,35,36], yield \(a_e^{\mathrm{SM}}\) as a series in powers of the fine structure constant \(\alpha \). On the other hand, we have precise measurements of \(a_e^{\mathrm{Exp}}\) such as [37]. In the past, such measurements were indeed used to infer values of \(\alpha \). On the contrary, measurements of atomic recoils [38] provide now more precise determinations of \(\alpha \), which give values of \(a_e^{\mathrm{SM}}\) such that

$$\begin{aligned} \delta a_e^{\mathrm{Exp,Cs}}=-(8.7\pm 3.6)\times 10^{-13}\, , \end{aligned}$$
(2)

from measurements with \(^{133}\mathrm{Cs}\) [39], and

$$\begin{aligned} \delta a_e^{\mathrm{Exp,Rb}}=(4.8\pm 3.0)\times 10^{-13}\, , \end{aligned}$$
(3)

from measurements with \(^{87}\mathrm{Rb}\) [40].

In reference [41] the possibility to explain the values of \(\delta a_\mu ^{\mathrm{Exp}}\) from the Muon g-2 Brookhaven experiment [1] together with \(\delta a_e^{\mathrm{Exp,Cs}}\) in Eq. (2) was successfully addressed within a subclass of Two Higgs Doublets Models (2HDMs) with general flavor conservation [42, 43]. This was achieved, of course, without conflicting with a large set of high and low energy constraints. The specific model considered, the so-called I-g\(\ell \)FC 2HDM is a 2HDM without tree level scalar flavor changing neutral couplings (SFCNC): in the quark sector it is a type I 2HDM while in the lepton sector it is a general flavor conserving model. The existence of these two anomalies has been addressed in a variety of scenarios, including models with extra Higgs doublets [44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62], models with other scalar extensions [63,64,65,66,67,68,69,70,71,72,73] and supersymmetric models [74,75,76,77,78]. There are also plenty of studies with other approaches such as leptoquarks, vector-like fermions or extra gauge bosons, among others [79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95].

The present work extends and improves several aspects of [41].

  • An improved numerical exploration of the parameter space shows that some unexpected regions of interest can be appropriately covered.

  • Some theoretical assumptions like the perturbativity limits on Yukawa couplings had a significant impact on the analysis and were not fully considered.

  • The latest Muon g-2 Fermilab result [2, 3] consolidates the need of NP brought by the previous Brookhaven result.

  • For \(a_e^{\mathrm{Exp}}\) the situation is rather unclear: within the present scenario, accommodating the values in Eq. (2) or in Eq. (3) may have non-trivial consequences in the model, since they differ in size and in sign.

  • The recent measurement of the W boson mass by the CDF collaboration [96], which disagrees with SM expectations [97], can also be addressed in this context.

All in all, we are entering an era of exclusion or discovery at the LHC and improved analyses of such potential NP hints are necessary.

The manuscript is organized as follows. In Sect. 2, the model is presented. Section 3 is devoted to a discussion of general constraints which apply regardless of \(\delta a_\ell \). The new contributions to \(\delta a_\ell \) are analysed in Sect. 4. The main aspects of the numerical analysis are introduced in Sect. 5. Next, Sect. 6 contains the results of the different analyses together with the corresponding discussions. Finally, the conclusions are presented in Sect. 7. We relegate to the appendices some aspects concerning different sections.

2 Model

The 2HDM is based on the SM gauge group with identical fermion matter contentFootnote 2 and an additional complex scalar doublet. Hence, we have \(\varPhi _{j}\) (\(j=1,2\)) and their corresponding C-conjugate fields defined as \({\tilde{\varPhi }}_{j}\equiv i\sigma _2\varPhi _{j}^*\,\), with opposite sign hypercharge.

The most general scalar potential of 2HDMs can be written as

$$\begin{aligned} \begin{aligned} {\mathscr {V}}(\varPhi _{1}, \varPhi _{2})&= \mu _{11}^2 \varPhi _{1}^\dagger \varPhi _{1} + \mu _{22}^2 \varPhi _{2}^\dagger \varPhi _{2} + (\mu _{12}^2 \varPhi _{1}^\dagger \varPhi _{2} + \text {H.c.}) \\&\quad + \lambda _1 (\varPhi _{1}^\dagger \varPhi _{1})^2 + \lambda _2 (\varPhi _{2}^\dagger \varPhi _{2})^2 \\&\quad + 2 \lambda _3 (\varPhi _{1}^\dagger \varPhi _{1})(\varPhi _{2}^\dagger \varPhi _{2}) + 2 \lambda _4 (\varPhi _{1}^\dagger \varPhi _{2})(\varPhi _{2}^\dagger \varPhi _{1}) \\&\quad + [\lambda _{5} (\varPhi _{1}^\dagger \varPhi _{2})^2 + \text {H.c.}] \\&\quad + [\lambda _{6} (\varPhi _{1}^\dagger \varPhi _{1}) (\varPhi _{1}^\dagger \varPhi _{2}) \\&\quad + \lambda _{7} (\varPhi _{2}^\dagger \varPhi _{2}) (\varPhi _{1}^\dagger \varPhi _{2}) + \text {H.c.}]\, , \end{aligned}\nonumber \\ \end{aligned}$$
(4)

with real \(\mu _{11}^2\), \(\mu _{22}^2\) and \(\lambda _i\) (\(i = 1\ \mathrm {to}\ 4\)), whereas \(\mu _{12}^2\) and \(\lambda _j\) (\(j = 5\ \mathrm {to}\ 7\)) are complex in general. We assume that \({\mathscr {V}}(\langle \varPhi _{1} \rangle ,\langle \varPhi _{2} \rangle )\) has an appropriate minimum at

$$\begin{aligned} \left\langle 0\big |\varPhi _{j}\big |0 \right\rangle =\frac{1}{\sqrt{2}}\begin{pmatrix}0\\ v_j e^{i\theta _j}\end{pmatrix}\, , \end{aligned}$$
(5)

being \(\theta _j\) and \(v_j\) (\(v_j \ge 0\)) real numbers. Taking this into account, the Higgs doublets can be parametrized around the vacuum as

$$\begin{aligned} \varPhi _{j} = e^{i\theta _j} \begin{pmatrix} \varphi _j^{+} \\ (v_j + \rho _j + i\eta _j)/\sqrt{2}\end{pmatrix}. \end{aligned}$$
(6)

IntroducingFootnote 3\(c_\beta \equiv \cos {\beta } \equiv v_1/v\), \(s_\beta \equiv \sin {\beta } \equiv v_2/v\), \(t_\beta \equiv \tan {\beta } = v_2/v_1\), with \(\beta \in [0;\pi /2]\) and \(v^2 = v_1^2 + v_2^2 = (\sqrt{2} G_F)^{-1} \simeq (246\ \mathrm {GeV})^2\), one can perform a global SU(2) rotation in the scalar space and express the scalar doublets in the so-called Higgs basis [98,99,100]

$$\begin{aligned}&\begin{pmatrix}H_1\\ H_2\end{pmatrix}={\mathscr {R}}_{\beta }\, \begin{pmatrix}e^{-i\theta _1}\varPhi _{1}\\ e^{-i\theta _2}\varPhi _{2}\end{pmatrix},\quad \text {with} \nonumber \\&\quad {\mathscr {R}}_{\beta }=\begin{pmatrix} c_\beta &{}\quad s_\beta \\ -s_\beta &{}\quad c_\beta \end{pmatrix}\quad \text {and} \quad {\mathscr {R}}_{\beta }^T= {\mathscr {R}}_{\beta }^{-1}, \end{aligned}$$
(7)

where only one linear combination of the scalar doublets, namely \(H_1\), has a non-zero vacuum expectation value (vev):

$$\begin{aligned} \langle H_1\rangle = \frac{v}{\sqrt{2}}\begin{pmatrix} 0 \\ 1 \end{pmatrix}, \quad \langle H_2\rangle = \begin{pmatrix} 0 \\ 0 \end{pmatrix}. \end{aligned}$$
(8)

The explicit degrees of freedom in this basis are defined by

$$\begin{aligned} H_1= \begin{pmatrix} G^{+} \\ \frac{v+{H}^0+iG^{0}}{\sqrt{2}} \end{pmatrix}, \quad H_2= \begin{pmatrix} \mathrm {H}^+\\ \frac{{R}^0+i{I}^0}{\sqrt{2}} \end{pmatrix}, \end{aligned}$$
(9)

where

$$\begin{aligned}&\begin{pmatrix} G^{+}\\ \mathrm {H}^+\end{pmatrix} = {\mathscr {R}}_{\beta }\begin{pmatrix} \varphi _1^{+}\\ \varphi _2^{+} \end{pmatrix}, \quad \begin{pmatrix} {H}^0\\ {R}^0\end{pmatrix} = {\mathscr {R}}_{\beta }\begin{pmatrix} \rho _1\\ \rho _2 \end{pmatrix}, \nonumber \\&\quad \begin{pmatrix} G^{0}\\ {I}^0\end{pmatrix} = {\mathscr {R}}_{\beta }\begin{pmatrix} \eta _1\\ \eta _2 \end{pmatrix}. \end{aligned}$$
(10)

As we can check, the would-be Goldstone bosons \(G^0\) and \(G^\pm \) get isolated as components of the first Higgs doublet. Likewise, we already identify two charged physical scalars \(\mathrm {H}^\pm \) and three neutral fields \(\{{H}^0,{R}^0,{I}^0\}\) that are not, in general, the mass eigenstates. The latter are determined by the scalar potential, which generates their mass matrix \({\mathscr {M}}_0^2\). This can be diagonalized by a \(3 \times 3\) real orthogonal transformation, \({\mathscr {R}}\), as

$$\begin{aligned} {\mathscr {R}}^T{\mathscr {M}}_0^2{\mathscr {R}}= \mathrm {diag}(m_{\mathrm {h}}^2,m_{\mathrm {H}}^2,m_{\mathrm {A}}^2), \quad {\mathscr {R}}^T= {\mathscr {R}}^{-1}, \end{aligned}$$
(11)

and thus the physical scalars \(\{\mathrm {h},\mathrm {H},\mathrm {A}\}\) are given by

$$\begin{aligned} \begin{pmatrix} \mathrm {h}\\ \mathrm {H}\\ \mathrm {A}\end{pmatrix} = {\mathscr {R}}^T\begin{pmatrix} {H}^0\\ {R}^0\\ {I}^0\end{pmatrix}. \end{aligned}$$
(12)

Neglecting CP violation in the scalar sector, one has

$$\begin{aligned} {\mathscr {R}}= \begin{pmatrix} s_{\alpha \beta }&{}\quad -c_{\alpha \beta }&{}\quad 0 \\ c_{\alpha \beta }&{}\quad s_{\alpha \beta }&{}\quad 0 \\ 0 &{}\quad 0 &{}\quad 1 \end{pmatrix}, \end{aligned}$$
(13)

where \(s_{\alpha \beta }\equiv \sin (\alpha + \beta )\) and \(c_{\alpha \beta }\equiv \cos (\alpha + \beta )\), with \(\pi /2 - \alpha \) being the mixing angle that parametrizes the change of basis from the fields in Eq. (6) to the mass eigenstates in Eq. (12). We should point out that different conventions for Eq. (13) can be found in the literature.

Regarding the Yukawa sector, it is extended to

$$\begin{aligned} \mathscr {L}_{\mathrm{Y}}= & {} -\bar{Q}_{L}^0\left( \varPhi _{1}Y_{d1}+\varPhi _{2}Y_{d2}\right) {d}_{R}^0\nonumber \\&-\bar{Q}_{L}^0\left( {\tilde{\varPhi }}_{1}Y_{u1}+{\tilde{\varPhi }}_{2}Y_{u2}\right) {u}_{R}^0 \nonumber \\&-\bar{L}_{L}^0\left( \varPhi _{1}Y_{\ell 1}+\varPhi _{2}Y_{\ell 2}\right) {\ell }_{R}^0 +\text {H.c.}\, , \end{aligned}$$
(14)

where the couplings \(Y_{dj}\), \(Y_{uj}\) and \(Y_{\ell j}\) \((j = 1,2)\) are \(3 \times 3\) complex matrices in flavor space. One should notice that there are only two flavor structures in the leptonic sector because we are not considering right-handed neutrinos. In the Higgs basis, the Yukawa Lagrangian takes the form

$$\begin{aligned} \mathscr {L}_{{\mathrm{Y}}}= & {} -\frac{\sqrt{2}}{v}\bar{Q}_{L}^0\left( H_1\mathrm{M}_{d}^{0}+H_2\mathrm{N}_{d}^{0}\right) {d}_{R}^0\nonumber \\&-\frac{\sqrt{2}}{v}\bar{Q}_{L}^0\left( {\tilde{H}}_1\mathrm{M}_{u}^{0}+{\tilde{H}}_2\mathrm{N}_{u}^{0}\right) {u}_{R}^0\nonumber \\&-\frac{\sqrt{2}}{v}\bar{L}_{L}^0\left( H_1\mathrm{M}_{\ell }^{0}+H_2\mathrm{N}_{\ell }^{0}\right) {\ell }_{R}^0+\text {H.c.}. \end{aligned}$$
(15)

It is then clear that the matrices \(\mathrm{M}_{f}^{0}\) (\(f = d, u, \ell \)) are the non-diagonal fermion mass matrices since they are coupled to the only Higgs doublet that acquires a non-vanishing vev, i.e., \(H_1\).

The model we are considering in the quark sector is defined by

$$\begin{aligned} Y_{d2} = d Y_{d1}\, , \quad Y_{u2} = d^{*} Y_{u1}\, , \end{aligned}$$
(16)

which is equivalent to

$$\begin{aligned} \mathrm{N}_{d}^{0}=t_\beta ^{-1}\mathrm{M}_{d}^{0}\, , \quad \mathrm{N}_{u}^{0}=t_\beta ^{-1}\mathrm{M}_{u}^{0}\, . \end{aligned}$$
(17)

In the leptonic sector, there exist two unitary matrices \(W_L\) and \(W_R\) such that both \(W_L^\dagger Y_{\ell i} W_R\) (\(i=1,2\)) get simultaneously diagonalized. It is well-known that the structure of the quark sector can be enforced through a \({\mathbb {Z}}_{2}\) symmetry, but this is not the case in the lepton sector. Nevertheless, as it is shown in Appendix A, the entire Yukawa structure is stable under one loop renormalization group evolution (RGE) and, therefore, the model is free from unwanted SFCNC.

Going to the fermion mass bases for our I-g\(\ell \)FC model – type I in the quark sector and general flavor conserving in the lepton sector – we get the relevant new Yukawa structures:

$$\begin{aligned} \mathscr {L}_{{\mathrm{Y}}}= & {} -\frac{\sqrt{2}}{v}\bar{Q}_{L}\left( H_1\mathrm{M}_{d}+H_2\mathrm{N}_{d}\right) {d}_{R}\nonumber \\&-\frac{\sqrt{2}}{v}\bar{Q}_{L}\left( {\tilde{H}}_1\mathrm{M}_{u}+{\tilde{H}}_2\mathrm{N}_{u}\right) {u}_{R}\nonumber \\&-\frac{\sqrt{2}}{v}\bar{L}_{L}\left( H_1\mathrm{M}_{\ell }+H_2\mathrm{N}_{\ell }\right) {\ell }_{R}+\text {H.c.}\, , \end{aligned}$$
(18)

with

$$\begin{aligned}&\mathrm{N}_{d}=t_\beta ^{-1}\mathrm{M}_{d}, \mathrm{N}_{u}=t_\beta ^{-1}\mathrm{M}_{u}, \mathrm{N}_{\ell }=\mathrm {diag}(n_{e},n_{\mu },n_{\tau })\, , \end{aligned}$$
(19)

and \(\mathrm{M}_{f}\) (\(f = u,d,\ell \)) the corresponding diagonal fermion mass matrices. Note that the quark couplings \(\mathrm{N}_{u}\) and \(\mathrm{N}_{d}\) are those from 2HDMs of type I or X. On the other hand, the matrices \(\mathrm{N}_{\ell }\) correspond to a general flavor conserving lepton sector. Therefore, they are diagonal, arbitrary and one loop stable under RGE, as it was shown in [43], meaning that they remain diagonal.

We must stress that it is the fact that \(n_{e}\) and \(n_{\mu }\) are completely independent what implements the desired decoupling between electron and muon NP couplings in order to have enough freedom to address the corresponding \((g-2)_\ell \) anomalies. We assume that these couplings are real, i.e., \(\mathrm {Im}(n_{\ell }) = 0\). This prevents us from dangerous contributions to electric dipole moments (EDMs), that are tightly constrained: \(|d_e| < 1.1 \times 10^{-29}\ \text {e}\cdot \text {cm}\) [101].

Furthermore, we consider an scalar potential shaped by a \({\mathbb {Z}}_{2}\) symmetry that is softly broken by the term \(\mu _{12}^2 \ne 0\). Hence, we have to take \(\lambda _6 = \lambda _7 = 0\) in Eq. (4). We also assume that there is no CP violation in the scalar sector, so Eq. (13) is fulfilled.

Fig. 1
figure 1

Oblique parameters: allowed regions in \(m_{\mathrm {A}}-m_{\mathrm {H}^\pm }\) vs. \(m_{\mathrm {H}}-m_{\mathrm {H}^\pm }\). Darker to lighter colors correspond to 2D-\(\varDelta \chi ^2\) 1, 2 and 3\(\sigma \) regions. The plot corresponds to \(m_{\mathrm {H}^\pm }=1\) TeV and scalar alignment

Under these assumptions, the flavor conserving Yukawa interactions of neutral scalars read

$$\begin{aligned} \begin{aligned} \mathscr {L}_{{\mathrm{N}}}&= -\frac{m_{u_j}}{v} \left( s_{\alpha \beta }+ c_{\alpha \beta }t_\beta ^{-1}\right) \mathrm {h}\, {\bar{u}}_j u_j \\&\quad -\frac{m_{d_j}}{v} \left( s_{\alpha \beta }+ c_{\alpha \beta }t_\beta ^{-1}\right) \mathrm {h}\, {\bar{d}}_j d_j \\&\quad -\frac{m_{\ell _j}}{v} \left( s_{\alpha \beta }+ c_{\alpha \beta }\frac{\text {Re}\left( n_{\ell _j}\right) }{m_{\ell _j}} \right) \mathrm {h}\, {\bar{\ell }}_j \ell _j \\&\quad -\frac{m_{u_j}}{v} \left( -c_{\alpha \beta }+ s_{\alpha \beta }t_\beta ^{-1}\right) \mathrm {H}\, {\bar{u}}_j u_j \\&\quad -\frac{m_{d_j}}{v} \left( -c_{\alpha \beta }+ s_{\alpha \beta }t_\beta ^{-1}\right) \mathrm {H}\, {\bar{d}}_j d_j \\&\quad -\frac{m_{\ell _j}}{v} \left( -c_{\alpha \beta }+ s_{\alpha \beta }\frac{\text {Re}\left( n_{\ell _j}\right) }{m_{\ell _j}} \right) \mathrm {H}\, {\bar{\ell }}_j \ell _j \\&\quad +i\frac{m_{u_j}}{v} t_\beta ^{-1}\mathrm {A}\, {\bar{u}}_j \gamma _5 u_j -i\frac{m_{d_j}}{v} t_\beta ^{-1}\mathrm {A}\, {\bar{d}}_j \gamma _5 d_j \\&\quad -i\frac{\text {Re}\left( n_{\ell _j}\right) }{v} \mathrm {A}\, {\bar{\ell }}_j \gamma _5 \ell _j\, , \\ \end{aligned} \end{aligned}$$
(20)

and those involving charged scalars are

$$\begin{aligned} \mathscr {L}_{{\mathrm{Ch}}}= & {} \frac{\mathrm {H}^-}{\sqrt{2}v}\, {\bar{d}}_i\, V_{ji}^{*} t_\beta ^{-1}\left[ (m_{u_j}-m_{d_i}) + (m_{u_j}+m_{d_i}) \gamma _5\right] u_j \nonumber \\&+ \frac{\mathrm {H}^+}{\sqrt{2}v}\, {\bar{u}}_j\, V_{ji} t_\beta ^{-1}\left[ (m_{u_j}-m_{d_i}) - (m_{u_j}+m_{d_i}) \gamma _5\right] d_i \nonumber \\&- \frac{\mathrm {H}^-}{\sqrt{2}v}\, {\bar{\ell }}_j\, \text {Re}\left( n_{\ell _j}\right) \left( 1-\gamma _5\right) \nu _j \nonumber \\&- \frac{\mathrm {H}^+}{\sqrt{2}v}\, {\bar{\nu }}_j\, \text {Re}\left( n_{\ell _j}\right) \left( 1+\gamma _5\right) \ell _j\, , \nonumber \\ \end{aligned}$$
(21)

with \(i,j = 1,2,3\) summing over generations. It is easy to check that \(\mathrm {h}\) presents the same couplings as the SM Higgs boson when we take the scalar alignment limit, i.e., \(s_{\alpha \beta }\rightarrow 1\).

3 General constraints

Before addressing the different contributions to the anomalous magnetic moments \(\delta a_\ell \), we discuss in this section some general constraints which are relevant in the scenario under consideration. By “general” we mean that they do not depend specifically on the values of \(\text {Re}\left( n_{e}\right) \), \(\text {Re}\left( n_{\mu }\right) \), \(\delta a_e\) and \(\delta a_\mu \). Furthermore, their effects can be understood in simple terms.

  • Alignment. The couplings of the scalar \(\mathrm {h}\), assumed to be the SM-Higgs-like particle with \(m_{\mathrm {h}}=125\) GeV, deviate from SM values through the scalar mixing in Eq. (13). Measurements of the signal strengths in the usual set of production mechanisms and decay channels impose \(c_{\alpha \beta }\ll 1\). Concerning the scalar sector, we are thus in the alignment limit.

  • \(\underline{\hbox {Oblique parameters and}\, M_W}\). Electroweak precision measurements constrain deviations in the oblique parameters S and T [97, 102]:

    $$\begin{aligned}&\varDelta S=0.00\pm 0.07,\quad \varDelta T=0.05\pm 0.06,\nonumber \\&\quad \rho =0.92 \ \text {(correlation)}. \end{aligned}$$
    (22)

    In 2HDMs, in the alignment limit mentioned above, one can observe that the corrections to S and T are kept under control when either \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\) or \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\), as shown in Fig. 1a. Recently, the CDF collaboration announced a measurement of the W boson mass which disagrees with SM expectations [96]. In fits of electroweak precision observables this disagreement can be translated into values of the oblique parameters \((\varDelta S,\varDelta T)\ne (0,0)\) [103, 104] (although fits including \(\varDelta U\) have also been considered, we focus on the case \(\varDelta U=0\), appropriate here). In order to “explain” the CDF \(M_W\) “anomaly” one can thus consider \((\varDelta S,\varDelta T)\) constraints from [103, 104] instead of Eq. (22). We can consider, in particular,

    1. (i)

      the “conservative scenario” in [103] which combines the CDF with previous measurements and gives

      $$\begin{aligned}&\varDelta S=0.086 \pm 0.077,\quad \varDelta T=0.177\pm 0.070,\nonumber \\&\quad \rho =0.89, \end{aligned}$$
      (23)
    2. (ii)

      the results in [104] which solely use the CDF measurement and give

      $$\begin{aligned}&\varDelta S=0.15 \pm 0.08,\quad \varDelta T=0.27\pm 0.06,\nonumber \\&\quad \rho =0.93. \end{aligned}$$
      (24)

    In the alignment limit, for \(m_{\mathrm {H}^\pm }=1\) TeV, Eqs. (23) and (24) give the allowed regions represented in Fig. 1b, c respectively. In sharp contrast with Fig. 1a, notice in Fig. 1b, c how near degeneracy of the three new scalars is excluded, and how even near degeneracies \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\) or \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\) are quite disfavored. Furthermore, notice that the \(1\sigma \) region (2D-\(\varDelta \chi ^2\le 2.23\)) does not appear in Fig. 1c: contrary to Eq. (23), with Eq. (24) one cannot obtain the minimum \(\chi ^2_{\mathrm{Min}}\) with \(m_{\mathrm {H}^\pm }=1\) TeV.

  • \(\underline{\mathrm {H}^\pm -\hbox {induced FCNC}}\). The charged scalar \(\mathrm {H}^\pm \) can contribute to \(\varDelta F=1\) and \(\varDelta F=2\) FCNC processes like \(b\rightarrow s\gamma \) and \(B_q\)\({{\bar{B}}}_q\) mixings (for example, through SM-like box diagrams for \(B_q\)\({{\bar{B}}}_q\) in which \(W^\pm \) are replaced with \(\mathrm {H}^\pm \)). The dominant contributions involve virtual top quarks as in the SM, with couplings including now \(t_\beta ^{-1}\) factors. Keeping those contributions within experimental bounds only allows, roughly, the colored region in Fig. 2. For each value of \(m_{\mathrm {H}^\pm }\) there is a lower bound on \(t_\beta \). See [105,106,107] for further details.

  • Scalar sector perturbativity. Additional constraints on scalar masses vs. \(t_\beta \) arise from perturbativity requirements on the quartic coefficients of the scalar potential and from perturbative unitarity of \(2\rightarrow 2\) scattering amplitudes [108,109,110,111,112,113,114]. With a \({\mathbb {Z}}_{2}\) symmetric potential, it is difficult to obtain masses above 1 TeV and values of \(t_\beta \) larger than 8. Larger values of the masses and larger values of \(t_\beta \) can be nevertheless obtained with the introduction of a soft symmetry breaking term \(\mu _{12}^2\ne 0\) in Eq. (4) [114, 115].

  • Gluon-gluon fusion production cross section. Let us consider the production cross section of \(\mathrm {H}\) and \(\mathrm {A}\) through the one loop gluon-gluon fusion process. In the scalar alignment limit, one can read from Eq. (20) that the same \(t_\beta ^{-1}\) factor applies to both pure scalar \(\mathrm {H}\) and pure pseudoscalar \(\mathrm {A}\) couplings with the top quark in the triangle loop:

    $$\begin{aligned}&\sigma (pp\rightarrow S)_{\mathrm{ggF}}\propto t_\beta ^{-2}|F_{S}(x)|^2,\ \nonumber \\&x=4m_t^2/m_S^2,\quad S=\mathrm {H},\mathrm {A}. \end{aligned}$$
    (25)

    The corresponding loop functions \(F_{\mathrm {H}}\) and \(F_{\mathrm {A}}\) [116,117,118,119,120,121] are different due to the scalar or pseudoscalar character:

    $$\begin{aligned} \begin{aligned}&F_{\mathrm {H}}(x)=-2x(1+(1-x)f(x)),\\&F_{\mathrm {A}}(x)=-2xf(x),\\&\quad f(x)= \left\{ \begin{aligned}&\arcsin ^2(1/\sqrt{x}),\quad&x\ge 1\\&-\frac{1}{4}\left( \ln \left( \frac{1+\sqrt{1-x}}{1-\sqrt{1-x}}\right) -i\pi \right) ^2,&x<1 \end{aligned}\right\} . \end{aligned}\nonumber \\ \end{aligned}$$
    (26)

    Figure 3 shows \(|F_{\mathrm {H}}(x)|^2\), \(|F_{\mathrm {A}}(x)|^2\) and the ratio \(|F_{\mathrm {A}}(x)|^2 /|F_{\mathrm {H}}(x)|^2\) as a function of the scalar mass. It is clear that the pseudoscalar \(\mathrm {A}\) has a larger gluon-gluon production cross section than the scalar \(\mathrm {H}\) for \(m_{\mathrm {A}}=m_{\mathrm {H}}\) (up to a factor of 6 for \(m_{\mathrm {A}}=m_{\mathrm {H}}=2m_t\)). Since dimuon searches \([pp]_{\mathrm{ggF}}\rightarrow S\rightarrow \mu ^+\mu ^-\) at the LHC can be rather constraining for scalar masses \(m_{\mathrm{S}}<1\) TeV, one can expect that in that low mass region \(m_{\mathrm {A}}> m_{\mathrm {H}}\). One could have worried about the validity of this expectation in case \(\text {Br}\left( \mathrm {A}\rightarrow \mu ^+\mu ^-\right) \ll \text {Br}\left( \mathrm {H}\rightarrow \mu ^+\mu ^-\right) \), but the only way to achieve a suppression of \(\text {Br}\left( \mathrm {A}\rightarrow \mu ^+\mu ^-\right) \) relative to \(\text {Br}\left( \mathrm {H}\rightarrow \mu ^+\mu ^-\right) \) is through the existence of \(\mathrm {A}\rightarrow \mathrm {H}Z\) decays, which are only available if \(m_{\mathrm {A}}> m_{\mathrm {H}}\), and thus cannot change that expectation.

  • \(\underline{e^+e^-\rightarrow \mu ^+\mu ^-\, \hbox {at LEP}}\). Sizable \(n_e\) and \(n_\mu \) are necessary ingredients for the contributions to \(a_e\) and \(a_\mu \) involving the new scalars \(\mathrm {H}\), \(\mathrm {A}\) and \(\mathrm {H}^\pm \). Data from LEP [122] on \(e^+e^-\rightarrow \mu ^+\mu ^-\) with \(\sqrt{s}\) up to 210 GeV are sensitive to s-channel \(\mathrm {H}\) and \(\mathrm {A}\) mediated contributions (contrary to the LHC gluon-gluon fusion process, being scalar or pseudoscalar does not change the sensitivity of LEP data). One can roughly expect that agreement with LEP data imposes \(m_{\mathrm {A}},m_{\mathrm {H}}\ge 210\) GeV.

Fig. 2
figure 2

\(\mathrm {H}^\pm \) FCNC: \(m_{\mathrm {H}^\pm }\) vs. \(t_\beta \) allowed region when contributions of \(\mathrm {H}^\pm \) to \(B_q\)\({{\bar{B}}}_q\) are below experimental uncertainty in \(\varDelta M_{B_q}\)

Fig. 3
figure 3

Loop functions controlling gluon-gluon production cross sections of scalars

4 Contributions to \(\delta a_\ell \)

The complete prediction of the anomalous magnetic moment \(a_\ell ^{\mathrm{Th}}\), \(\ell =e,\mu \), is

$$\begin{aligned} a_{\ell }^{\mathrm{Th}}=a_{\ell }^{\mathrm{SM}}+\delta a_{\ell }\, , \end{aligned}$$
(27)

where \(a_{\ell }^{\mathrm{SM}}\) is the SM contribution and \(\delta a_\ell \) the NP correction. The anomalies in Eqs. (1)–(3) are “solved” for \(\delta a_e=\delta a_e^{\mathrm{Exp}}\) and \(\delta a_\mu =\delta a_\mu ^{\mathrm{Exp}}\). We introduce for convenience \(\varDelta _\ell \) such that

$$\begin{aligned} \delta a_\ell =K_\ell \,\varDelta _\ell ,\qquad K_\ell =\frac{1}{8\pi ^2}\left( \frac{m_\ell }{v}\right) ^2. \end{aligned}$$
(28)

For \(\delta a_\mu \) one needs

$$\begin{aligned} \varDelta _\mu \simeq 1\, , \end{aligned}$$
(29)

while for \(\delta a_e\) one needs

$$\begin{aligned} \varDelta _e^{\mathrm{Cs}}\simeq -16\, ,\qquad \varDelta _e^{\mathrm{Rb}}\simeq 9, \end{aligned}$$
(30)

where the superscript corresponds to the different values in Eqs. (2) and (3).

In the model considered here, it is well known that both one loop and two loop (of Barr–Zee type) contributions can be dominant. In this section we analyse both types of contributions in the scalar alignment limit \(s_{\alpha \beta }\rightarrow 1\) and keeping only leading terms in \(\frac{m_\ell ^2}{m_{S}^2}\), \({S}=\mathrm {H},\mathrm {A},\mathrm {H}^\pm \). Full results, used for instance in the numerical analyses, can be found in Appendix B.

4.1 One loop contributions to \(\delta a_\ell \)

The one loop result \(\varDelta _\ell ^{(1)}\) has contributions from \(\mathrm {H}\), \(\mathrm {A}\) and \(\mathrm {H}^\pm \). With the approximations mentioned above and the couplings in Eqs. (20) and (21), we have

$$\begin{aligned} \varDelta _\ell ^{(1)}\simeq |n_{\ell }|^2\left( \frac{I_{\ell \mathrm {H}}}{m_{\mathrm {H}}^2}-\frac{I_{\ell \mathrm {A}}-2/3}{m_{\mathrm {A}}^2}-\frac{1}{6m_{\mathrm {H}^\pm }^2}\right) \, , \end{aligned}$$
(31)

where

$$\begin{aligned} I_{\ell S}=-\frac{7}{6}-2\ln \left( \frac{m_\ell }{m_S}\right) . \end{aligned}$$
(32)

The range of interest in our analyses will be \(m_S\in [0.2;2.5]\) TeV, in which case

$$\begin{aligned} I_{\mu S}\in [13.9;18.9], \end{aligned}$$
(33)

while

$$\begin{aligned} I_{eS}=I_{\mu S}+2\ln \left( \frac{m_\mu }{m_e}\right) =I_{\mu S}+10.7. \end{aligned}$$
(34)

In Eq. (31), the \(\mathrm {H}\) contribution is positive, the \(\mathrm {A}\) contribution is negative and the \(\mathrm {H}^\pm \) contribution is negligible. One can then anticipate the following.

  • The muon anomaly \(\varDelta _\mu \simeq 1\) can only be explained with the one loop \(\mathrm {H}\) contribution and provided

    $$\begin{aligned} 1 \simeq \frac{|n_{\mu }|^2}{m_{\mathrm {H}}^2} I_{\mu \mathrm {H}}\ \Rightarrow \ |n_{\mu }|\sim \frac{1}{4}m_{\mathrm {H}}. \end{aligned}$$
    (35)

    Considering \(|n_{\mu }|<250\) GeV, a priori there could be a one loop explanation of \(\delta a_\mu \) for \(m_{\mathrm {H}}< 1\) TeV. Since the \(\mathrm {A}\) contribution has opposite sign, if \(m_{\mathrm {A}}\sim m_{\mathrm {H}}\) a substantial cancellation would occur. As discussed in Sect. 3, it is precisely for light \(\mathrm {H}\) that one expects \(m_{\mathrm {A}}>m_{\mathrm {H}}\), in which case that cancellation is largely avoided and a one loop \(\mathrm {H}\) explanation viable. For heavier \(m_{\mathrm {H}}\), the muon anomaly needs other contributions.

  • For the electron \(\mathrm{Cs}\) anomaly, \(\varDelta _e^{\mathrm{Cs}}\simeq -16\) can only be explained with the one loop \(\mathrm {A}\) contribution provided

    $$\begin{aligned} -16\simeq -\frac{|n_{e}|^2}{m_{\mathrm {A}}^2} I_{e\mathrm {A}}\ \Rightarrow \ |n_{e}|\sim \frac{4}{5}m_{\mathrm {A}}. \end{aligned}$$
    (36)

    For \(|n_{e}|<250\) GeV, this would require the pseudoscalar \(\mathrm {A}\) to be rather light, \(m_{\mathrm {A}}<300\) GeV. On the other hand, \(m_{\mathrm {A}}>200\) GeV would require \(|n_{e}|>160\) GeV: besides perturbativity concerns, such values of \(|n_e|\) might be hard to reconcile with other constraints. More importantly, since we expect \(m_{\mathrm {H}}<m_{\mathrm {A}}\) for light \(\mathrm {A}\), we also expect a sizable cancellation among \(\mathrm {H}\) and \(\mathrm {A}\) contributions. From this simple analysis, obtaining \(\varDelta _e^{\mathrm{Cs}}\simeq -16\) with one loop contributions does not appear to be feasible.

  • For the electron \(\mathrm{Rb}\) anomaly, \(\varDelta _e^{\mathrm{Rb}}\simeq 9\) can only be explained with the one loop \(\mathrm {H}\) contribution and provided

    $$\begin{aligned} 9\simeq \frac{|n_{e}|^2}{m_{\mathrm {H}}^2} I_{e\mathrm {H}}\ \Rightarrow \ |n_{e}|\sim \frac{3}{5}m_{\mathrm {H}}. \end{aligned}$$
    (37)

    For \(m_{\mathrm {H}}>200\) GeV, this would require \(|n_{e}|>120\) GeV. If the same concerns on the values of \(|n_{e}|\) mentioned for \(\varDelta _e^{\mathrm{Cs}}\simeq -16\) apply here, obtaining \(\varDelta _e^{\mathrm{Rb}}\simeq 9\) does not seem to be feasible neither; otherwise \(\varDelta _e^{\mathrm{Rb}}\simeq 9\) would be “easier” to accommodate with one loop contributions than \(\varDelta _e^{\mathrm{Cs}}\simeq -16\) because of the sign difference and the smaller absolute value.

4.2 Two loop contributions to \(\delta a_\ell \)

The dominant two loop contributions are the Barr–Zee ones. Diagrammatically they correspond to contributions where a closed fermion loop is attached to the external lepton through two propagators: one photon and one of the new scalars \(\mathrm {H}\), \(\mathrm {A}\). In the scalar alignment limit,

$$\begin{aligned} \varDelta _\ell ^{(2)}=-\frac{2\alpha }{\pi }\,\frac{\text {Re}\left( n_{\ell }\right) }{m_\ell }\,F. \end{aligned}$$
(38)

It is important to notice that these contributions are linear in \(n_{\ell }\). Detailed expressions are provided in Appendix B. In Eq. (38) we have

$$\begin{aligned} \frac{2\alpha }{\pi m_e}\simeq 9.1\,\text {GeV}^{-1},\qquad \frac{2\alpha }{\pi m_\mu }\simeq 0.044\,\text {GeV}^{-1}. \end{aligned}$$
(39)

The function F depends on the masses of the fermions in the closed loop, their couplings to \(\mathrm {H}\) and \(\mathrm {A}\), and on \(m_{\mathrm {A}}\) and \(m_{\mathrm {H}}\). Considering the dominant contributions from top and bottom quarks, and also from tau and muon leptons since \(n_{\tau }\) and \(n_{\mu }\) are free parameters,

$$\begin{aligned} F= & {} \frac{t_\beta ^{-1}}{3}\left[ 4(f_{t\mathrm {H}}+g_{t\mathrm {A}})+(f_{b\mathrm {H}}-g_{b\mathrm {A}})\right] \nonumber \\&+\frac{\text {Re}\left( n_{\tau }\right) }{m_\tau }(f_{\tau \mathrm {H}}-g_{\tau \mathrm {A}})+\frac{\text {Re}\left( n_{\mu }\right) }{m_\mu }(f_{\mu \mathrm {H}}-g_{\mu \mathrm {A}}), \end{aligned}$$
(40)

with

$$\begin{aligned} f_{xS}=f\left( \frac{m^2_{x}}{m^2_{S}}\right) ,\quad g_{xS}=g\left( \frac{m^2_{x}}{m^2_{S}}\right) . \end{aligned}$$
(41)

The functions f(z) and g(z) are defined in Appendix B. It is to be noticed that (i) \(f(z)\sim g(z)\) in the range of interest, (ii) larger values correspond to heavier fermions, (iii) for the top quark loop, f and g vary between 0.08 and 1 in the relevant range of scalar masses, \(m_S\in [0.2;2.5]\) TeV.

  • If the electron anomaly is to be obtained through the two loop contributions,

    $$\begin{aligned} \varDelta _e\simeq -9.1\,F\,\text {Re}\left( n_{e}\right) /\text {GeV}, \end{aligned}$$
    (42)

    and thus

    $$\begin{aligned} \begin{aligned}&\text {from } \varDelta _e^{\mathrm{Cs}},\quad \text {Re}\left( n_{e}\right) \,F\simeq 1.8\,\text {GeV},\\&\text {from } \varDelta _e^{\mathrm{Rb}},\quad \text {Re}\left( n_{e}\right) \,F\simeq -1.0\,\text {GeV}. \end{aligned} \end{aligned}$$
    (43)

    The sign and the magnitude of F is fixed by the \(\text {Re}\left( n_{e}\right) \) value to fix \(\delta a_e\).

  • For \(m_{\mathrm {H}}>1\) TeV, two loop contributions are necessary to explain the muon anomaly, in which case

    $$\begin{aligned} \varDelta _\mu&\simeq -0.044\,F\,\text {Re}\left( n_{\mu }\right) /\text {GeV}\ \Rightarrow \ \text {Re}\left( n_{\mu }\right) \,\nonumber \\ F&\simeq -23\,\text {GeV}. \end{aligned}$$
    (44)

    If follows that, for \(m_{\mathrm {H}}>1\) TeV,

    $$\begin{aligned} \begin{aligned}&\text {for } \varDelta _e^{\mathrm{Cs}}\text { and }\varDelta _\mu ,\quad \text {Re}\left( n_{\mu }\right) \sim -13\text {Re}\left( n_{e}\right) ,\\&\text {for } \varDelta _e^{\mathrm{Rb}}\text { and }\varDelta _\mu ,\quad \text {Re}\left( n_{\mu }\right) \sim 23\text {Re}\left( n_{e}\right) . \end{aligned} \end{aligned}$$
    (45)

These correlations show that, in the present framework, the independence of \(n_{e}\) and \(n_{\mu }\) is essential to explain the different sign of \(\varDelta _e^{\mathrm{Cs}}\) and \(\varDelta _\mu \). This sign difference is challenging for many scenarios addressing simultaneously both anomalies. In this sense, addressing \(\varDelta _e^{\mathrm{Rb}}\) and \(\varDelta _\mu \) is less challenging.

5 Analysis

In Sect. 3 we have discussed some general constraints that apply without regard to the values of \(n_{e}\) and \(n_{\mu }\) of interest to reproduce the \(\delta a_\ell \) anomalies; in Sect. 4 we have explored the obtention of the \(\delta a_\ell \) anomalies through one and two loop contributions. It is now time to present the main aspects of our detailed numerical analyses. The goal of the numerical analyses is to explore the parameter space of the model and map the different regions where a chosen set of relevant constraints is satisfied and the \(\delta a_\ell \) anomalies are explained in terms of the new contributions. The independent parameters of the model are \(\{t_\beta ,m_{\mathrm {H}},m_{\mathrm {A}},m_{\mathrm {H}^\pm },\mu _{12}^2,c_{\alpha \beta },\text {Re}\left( n_{e}\right) \!,\text {Re}\left( n_{\mu }\right) \!,\text {Re}\left( n_{\tau }\right) \}\): \(\{t_\beta , m_{\mathrm {H}},m_{\mathrm {A}},m_{\mathrm {H}^\pm },\mu _{12}^2,c_{\alpha \beta }\}\) control the scalar sector (together with v and \(m_{\mathrm {h}}\)) while \(\{\text {Re}\left( n_{e}\right) ,\text {Re}\left( n_{\mu }\right) ,\text {Re}\left( n_{\tau }\right) \}\) give the lepton Yukawa couplings (quark Yukawa couplings are fixed by \(t_\beta \)). The set of relevant constraints includes the following.

  • Boundedness from below of the scalar potential [123], perturbativity of quartic couplings and perturbative unitarity of high energy \(2\rightarrow 2\) scattering in the scalar sector [110].

  • Corrections to the oblique parameters S and T in agreement with electroweak precision data [97, 102].

  • “Production \(\times \) decay” predictions for \(\mathrm {h}\) in agreement with Higgs signal strengths [124,125,126,127,128,129,130,131,132,133,134,135,136]. As already mentioned, this constraint forces the alignment limit in the scalar sector: in the analyses of Sect. 6 one obtains indeed \(c_{\alpha \beta }<3\times 10^{-3}\).

  • Lepton flavor universality in leptonic and semileptonic decays [97, 137, 138].

  • \(b\rightarrow s\gamma \) and \(B_q^0\)\({{\bar{B}}}_q^0\) data [97, 105, 106].

  • \(e^+e^-\rightarrow \mu ^+\mu ^-\), \(\tau ^+\tau ^-\) data from LEP (with center of mass energies up to 210 GeV) [122].

  • LHC searches: resonant \(pp\rightarrow S\rightarrow \mu ^+\mu ^-,\tau ^+\tau ^-\) searches with gluon-gluon fusion \(pp\rightarrow S\) production [139,140,141,142,143] and \(\mathrm {H}^\pm \) searches in \(pp\rightarrow \mathrm {H}^\pm tb\), \(\mathrm {H}^\pm \rightarrow \tau \nu ,tb\) [144,145,146,147].

For additional details on the different constraints we refer to [41]. The constraints are typically modelled with a gaussian likelihood or an equivalent \(\chi ^2\) term, the overall likelihood is sampled over parameter space using Markov chain Monte Carlo techniques in order to obtain the regions where (best) agreement with the constraints is obtained. There are two final aspects of central importance which require a specific discussion: (i) how are the anomalies included in the analyses, (ii) what ranges are considered for the \(n_\ell \) parameters.

Concerning the \(a_\ell \) anomalies, the situation for \(\delta a_\mu ^{\mathrm{Exp}}\) is clear: one should consider Eq. (1). On the contrary, for \(\delta a_e^{\mathrm{Exp}}\) the situation is not settled: we have Eqs. (2) and (3), which are rather incompatible. In order to have a complete picture, we analyse both cases separately. Furthermore, we also consider two additional possibilities concerning \(\delta a_e^{\mathrm{Exp}}\):

  • despite the marginal compatibility of \(\delta a_e^{\mathrm{Exp,Cs}}\) and \(\delta a_e^{\mathrm{Exp,Rb}}\), we combine them into

    $$\begin{aligned} \delta a_e^{\mathrm{Exp,Avg}}=-(2.0\pm 2.2)\times 10^{-13}, \end{aligned}$$
    (46)

    which has the same sign as \(\delta a_e^{\mathrm{Exp,Cs}}\), i.e. opposite to \(\delta a_\mu ^{\mathrm{Exp}}\), but a size close to 4 times smaller;

  • a conservative approach in which we only assume that \(|\delta a_e|\le 20\times 10^{-13}\). Rather than targeting a specific value, this analysis may help to single out regions of parameter space where one cannot reproduce \(\delta a_\mu ^{\mathrm{Exp}}\) together with any value of \(\delta a_e\) compatible with \(\delta a_e^{\mathrm{Exp,Cs}}\) or \(\delta a_e^{\mathrm{Exp,Rb}}\).

We will refer to these separate analyses as “\(a_e^{\mathrm{Cs}}\)”, “\(a_e^{\mathrm{Rb}}\)”, “\(a_e^{\mathrm{Avg}}\)”, “\(a_e^{\mathrm{Bound}}\)”. For their implementation in the analyses, we assign a joint \(\chi ^2\) contribution (corresponding to a gaussian factor in the likelihood)

$$\begin{aligned} \chi ^2_{g-2}(\delta a_e,\delta a_\mu )=\left( \frac{\delta a_e-c_e}{s_e}\right) ^2+\left( \frac{\delta a_\mu -c_\mu }{s_\mu }\right) ^2, \end{aligned}$$
(47)

where \(c_\ell \) is the experimental central value and \(s_\ell \) is the experimental uncertainty divided by 4. The scope of this choice – dividing the experimental uncertainty by 4 instead of simply using the experimental uncertainty – is to show clearly that the model can reproduce easily and simultaneously both the muon and the electron anomalies, and to guarantee that we are definitely reproducing a sizable deviation from the SM both in \(a_\mu \) and in all cases for \(a_e\), except the “\(a_e^{\mathrm{Bound}}\)” analysis where there is no \(\delta a_e\) term in Eq. (47) and \(|\delta a_e|\le 20\times 10^{-13}\) is imposed. As a summary, all four selected cases of \(\delta a_\mu ^{\mathrm{Exp}}\) vs. \(\delta a_e^{\mathrm{Exp}}\) are represented in Fig. 4.

Fig. 4
figure 4

Allowed \(\delta a_\mu ^{\mathrm{Exp}}\) vs. \(\delta a_e^{\mathrm{Exp}}\) regions in the different analyses

The different colored regions in Fig. 4 represent three contours in the joint \(\varDelta \chi ^2 = \chi ^2 - \chi ^2_{\mathrm{{Min}}}\) considering only Eq. (47). In a 2D-\(\varDelta \chi ^2\) distribution they correspond, darker to lighter, to 1, 2 and 3\(\sigma \) regions with 68.2% C.L., 95.4% C.L. and 99.7% C.L., respectively. The same color coding is used in the figures below illustrating the final results of the analyses, where all observables and constraints have been included.

Fig. 5
figure 5

Allowed regions: relevant correlations involving \(m_{\mathrm {H}}\) with \(|n_{\ell }|\le 100\) GeV

Finally, in [41], only \(|n_{\ell }|\le 100\) GeV were considered. Although for \(|n_{e}|,|n_{\mu }|\sim 100\) GeV lepton couplings to the new scalars are hugely enhanced with respect to \(\mathrm {h}\) couplings, it is true that \(\frac{|n_{\ell }|}{v}\sim \frac{100\text { GeV}}{v}\sim 0.4\) does not appear to pose a perturbativity challenge. In fact, the one loop correction to the imaginary part of the \(m_{\mathrm {H}}\) mass is controlled by \(\varGamma (\mathrm {H}\rightarrow \ell {{\bar{\ell }}})\) and the relevant ratio is \(\frac{\varGamma }{m_{\mathrm {H}}}=\frac{1}{8\pi }\frac{|n_{\ell }|^2}{v_{}^2}\), therefore arriving to \(|n_{\ell }|=v\sim 250\) GeV represents one loop corrections at the 4% level. For this reason the analyses have been done with \(|n_{\ell }|\le 250\) GeV \(\sim v\); furthermore, the analysis “\(a_e^{\mathrm{Cs}}\)” has been conducted both with \(|n_{\ell }|\le 100\) GeV (since this case is the closest one to [41]) and with \(|n_{\ell }|\le 250\) GeV.

6 Results

In the next subsections we discuss the most relevant results of the analyses done following the lines of the previous section. In Sect. 6.1 we consider the scenario “\(a_e^{\mathrm{Cs}}\)” when \(|n_{\ell }|\le 100\) GeV is imposed. The implications of changing this last assumption to \(|n_{\ell }|\le 250\) GeV are addressed in Sect. 6.2. The implications of the different assumptions for the electron anomaly, that is scenarios “\(a_e^{\mathrm{Rb}}\)”, “\(a_e^{\mathrm{Avg}}\)” and “\(a_e^{\mathrm{Bound}}\)” are explored in Sect. 6.3. The impact of the recent measurement of \(M_W\) by the CDF collaboration is considered in Sect. 6.4. Finally, to further illustrate these discussions, a few complete example cases are shown in Sect. 6.5.

6.1 \(|n_{\ell }|\le 100\) GeV

Here we present the results of the analysis “\(a_e^{\mathrm{Cs}}\)” with the perturbativity constraint \(|n_{\ell }|\le 100\) GeV. This serves to revisit the main results of [41] and as a reference for the analysis with \(|n_{\ell }|\le 250\) GeV addressed in the following subsection.

The perturbativity constraint limits the possibility of explaining \(\delta a_e^{\mathrm{Exp,Cs}}\) via the one loop contribution, since it requires \(m_{\mathrm {A}}\le 125\) GeV for \(|n_{e}|\le 100\) GeV (see Eq. (36)) which is not allowed by \(e^+e^-\rightarrow \mu ^+\mu ^-\) LEP data. On that respect, lepton flavor universality constraints also limit the possibility of a one loop explanation for the electron anomaly, as discussed later. This leaves us with two scenarios, one where both anomalies are explained via the two loop contribution, following the scaling law in Eq. (45), and another where the muon anomaly is one loop dominated while the electron one is still generated at two loops.

In Fig. 5a the allowed regions for \(\text {Re}\left( n_{\mu }\right) \) are presented as a function of \(m_{\mathrm {H}}\). Three disjoint regions in the scalar mass can be seen: two in the 200–400 GeV range and the other above 1.2 TeV. The low mass regions belong to the scenario where the muon anomaly is obtained through the one loop contribution in agreement with the relation in Eq. (35). Note that this contribution depends on the absolute value of the coupling, so both signs are allowed for \(\text {Re}\left( n_{\mu }\right) \). In the large mass region both leptonic anomalies are two loop dominated.

Figure 5b shows \(m_{\mathrm {H}}\) vs. \(t_\beta \). It contains two separate allowed regions again: in the \(t_\beta \sim 1\) regime only scalar masses above 1.2 TeV are allowed; conversely for \(t_\beta \) larger than 10, \(m_{\mathrm {H}}\) lies in the 200–400 GeV interval.

To complement the previous two plots, in Fig. 5c the relation between the masses \(m_{\mathrm {H}}\) and \(m_{\mathrm {H}^\pm }\) is shown. In the low mass region we can clearly distinguish two scenarios. One where \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\) and another where \(m_{\mathrm {H}^\pm }>m_{\mathrm {H}}\); in the latter, \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\).The degeneracy of \(\mathrm {H}^\pm \) with either \(\mathrm {H}\) or \(\mathrm {A}\) arises from the oblique parameters constraint, as mentioned in Sect. 3. In the large mass region the mass differences do not exceed \(\pm 300\) GeV.

Fig. 6
figure 6

\(\sigma (pp\rightarrow S)_{\mathrm{[ggF]}}\times \text {Br}\left( S\rightarrow \mu ^+\mu ^-\right) \) vs. \(m_S\) allowed regions with \(|n_{\ell }|\le 100\) GeV

Figure 6 illustrates the allowed regions for the resonant process \([pp]_{\mathrm{{ggF}}} \rightarrow S \rightarrow \mu ^{+}\mu ^{-}\) with respect to the scalar mass \(m_S\) for \(S = \mathrm {H}, \mathrm {A}\). The black line corresponds to the limit observed by CMS [140]. Although LHC direct searches are already constraining the allowed regions, there is ample room for extra scalars that can explain both \(g-2\) anomalies simultaneously.

Fig. 7
figure 7

Allowed regions for \(\text {Re}\left( n_{e}\right) \) with \(|n_{\ell }|\le 100\) GeV

Let us now discuss some results concerning \(\text {Re}\left( n_{e}\right) \) and \(\text {Re}\left( n_{\tau }\right) \). With a two loop explanation of the electron anomaly, it follows from Eqs. (38) and (40) (see Appendix B for further details) that one could have expected that both the coupling \(\text {Re}\left( n_{e}\right) \) and the deviation \(\delta a_{e}\) have opposite sign: this is confirmed in Fig. 7 in the 1\(\sigma \) region. However, this figure also contains regions where \(\text {Re}\left( n_{e}\right) \) is negative. This behavior might be understood by analysing with some detail the two loop contribution to \(\delta a_e\): F in Eq. (40) can be decomposed as \(F= F_q+F_\tau +F_\mu \), where \(F_f\) is the contribution with fermion f running in the closed loop. One can estimate the importance of the different contributions for different \(t_\beta \), \(m_{\mathrm {H}}\) and \(m_{\mathrm {A}}\) ranges.

  • For \(t_\beta \sim 1\) and \(m_{\mathrm {H}},m_{\mathrm {A}}>1.2\) TeV, Eq. (40) gives

    $$\begin{aligned}&F_q>0.18,\quad \left| F_\tau \right|<0.13\times 10^{-4}\frac{\text {Re}\left( n_{\tau }\right) }{1\text { GeV}},\nonumber \\&\left| F_\mu \right| <0.12\times 10^{-5}\frac{\text {Re}\left( n_{\mu }\right) }{1\text { GeV}}. \end{aligned}$$
    (48)

    It is clear that in this region the quark-induced contribution \(F_t\) is (i) necessarily dominant and (ii) it requires \(\text {Re}\left( n_{e}\right) \sim 4{-}10\) GeV, as Fig. 7a illustrates, in order to reproduce \(\varDelta _e^{\mathrm{Cs}}\simeq -16\).

  • For \(t_\beta >10\), \(m_{\mathrm {H}}\in [200, 400]\) GeV and \(m_{\mathrm {A}}\in [400,1000]\) GeV, Eq. (40) gives

    $$\begin{aligned}&F_q<0.18,\quad F_\tau \in [0.02, 0.15]\times 10^{-2}\frac{\text {Re}\left( n_{\tau }\right) }{1\text { GeV}}, \nonumber \\&F_\mu \in [0.04, 0.25]\times 10^{-3}\frac{\text {Re}\left( n_{\mu }\right) }{1\text { GeV}}. \end{aligned}$$
    (49)

    In this case, large values of \(\text {Re}\left( n_{\tau }\right) \simeq \pm 100\) GeV give \(\tau \)-induced contributions at the same level of, or even larger than, the quark-induced contribution. This occurs despite some cancellation among the \(\tau \mathrm {H}\) and \(\tau \mathrm {A}\) contributions in Eq. (40). This scenario would require \(\text {Re}\left( n_{e}\right) \lesssim -15\) GeV or \(\text {Re}\left( n_{e}\right) \gtrsim 7\) GeV, as shown in Fig. 7a, in order to reproduce \(\varDelta _e^{\mathrm{Cs}}\simeq -16\).

From this simple estimates one can conclude that, besides the expected regions where \(\delta a_e\) arises from quark-induced two loop contributions, regions where the \(\tau \)-induced contributions have an important role might be present. For this to occur, one might expect some peculiarities: besides light \(\mathrm {H}\) and large \(t_\beta \), large values of both \(|\text {Re}\left( n_{\tau }\right) |\) and \(|\text {Re}\left( n_{e}\right) |\), with \(\text {Re}\left( n_{\tau }\right) \) and \(\text {Re}\left( n_{e}\right) \) having the same sign, are required. Contrary to the case with dominating quark-induced contributions, one might then have allowed regions where \(\text {Re}\left( n_{e}\right) <0\). This is illustrated in Fig. 7a, b where one can observe how allowed \(\text {Re}\left( n_{e}\right) <0\) only appear for a light \(\mathrm {H}\), and how the regions with large \(\pm \text {Re}\left( n_{e}\right) \) correspond to large \(\pm \text {Re}\left( n_{\tau }\right) \).

Fig. 8
figure 8

Allowed regions: relevant correlations involving \(m_{\mathrm {H}}\) with \(|n_{\ell }|\le 250\) GeV

To close this subsection, it is worth analysing in detail the role of the lepton flavor universality constraints mentioned in Sect. 5. As justified later, we focus on observables involving only \(\mu \)’s and e’s. For the ratios

$$\begin{aligned} R_{\mu e}^P=\frac{\varGamma (P^+\rightarrow \mu ^+\nu )}{\varGamma (P^+\rightarrow \mu ^+\nu )_{\mathrm{SM}}}\frac{\varGamma (P^+\rightarrow e^+\nu )_{\mathrm{SM}}}{\varGamma (P^+\rightarrow e^+\nu )}, \end{aligned}$$
(50)

the current constraints are [97]

$$\begin{aligned} R_{\mu e}^\pi= & {} 1+(4.1\pm 3.3)\times 10^{-3},\nonumber \\ R_{\mu e}^K= & {} 1-(4.8\pm 4.7)\times 10^{-3}. \end{aligned}$$
(51)

In the present scenario,

$$\begin{aligned} R_{\mu e}^P=\frac{|1-\varDelta _\mu ^P|^2}{|1-\varDelta _e^P|^2},|1-\varDelta _\ell ^P|^2=\left| 1-\frac{M_P^2}{t_\beta m_{\mathrm {H}^\pm }^2}\frac{\text {Re}\left( n_{\ell }\right) }{m_\ell }\right| ^2, \end{aligned}$$
(52)

and thus, for \(\varDelta _\ell ^P\ll 1\),

$$\begin{aligned} R_{\mu e}^P\simeq 1+2\frac{M_P^2}{t_\beta m_{\mathrm {H}^\pm }^2}\left( \frac{\text {Re}\left( n_{e}\right) }{m_e}-\frac{\text {Re}\left( n_{\mu }\right) }{m_\mu }\right) . \end{aligned}$$
(53)

The presence of \(M_P^2\) and the lepton masses allows us to concentrate on \(R_{\mu e}^K\) and neglect the \(n_{\mu }\) contribution. Therefore from Eq. (51) we get the constraint

$$\begin{aligned} \text {Re}\left( n_{e}\right) < 5\frac{t_\beta m_{\mathrm {H}^\pm }^2}{1\text { TeV}^2}\text { GeV}. \end{aligned}$$
(54)

Then,

  • for \(t_\beta \simeq 1\) and \(m_{\mathrm {H}^\pm }\simeq 2\) TeV, \(\text {Re}\left( n_{e}\right) < 20\) GeV,

  • while for \(t_\beta \simeq 10^2\) and \(m_{\mathrm {H}^\pm }\simeq 0.5\) TeV, \(\text {Re}\left( n_{e}\right) < 125\) GeV.

From muon decay constraints on the \(\mathrm {H}^\pm \) mediated contributions we also have a \(t_\beta \) independent constraint (since the process is purely leptonic) [97, 107]:

$$\begin{aligned} \left| \frac{n_{e} n_{\mu }}{m_{\mathrm {H}^\pm }^2}\right| <0.035. \end{aligned}$$
(55)

This constraint is relevant for the low mass region: for \(\text {Re}\left( n_{\mu }\right) \simeq 100\) GeV, we can rewrite

$$\begin{aligned} |n_{e}|< 87\left( \frac{m_{\mathrm {H}^\pm }}{0.5\text { TeV}}\right) ^2\text { GeV}, \end{aligned}$$
(56)

which is more restrictive than the bound from \(R_{\mu e}^K\) above. Concerning other observables involving \(\tau \) leptons, semileptonic processes are not sensitive to \(n_{\tau }\) due to \(\frac{m_e}{m_\tau }\) and \(\frac{m_\mu }{m_\tau }\) suppressions, while purely leptonic decays have looser bounds than Eq. (55).

This simple numerical exercise confirms that \(\delta a_e^{\mathrm{Exp}}\) cannot be explained through one loop contributions.

Fig. 9
figure 9

Allowed regions for scalar masses with \(|n_{\ell }|\le 250\) GeV

Fig. 10
figure 10

\(\sigma (pp\rightarrow S)_{\mathrm{[ggF]}}\times \text {Br}\left( S\rightarrow \mu ^+\mu ^-\right) \) vs. \(m_S\) allowed regions with \(|n_{\ell }|\le 250\) GeV

6.2 \(|n_{\ell }|\le 250\) GeV

As previously motivated, perturbativity bounds on the Yukawa couplings should be studied in detail. Here we explore higher scales in \(n_{\ell }\), namely changing from \(|n_{\ell }| \le 100\ \mathrm {GeV}\) to \(|n_{\ell }| \le 250\ \mathrm {GeV}\) while maintaining the same constraints of the previous section. Conversely to what one would naively expect, it is not just the allowed regions in the different \(n_{\ell }\) that might change, but it has direct consequences on other physical observables such as the scalar masses and \(t_\beta \), among others.

Figure 8a shows results for \(\text {Re}\left( n_{\mu }\right) \) vs. \(m_{\mathrm {H}}\). It is clear that the allowed regions in parameter space are notably enlarged with respect to those in Fig. 5a, which are completely embedded in the ones of this new analysis, as one could have expected. On that respect, one may realize of the appearance of a new set of intermediate values for the scalar mass, \(m_{\mathrm {H}}\in [0.4; 1.2]\ \mathrm {TeV}\), when increasing our perturbativity upper bound. It can be easily understood by tracing an horizontal line at \(\text {Re}\left( n_{\mu }\right) _0 = -100\) GeV: we eliminate the blue region “bridge” connecting the low and high mass solutions. Therefore, this new range of scalar masses requires large values of \(|\text {Re}\left( n_{\mu }\right) |\).

To fully characterize the impact of perturbativity on the allowed parameter space, Fig. 8b illustrates the scalar mass \(m_{\mathrm {H}}\) in terms of \(t_\beta \). Taking into account the appearance of new intermediate solutions in \(m_{\mathrm {H}}\), one could expect that this behavior is translated into \(t_\beta \). As Fig. 8b corroborates, new values of \(t_\beta \), roughly in the range \(1 \le t_\beta \le 10\), are allowed when changing the perturbativity requirement from \(|n_{\ell }| \le 100\) GeV to \(|n_{\ell }| \le 250\) GeV. Furthermore, one may also notice by comparing with Fig. 5b that the top blue region for large \(m_{\mathrm {H}}\) becomes wider, around a factor 2.5 in \(t_\beta \) for each value of the scalar mass.

Figure 8c shows correlations among the scalar masses \(m_{\mathrm {H}}\) and \(m_{\mathrm {H}^\pm }\). Concerning the low mass regions where \(\mathrm {H}^\pm \) is degenerate either with \(\mathrm {H}\) or \(\mathrm {A}\), already mentioned in Fig. 5c, it can be observed that enlarging perturbativity bounds pushes the upper limit of these regions in such a way that \(m_{\mathrm {H}^\pm }\in [0.2; 1.2]\) TeV for \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\) and \(m_{\mathrm {H}^\pm }\in [0.4; 1.2]\) TeV for \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\), to a high degree of accuracy. Figure 9a, b complete the results for the scalar masses. For instance, it is still true that \(m_{\mathrm {A}}> m_{\mathrm {H}}\) in the low mass region, according to the general constraints presented in Sect. 3.

On the other hand, Fig. 10 shows the resonant process \([pp]_{\mathrm{{ggF}}} \rightarrow S \rightarrow \mu ^{+}\mu ^{-}\) as a function of the scalar mass \(m_S\) for \(S = \mathrm {H}, \mathrm {A}\), which acquire an important role since we may be entering an era of exclusion or discovery at the LHC. As disclosed above, the existence of an intermediate set of solutions, \(m_{\mathrm {H}}\in [0.4; 1.2]\) TeV and \(m_{\mathrm {A}}\in [0.9; 1.2]\) TeV, opens the possibility to detect a sizeable signal in that range of scalar masses that was not contemplated in Fig. 6. Moreover, it is clear that increasing \(|n_{\ell }|\) up to \(|n_{\ell }| \le 250\) GeV modifies our expectations for \(\text {Br}\left( S\rightarrow \ell ^{+}\ell ^{-}\right) \) and, in particular, enlarges the allowed parameter space, as one can easily check.

Fig. 11
figure 11

Allowed regions for \(\text {Re}\left( n_{e}\right) \) with \(|n_{\ell }|\le 250\) GeV

Finally, we should stress some aspects concerning \(\text {Re}\left( n_{e}\right) \) and \(\text {Re}\left( n_{\mu }\right) \) from Fig. 11. In spite of increasing our perturbativity bound up to \(|n_{\ell }| \le 250\) GeV, it still seems difficult to obtain a one loop explanation for the electron anomaly since it requires quite large couplings, namely \(|n_{e}| > 160\) GeV in the Cs case. Figure 11b shows that \(|n_{e}| < 150\) GeV in the relevant range of scalar masses, thus indicating that \(\delta a_{e}^{\mathrm{Exp,Cs}}\) is mainly explained at two loops. This agrees with the discussion on universality constraints closing Sect. 6.1.

As it was already explained in the discussion of Fig. 7a, now in Fig. 11b and for large scalar masses, one can easily check that the electron coupling is positive and lies in the range \(\text {Re}\left( n_{e}\right) \sim 4{-}20\) GeV. Furthermore, according to Eq. (45), there exists a linear relation between \(\text {Re}\left( n_{\mu }\right) \) and \(\text {Re}\left( n_{e}\right) \) for \(m_{\mathrm {H}}> 1.2\) TeV, which implies that they have opposite sign in the Cs case and therefore \(\text {Re}\left( n_{\mu }\right) \) should be negative in this region. The region \(\text {Re}\left( n_{\mu }\right) =-13\text {Re}\left( n_{e}\right) \) can be seen in the lower part of Fig. 11a inside the 1\(\sigma \) region as it should. Departure from this straight line introduces an important one loop contribution to the muon anomaly lowering also the scalar masses ranges. On the other hand, for light scalar masses, \(\text {Re}\left( n_{e}\right) \) might be either positive or negative by the same arguments discussed in Sect. 6.1. It is also important to recall that, in this low mass region, \(\varDelta _\mu \) receives dominant one loop contributions and thus \(\text {Re}\left( n_{\mu }\right) \) could naturally appear with both signs. From Fig. 11a, one may notice as well that \(|n_{\mu }|\) is in general larger than \(|n_{e}|\) in the whole parameter space.

6.3 Different \(\delta a_e\)

As commented in Sect. 5, the situation concerning \(\delta a_e^{\mathrm{{Exp}}}\) is to some extent unclear. In this section we discuss the implications of different assumptions for the value of \(\delta a_e^{\mathrm{Exp}}\), that is, in terms of the model, the implications of requiring different values of the new contributions \(\delta a_e\). The ultimate answer is definitely provided by repeating detailed numerical analyses under the different assumptions \(\delta a_e^{\mathrm{Exp}}\). However, one can anticipate part of the answer with simple considerations. As analysed in Sect. 4, \(\delta a_e\) arises from two loop contributions proportional to \(\text {Re}\left( n_{e}\right) \): this fact, together with the results of Sect. 6.2 corresponding to \(\delta a_e\simeq \delta a_e^{\mathrm{Exp,Cs}}\), can give us a first insight. Consider for example an allowed point in parameter space (i.e. a point respecting all imposed constraints) which gives \(\delta a_e\simeq \delta a_e^{\mathrm{Exp,Cs}}\). This point has a certain \(\text {Re}\left( n_{e}\right) =\text {Re}\left( n_{e}\right) _{\mathrm{Cs}}\); it is straightforward that changing \(\text {Re}\left( n_{e}\right) \mapsto \text {Re}\left( n_{e}\right) '=\text {Re}\left( n_{e}\right) _{\mathrm{Rb}}=\text {Re}\left( n_{e}\right) _{\mathrm{Cs}}\times \frac{\delta a_e^{\mathrm{Exp,Rb}}}{\delta a_e^{\mathrm{Exp,Cs}}}\) and no other parameter, one would obtain \(\delta a_e\simeq \delta a_e^{\mathrm{Exp,Rb}}\). The question is, of course, if such a change in \(\text {Re}\left( n_{e}\right) \) alone still gives an allowed point. On that respect, one needs to analyse which observables constrain \(\text {Re}\left( n_{e}\right) \) and how those constraints work. These are the ones related to lepton flavor universality in leptonic decays \(\ell _i\rightarrow \ell _j\nu {{\bar{\nu }}}\) and in semileptonic decays involving kaons and pions, analysed in Sect. 6.1. In particular, attending to \(\delta a_e^{\mathrm{Exp,Cs}}\), \(\delta a_e^{\mathrm{Exp,Rb}}\) and \(\delta a_e^{\mathrm{Exp,Avg}}\) in Eqs. (2), (3) and (46), one is interested in the effect on those constraints of

$$\begin{aligned} \text {Re}\left( n_{e}\right)&=\text {Re}\left( n_{e}\right) _{\mathrm{Cs}}\mapsto \text {Re}\left( n_{e}\right) '\nonumber \\&=\text {Re}\left( n_{e}\right) _{\mathrm{Rb}}=-0.55\text {Re}\left( n_{e}\right) _{\mathrm{Cs}},\nonumber \\ \text {Re}\left( n_{e}\right)&=\text {Re}\left( n_{e}\right) _{\mathrm{Cs}}\mapsto \text {Re}\left( n_{e}\right) '\nonumber \\&=\text {Re}\left( n_{e}\right) _{\mathrm{Avg}}=0.23\text {Re}\left( n_{e}\right) _{\mathrm{Cs}}, \end{aligned}$$
(57)
Fig. 12
figure 12

\(\text {Re}\left( n_{e}\right) \) correlations

Fig. 13
figure 13

\(\text {Re}\left( n_{\mu }\right) \) vs. \(m_{\mathrm {H}}\) correlations

when no other parameter is changed. There are two different aspects:

  1. 1.

    since \(|\text {Re}\left( n_{e}\right) _{\mathrm{Avg}}|,|\text {Re}\left( n_{e}\right) _{\mathrm{Rb}}|<|\text {Re}\left( n_{e}\right) _{\mathrm{Cs}}|\), the constraint on \(|n_{e}|\) from \(\mu \rightarrow e\nu {{\bar{\nu }}}\) decays in Eq. (56) is necessarily less restrictive when Eq. (57) is considered;

  2. 2.

    besides the uncertainty in \(R^K_{\mu e}\) in Eq. (51), as discussed previously, there is a “sign” question concerning the deviation, at the same \(\simeq 5\times 10^{-3}\) level of the uncertainty, from \(R^K_{\mu e}=1\). In order to obtain \(\delta a_e\simeq \delta a_e^{\mathrm{Exp,Cs}}<0\), the expectation is \(\text {Re}\left( n_{e}\right) >0\), and that produces \(R^K_{\mu e}-1>0\) in Eq. (53), which goes “in the wrong direction”. For both cases in Eq. (57), that problem is alleviated.

It is then clear that the analysis with \(\delta a_e^{\mathrm{Exp,Cs}}\) is somehow a “worst case” scenario in terms of the dependence of the constraints on \(\text {Re}\left( n_{e}\right) \): besides the naive mapping of allowed regions expected from Eq. (57), one might then expect larger allowed regions not only for \(\text {Re}\left( n_{e}\right) \) but also for other quantities of interest. As mentioned in Sect. 5, we also perform an analysis where \(|\delta a_e|\le 20\times 10^{-13}\) is imposed (instead of requiring some specific value, as summarized in Fig. 4). This serves a double purpose: identifying which allowed regions are necessary in order to obtain an appropriate \(\delta a_\mu \) without regard to \(\delta a_e\), and identifying which regions are absolutely excluded for any value of \(\delta a_e\) reasonably compatible with \(\delta a_e^{\mathrm{Exp,Cs}}\) or \(\delta a_e^{\mathrm{Exp,Rb}}\) that one could consider.

Fig. 14
figure 14

\(\delta a_e\) vs. \(\text {Re}\left( n_{e}\right) \)

Fig. 15
figure 15

Results with \((\varDelta S,\varDelta T)\) in Eq. (23), “conservative” case in [103]

In Figs. 1213 and 14, the color coding follows Fig. 4. Figure 12 shows \(\text {Re}\left( n_{\mu }\right) \) vs. \(\text {Re}\left( n_{e}\right) \) and \(\text {Re}\left( n_{e}\right) \) vs. \(m_{\mathrm {H}}\) allowed regions: comparison with Fig. 11a, b confirms the simple expectations of the previous discussion in terms of the position of the allowed regions and their extension. The same applies to Fig. 13, which shows \(\text {Re}\left( n_{\mu }\right) \) vs. \(m_{\mathrm {H}}\) (to be compared with Fig. 8a). In particular it is clear from Fig. 13c that once \(\delta a_\mu \simeq \delta a_\mu ^{\mathrm{Exp}}\) is imposed, the allowed regions for some parameters (besides \(\text {Re}\left( n_{\mu }\right) \), obviously) are coarsely determined and the sensitivity of the analysis on the requirement for \(\delta a_e\) only concerns a finer level of detail.

Fig. 16
figure 16

Results with \((\varDelta S,\varDelta T)\) in Eq. (24), from [104]

There is a final point that the analysis with \(|\delta a_e|\le 20\times 10^{-13}\) confirms. Figure 14 shows \(\delta a_e\) vs. \(\text {Re}\left( n_{e}\right) \): under the simple expectations for the two loop contributions discussed in Sect. 4.2, one would have \(\text {Re}\left( n_{e}\right) \times \delta a_e<0\). Besides that expected region, one can observe smaller allowed regions where \(\text {Re}\left( n_{e}\right) \times \delta a_e>0\): they correspond to the unexpected situation in which the two loop contributions are dominated by virtual \(\tau \)’s in the fermion loop, and furthermore it is clear that the values of \(\delta a_e\) that can be obtained in this manner are more restricted, with \(|\delta a_e|<10^{-12}\).

6.4 The CDF \(M_W\) anomaly

As mentioned in Sect. 3, one can use deviations from the SM in the oblique parameters \((\varDelta S,\varDelta T)\ne (0,0)\) in order to “explain” the CDF measurement of \(M_W\) in [96]: this subsection is devoted to that “explanation”. Figures 15 and 16 show results analogous to the ones in Sect. 6.2 – which use Eq. (22) – except for a different \((\varDelta S,\varDelta T)\) constraint. Figure 15 is obtained with Eq. (23) (the “conservative” average of [103]) and Fig. 16 is obtained with Eq. (24) (the results in [104]). The coloring of the allowed regions corresponds, darker to lighter, to \(1,2,3\sigma \) levels of a 2D-\(\varDelta \chi ^2\). For \(\varDelta \chi ^2\) we use the \(\chi ^2_{\mathrm{Min}}\) value of the analysis in Sect. 6.2 (that is, with the constraint in Eq. (22) for \(\varDelta S\), \(\varDelta T\)). A few comments are in order.

  • Besides the absence of degeneracies \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\) or \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\), masses of the new scalars larger than 2 TeV are more difficult to accommodate. This can be understood attending to the clash between the mass differences discussed in Sect. 3 that Eqs. (23) or (24) require, and the need of near degenerate scalars that the perturbativity requirements on the scalar potential impose for new scalar masses much larger than \(v_{}\).

  • Overall agreement with the imposed constraints is worse in several regions in Figs. 15 and 16 than it was in the analyses of Sect. 6.2 (Figs. 89 and 11). This is more dramatic in Fig. 16, where the agreement with constraints is worse than in Fig. 15 to the point that several regions are beyond the represented contour levels.

Despite these changes, the main characteristics of the allowed regions discussed in the previous sections still apply and are clearly identified in both Figs. 15 and 16.

Fig. 17
figure 17

\(\varDelta S\) vs. \(\varDelta T\), 1, 2 and 3\(\sigma \) 2D-\(\varDelta \chi ^2\) contours from the imposed constraint are shown together with the allowed regions

Table 1 Example points, masses and \(\text {Re}\left( n_{\ell }\right) \)’s in GeV
Table 2 Example points, \(\delta a_e\) values; columns 3 to 9 show the relative contributions of the different one and two loop terms to the value of \(\delta a_e\) in the second column

Finally, since the oblique parameters S and T play an important role, Fig. 17 shows allowed regions for \(\varDelta S\) vs. \(\varDelta T\) in the two scenarios considered for the CDF \(M_W\) “explanation”, together with the imposed \((\varDelta T,\varDelta S)\) constraint in each case. As anticipated, the constraint in Eq. (24) appears to be more difficult to accommodate than the constraint in Eq. (23). In fact, despite the different position of the ellipses corresponding to the \((\varDelta T,\varDelta S)\) constraints in Fig. 17a, b, the allowed regions are quite similar in both cases, that is, the model appears to be unable to accommodate values \(\varDelta T>0.22\) together with \(\varDelta S>0.02\). Other possible explanations of the CDF \(M_W\) anomaly have been addressed in [148,149,150,151,152,153,154,155].

Table 3 Example points, \(\delta a_\mu \) values; columns 3 to 9 show the relative contributions of the different one and two loop terms to the value of \(\delta a_\mu \) in the second column

6.5 Example points

In this section, some example points of the allowed parameter space are presented in order to specify the behavior pointed out in the previous plots. For the sake of clarity, we only focus on the analysis with “\(a_e^{\mathrm{Cs}}\)” concerning the electron anomaly: other cases do not change substantially beyond the differences already mentioned in Sect. 6.3.

From Table 1, it is clear that points 1–2 correspond to the solution with small values of \(t_\beta \) and large scalar masses: all scalars are above \(1.2\ \mathrm {TeV}\) and their mass differences do not exceed \(\pm 200\) GeV. In this region, both anomalies are explained at two loops through the top quark terms, as one can easily check in Tables 2 and 3, where all loop contributions are normalized to the total \(\delta a_{\ell }\) in such a way that their sum must be 1. One may also notice that \(\delta a_{\mu }\) receives a subdominant one loop contribution. The lepton couplings \(\text {Re}\left( n_{e}\right) \) and \(\text {Re}\left( n_{\mu }\right) \) have opposite sign and they roughly satisfy the linear relation in Eq. (45) for the Cs case.

Regarding the appearance of the intermediate values of the scalar masses and \(t_\beta \) previously commented in Sect. 6.2, our point 3 gives a perfect example of that behavior. It is important to realize that large values of \(|\text {Re}\left( n_{\mu }\right) |\) are required in this region; in fact, they are almost reaching the perturbativity upper bound \(|n_{\ell }| \le 250\) GeV. On the other hand, although the top dominance still holds at two loops in the electron anomaly, the corresponding tau contributions begin to play a relevant role. This trend will continue as \(t_\beta \) grows and the quark contributions are more suppressed.

Finally, points 4–9 belong to the low mass region corresponding to a wide range of \(t_\beta \gg 1\) values. As we have stressed before, two possible scenarios arise: one where \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\) (points 4–6) and another where \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\) (points 7–9). In all cases, the scalar masses are below 1 TeV and \(m_{\mathrm {A}}> m_{\mathrm {H}}\), as anticipated. Taking into account the large values of \(t_\beta \), the two loop contribution that dominates \(\delta a_e\) is generated by the tau loop. This confirms our expectation for \(\text {Re}\left( n_{e}\right) \): its sign is not fixed and it could be either positive or negative (point 9). Furthermore, in this region the muon anomaly is clearly one loop dominated, albeit there exists a subdominant contribution from the tau loop as well. This in turn means that \(\text {Re}\left( n_{\mu }\right) \) can take both signs, as one can easily check.

For completeness, the last two points have been included to give an example of the allowed parameter space in Sect. 6.4 considering the CDF \(M_W\) anomaly. It is clear that point 10 mimics the behavior of points 1–2, while point 11 presents the same features as points 4–6.

7 Conclusions

The experimental determinations of the muon and the electron anomalous magnetic moment point towards the necessity of lepton flavor non-universal New Physics. Aiming to address both leptonic anomalies simultaneously, we have considered a type I or type X 2HDM with a general flavor conserving lepton sector, one loop stable under renormalization, in which the new Yukawa couplings are completely decoupled from lepton mass proportionality. The latter turns out to be crucial in order to reproduce the \(g-2\) muon anomaly together with the different scenarios one can consider for the \(g-2\) electron anomaly, related to the Cs and/or the Rb recoil measurements of the fine structure constant. A thorough analysis of the parameter space of the model has been performed including all relevant theoretical and experimental constraints. The results show that the muon anomaly receives dominant one loop contributions for light new scalar masses in the 0.2–1.0 TeV range together with a significant hierarchy in the vacuum expectation values of the scalars, that is \(t_\beta \gg 1\), while two loop Barr–Zee diagrams are also needed for heavy new scalars with masses above 1.2 TeV together with \(t_\beta \sim 1\). On the other hand, the electron anomaly receives dominant two loop contributions in the whole range of scalar masses. Furthermore, we have analysed how the perturbativity assumptions on the lepton Yukawa couplings have direct impact on relevant physical observables: intermediate values of the scalar masses and \(t_\beta \) only arise when the perturbativity upper bound on \(n_{\ell }\) reaches the electroweak scale. This might be relevant since we are entering an era of exclusion or discovery at the LHC, so that the allowed parameter space of the model must be fully scrutinized. The disagreement between the recent CDF measurement of \(M_W\) and the SM expectations for electroweak precision results can be translated into deviations \((\varDelta S,\varDelta T)\ne (0,0)\) of the oblique parameters. We have considered two different scenarios for \((\varDelta S,\varDelta T)\) values which “explain” the CDF disagreement. Both scenarios require a scalar spectrum where near degeneracies \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {H}}\) or \(m_{\mathrm {H}^\pm }\simeq m_{\mathrm {A}}\) are now disfavored, and where masses larger than 2 TeV are more difficult to accommodate. However, concerning the \(n_{\ell }\) couplings and \(t_\beta \), the allowed regions have the same characteristics as in the analyses compatible with \((\varDelta S,\varDelta T)=(0,0)\).