Impact of Fermionic Electroweak Multiplet Dark Matter on Vacuum Stability with One-loop Matching

We investigate the effect of fermionic electroweak multiplet dark matter models on the stability of the electroweak vacuum using two-loop renormalization group equations (RGEs) and one-loop matching conditions. Such a treatment is crucial to obtain reliable conclusions, compared with one-loop RGEs and tree-level matching conditions. In addition, we find that the requirement of perturbativity up to the Planck scale would give strong and almost mass-independent constraints on the Yukawa couplings in the dark sector. We also evaluate these models via the idea of finite naturalness for the Higgs mass fine-tuning issue.


I. INTRODUCTION
The discovery of the Higgs boson in 2012 at the Large Hadron Collider (LHC) [1,2] confirms the particle content of the standard model (SM) and the validity of the Brout-Englert-Higgs mechanism. With present experimental values of SM parameters, if the SM is valid up to the Planck scale (1.2 × 10 19 GeV), a deeper minimum would appear at ∼ 10 17 GeV, indicating the electroweak (EW) vacuum is metastable [3][4][5][6]. Nevertheless, such a situation could be modified by new physics beyond the standard model (BSM) at scales lower than the Planck scale. Therefore, the requirement of a stable or metastable EW vacuum may strongly constrain BSM models.
On the other hand, the existence of dark matter (DM) has been established by solid astrophysical and cosmological observations. This undoubtedly indicates that BSM physics must exist. Among various DM candidates proposed, weakly interacting massive particles (WIMPs) are very compelling and have been widely studied. WIMP models can be easily constructed by extending the SM with new electroweak multiplets, such as minimal dark matter models [7][8][9][10][11][12][13][14], and other models which contain more than one SU(2) L multiplets . In this paper, we focus on a class of fermionic electroweak multiplet dark matter (FEMDM) models which involve a dark sector with more than one fermionic SU(2) L multiplets.
Specifically, we take the following three models as illuminating examples: • Singlet-doublet fermionic dark matter (SDFDM) model: the dark sector contains one singlet Weyl spinor and two doublet Weyl spinors; • Doublet-triplet fermionic dark matter (DTFDM) model: the dark sector contains one triplet Weyl spinor and two doublet Weyl spinors; • Triplet-quadruplet fermionic dark matter (TQFDM) model: the dark sector contains one triplet Weyl spinor and two quadruplet Weyl spinors; After the electroweak symmetry breaking (EWSB), these multiplets can mix with each other through Yukawa couplings, and the mass eigenstates include neutral Majorana fermions χ 0 i , singly charged fermions χ ± i , and (if in the TQFDM model) a doubly charged fermion χ ±± . By imposing a discrete Z 2 symmetry, the lightest neutral fermion χ 0 1 is stable, serving as a DM candidate. All these additional EW multiplets can alter the high energy behaviors of the running couplings through the renormalization group equations (RGEs). For instance, given the contributions of the new physical states, the quartic couplings λ in the Higgs potential might stay positive up to the Planck scale. In this case, these models render a stable EW vacuum up to the Planck scale. However, there could be an opposite effect if the new couplings are too large. In such a case, Landau poles would appear and render the breakdown of the theory. Thus, by investigating the conditions for EW vacuum stability and Landau poles, the parameter spaces of these FEMDM models could be constrained.
In addition, the mass term of the Higgs doublet could receive loop corrections from Yukawa couplings with the new states. These corrections are generally proportional to mass squares of the new particles, so they would give rise to a naturalness problem if the new particles are too heavy. In this paper we will adopt an idea called finite naturalness [42] to evaluate such a effect.
Note that some papers in the literature have studied the above effects utilizing one-loop or two-loop RGEs, but they concentrate on different models, such as singlet extensions [4,40,[43][44][45][46], triplet extensions [3,6], two Higgs doublet models [47][48][49][50][51][52], and so on. Furthermore, when it comes to the initial values of running parameters, only the tree level matching is considered in these works. Nonetheless, it is well known that the quartic coupling λ almost vanishes at high energy scales in the SM, and hence the next-to-next-to-leading-order (NNLO) corrections to λ are important to determine the fate of the EW vacuum [5]. Therefore, the tree-level matching seems not sufficient to accurately determine the initial values of running parameters when considering corrections to λ from new physics [53,54]. The loop contributions from the dark sector deserve accurate calculations for studying the vacuum stability problem. In our calculations, we utilize three-loop RGEs and two-loop matching for the SM sector [5], as well as two-loop RGEs and one-loop matching for the dark sector. This paper is outlined as follows. In Sec. II we give a brief introduction of our strategy for the one-loop matching of the dark sector. In Sec. III we study the RGE running of dimensionless couplings in the SDFDM model and the effects on the EW vacuum. Using the same method, the results for the DTFDM and TQFDM models are demonstrated in Sec. IV and Sec. V. The conclusions are given in Sec. VI. Appendix A gives the explicit expressions for two-loop RGEs in the FEMDM models.

II. MATCHING AND RUNNING
To study the evolution of a theory from a low energy scale to a high energy scale, two ingredients are necessary: • The RGEs of all the running parameters.
• The initial values of these parameters at the low energy scale where the evolution starts.
The first ingredient involves the calculations of β-functions for the given theory; the second ingredient concerns the matching conditions between the running parameters and observables. In this paper, we always carry out the loop calculations in the MS scheme, because in this scheme all the parameters have gauge-invariant RGEs [5,55].
The β-function describing the evolution of a given parameter ξ can be defined as where µ is the energy scale. By expanding β(ξ) in a perturbative series, we have where M W 0 is the bare mass of the W boson, A W W is the W self-energy, V W is the vertex contribution in the muon decay process, B W is the box contribution, E is a term due to the renormalization of external legs, and M is a mixed contribution due to a product of different objects among V W , A W W , B W , and E. All quantities in Eq. (15) are computed at zero external momentum. Thus the one-loop term of ∆r 0 is given by where we have used M (1) = 0 [5] to get ride of the last term in Eq. (15). With Eq. (5), we can get the one-loop relations between MS parameters and physical observables as follows: Note that here we only list the most important results for the one-loop matching, and we will use these formulas to calculate the one-loop corrections of the FEMDM models to the SM fundamental parameters. More detailed derivations and even two-loop matching conditions can be found in Refs. [5,67,68].

A. Model Details
In the SDFDM model [15-18, 20, 23, 27-30, 36], the dark sector involves colorless Weyl fermions S, D 1 , and D 2 obeying the following (SU(2) L , U(1) Y ) gauge transformations: The signs of the hypercharges of the two doublets are opposite, making sure that the model is anomaly free. The gauge-invariant Lagrangian is given by where L SM is the SM Lagrangian and In Eq. (21), m S and m D are the masses of the singlet and doublets, y 1 and y 2 are Yukawa couplings in the dark sector, and D µ = ∂ µ − igW a µ τ a − ig Y B µ is the covariant derivative with τ a being the generators for the corresponding SU(2) L representations. H is the Higgs doublet andH ≡ iσ 2 H * . The constants b ij , c ij , and d ij render the gauge invariance of the mass and Yukawa terms. They can be decoded from Clebsch-Gordan (CG) coefficients multiplied by an arbitrary normalization factor. The nonzero values for them are After the EWSB, the Higgs doublet obtains a vacuum expectation value v, resulting in mixing between the singlet and doublets because of the Yukawa coupling terms. It is instructive to rewrite the Lagrangian with mass eigenstates instead of gauge eigenstates, and reform the interaction terms with four-component spinors. Nonetheless, in this paper we will not show all these details, which can be found in Ref. [36].

B. RGEs and One-loop MS Parameters
The evolution of various dimensionless couplings in the SDFDM model with an energy scale µ is performed in this subsection. We use PyR@TE 2 [60] to carry out auxiliary calculations, and obtain RGEs at the two-loop level. The full β-functions consist of two parts: where β SM is the contribution of the SM sector, and β SD is the contribution of the dark sector. β SM of SM couplings at the three-loop level in the MS scheme can be found in Ref. [5]. β SD at the two-loop level are listed in Appendix A 1.
In Sec. II we have shown how to calculate one-loop MS parameters, and now we can use Eqs. (11)- (18) to calculate the initial values of the running couplings. Instead of demonstrating all the details of the loop calculations, here we give a brief description, and show some important results for some benchmark points (BMPs).
In Eqs. (11)- (18) we can see that there are several quantities needed for a given theory: tadpole diagrams of the Higgs boson iT , mass counterterms of the W , Z, and Higgs bosons, and ∆r 0 . Here we denote the mass counterterms as All these self-energies are computed with on-shell external particles, and the corresponding Feynman diagrams can be found in Fig. 1. Initial values for RGE running All the loop diagrams can be expressed via Passarino-Veltman (PV) functions, whose numerical values are obtained using LoopTools 2.13 [69]. All the calculations are stick to the MS scheme with the renormalization scale setting at M t . In Table II, we present the values of the fundamental parameters for SM and SM+SDFDM at µ = M t . Here we adopt two BMPs for parameters in the SDFDM model: BMP1 with y 1 = y 2 = 0.1, m S = m D = 1000 GeV, and BMP2 with y 1 = y 2 = 0.6, m S = m D = 1000 GeV. Yukawa couplings in BMP2 are larger than those in BMP1.
With the β-functions listed in the Appendix A 1, we solve the RGEs to give the evolution of these parameters along with the energy scale µ. In Fig. 2 we demonstrate the evolution of the Higgs quartic coupling λ (solid blue line) and the Yukawa couplings y 1 (solid green line) for the two BMPs. As the one-loop β-functions of λ and y 1 can sufficient to give a better understanding of the running behaviors, we write down their expression here: In Fig. 2(a) for BMP1, we find that the running value of y 1 is almost invariant up to the Planck Differences between tree-level and one-loop matching λ| µ=Mt λ min log 10 (µ/GeV)| λ=0 log 10 (µ/GeV)| λ=λ min log 10 (P 0 ) Differences between tree-level (2,T) and one-loop (2,1) matching for SM+SDFDM with BMP1 (y 1 = y 2 = 0.1, m S = m D = 1000 GeV) are presented. The second column list values in the SM. λ min is the minimal value of the running λ. P 0 in the last column represents the EW vacuum decay probability, which is discussed in Subsection III C.
scale. This is because its β-function is proportional to its value, and when y 1 and y 2 are small the running gets suppressed. In addition, because of small y 1 and y 2 , the contributions of the dark sector to the running of λ are also insignificant. Thus the running values of λ for SM (solid red line) and SM+SDFDM (solid blue line) are very close, and even coincide with each other at low energy scales. In addition, the differences between tree-level (dashed blue line) and one-loop (solid blue line) matching is also demonstrated. More accurate numerical analysis of such differences are given in Table III. We find that the matching conditions are important for the minimal running value λ min , and this would obviously influence on the decay probability of the EW vacuum.
demand that the theory remains perturbative up to the Planck scale, leading to strong constraints on the Yukawa couplings y 1 and y 2 .
In the above calculation, we have set the matching scale at Q = M t . In addition, we would like to investigate the effects of the matching scale Q on the evolution of λ. This can be realized through the following steps: firstly we carry out the coupling running with SM RGEs to the energy scale Q, and then perform the one-loop matching and carry out the running with SM+SDFDM RGEs up to higher scales. In Table IV we have compared the differences in the evolution of λ for BMP1 (y 1 = y 2 = 0.1, m S = m D = 1000 GeV) with Q setting at M t , 300 GeV, 500 GeV, and 1000 GeV, respectively. The values of λ at 10 3 GeV, 10 5 GeV, 10 10 GeV, and 10 16.5 GeV for each value of Q are presented. For a low energy scale, say, µ = 10 3 GeV, we find that the effects of Q varying from M t to 10 3 GeV are at the order of 0.1%, but such differences will be amplified when λ runs to high energy scales. For instance, at µ = 10 16.5 GeV we find that the difference is approximately the same level with that caused by the choice of tree-level or one-loop matching. Note that the matching scale dependence will decrease if higher order calculation and matching get involved [54]. In the following analysis, we will just set Q = M t .
There is another conception called finite naturalness [42] we would like to introduce for evaluating the FEMDM models. The idea is that we should ignore uncomputable quadratic divergences, so that the Higgs mass is naturally small as long as there are no heavier particles that give large finite contributions to the Higgs mass. In this sense, the fine-tuning ∆ at one-loop level in the MS scheme can be defined as and a smaller ∆ means that the theory is more natural. For example, in the SM we have is the β-function of the quadratic coefficient of the Higgs potential m 2 at one-loop level [42]. According to Eq. (30), we can get ∆ = 0.133 for the SM with the renormalization scale setting at M t . This is a small number, and thus we can say that the SM satisfies finite naturalness.
In Fig. 3 we demonstrate the contours of ∆(µ = M t ) at the m S -m D plane for two sets of parameters: y 1 = 0.1, y 2 = 0.2 (red lines) and y 1 = 0.4, y 2 = 0.2 (blue lines). For each parameter set, the contours that corresponding to ∆ = −1 (dashed lines) and ∆ = −10 (solid lines) are presented. We find that the contours for fixed ∆ shrink as the Yukawa couplings increase, because the dark sector contributions to the Higgs mass correction δ (1) m 2 OS also increase. Therefore, if we demand a small fine-tuning, say, |∆| < 1, there will be upper bounds for the masses of dark sector particles.

C. Tunneling Probability and Phase Diagrams
The present experimental values of M h and M t indicate that the SM EW vacuum might be a false vacuum, which may decay to the true vacuum through quantum tunneling. The present vacuum decay probability P 0 is expressed as [5,70,71] where H 0 = 67.4 km sec −1 Mpc is the present Hubble rate, and S(Λ B ) is the action of bounce of size R = Λ −1 B , given by In practice Λ B is roughly determined by the condition β λ (Λ B ) = 0, and at that energy scale the Higgs quartic coupling λ achieves its minimum value λ min . Note that if Λ B > M Pl , we can only get , nor gravitational corrections to the tunneling rate [72].
Using the initial parameter values given in Sec. II, we obtain λ(Λ B ) and its corresponding energy scale Λ B through analyzing the evolution of λ. Therefore, with Eqs. (31) and (32) we can calculate the tunneling probability of the EW vacuum. In this paper, we classify different states of the EW vacuum using the following conventions.
• Stable: λ > 0 for µ < M Pl ; • Metastable: λ(Λ B ) < 0 and P 0 < 1; • Unstable: λ(Λ B ) < 0 and P 0 > 1; • Non-perturbative: |λ| > 4π before the Planck scale 1 . Now we can discuss the status of the EW vacuum in the SDFDM model. In Fig. 4(a), we display the non-perturbative region (red), the unstable region (yellow), and the metastable region (green) in the y 1 -y 2 plane with fixed mass parameters of m S = m D = 1000 GeV. A similar plot in the Y -M plane is demonstrated in Fig. 4(b), where Y = y 1 = y 2 and M = m S = m D . We can see that the non-perturbative region is almost independent of the masses of the singlet and doublets, while the metastable region shows a little dependency on these mass parameters. There are two reasons for this: (1) the β-functions in the MS scheme is mass-independent, meaning that the mass parameter do not directly enter the β-functions; (2) the mass parameters can only affect the initial values of running parameters through the loop matching conditions, but such a effect is only at one percent level or even smaller (see Table II and Fig. 5).
The effect of such tiny contributions on the evolution of the Yukawa couplings y 1 and y 2 , or on the position of the Landau poles, is negligible, but it indeed affects the decay rate of the EW vacuum (see Table. III). We conclude that the requirement of perturbativity gives a strong and almost mass-independent constraint on the Yukawa couplings, roughly |Y | 0.5, and an even stronger upper bound, roughly |Y | 0.3, can be obtained by demanding the metastability of EW vacuum.
Moreover, in Fig. 6 we present the phase diagram in the M h -M t plane for SM+SDFDM with two BMPs. Similar diagrams for SM are give in Refs. [5,67], where the SM vacuum stability up to the Planck scale is excluded at 2.8σ. In Fig. 6(a), we set y 1 = 0.1, y 2 = 0.2, m S = 300 GeV, m D = 500 GeV, and find that the vacuum stability up to the Planck scale in the SDFDM model is excluded at ∼ 2.0σ. Therefore, with relatively small y 1 and y 2 the EW vacuum is more stable than the SM case. This can also be concluded from the last column of Table III. However, in Fig. 6(b) with y 1 = 0.4, y 2 = 0.2, m S = m D = 500 GeV, we can see that the EW vacuum is more likely unstable for large y 1 and y 2 . Note that these results can also be further understood through Fig. 4(a), as the status of the EW vacuum is almost independent of the mass parameters in the dark sector. In Fig. 4(a), the point (y 1 , y 2 ) = (0.1, 0.2) locates in the metastable region, while the point (y 1 , y 2 ) = (0.4, 0.2) locates at the junction of the metastable and unstable regions.

IV. DTFDM MODEL
In the dark sector of the DTFDM model [21,29,36,37], we introduce two SU(2) L Weyl doublets and one SU(2) L Weyl triplet, which obey the following (SU(2) L , U(1) Y ) gauge transformations: The gauge invariant Lagrangian is with There are four independent parameters: m T , m D , y 1 , and y 2 . After the EWSB, the triplet and doublets mix with each other. More details can be found in Ref. [36]. The two-loop β-functions in the MS scheme are listed in Appendix A 2. Using the same strategy as in Sec. III we can study the influences of the dark sector on the stability of the EW vacuum. Analogue results are presented in Fig. 7. Comparing with Fig. 4 in the SDFDM model, there are two obvious differences: (1) the contours in Fig. 7(a) are rotated by roughly 45 degrees compared to the SDFDM case; (2) for the DTFDM model, the absolute stable region appears, and in Fig. 7(b) the metastable region have slightly larger dependence on the mass parameters of the multiplets.
The first difference is caused by the β-functions of Yukawa couplings y 1 and y 2 , and here we list the one-loop β-function for y 1 : Comparing with Eq. (28), we find that the contributions of the terms proportional to y 2 2 are smaller, and hence the evolution of y 1 is basically independent of y 2 .
The second difference is because the loop corrections to the initial parameter values in the DTFDM model are larger than those in the SDFDM model. This can be observed in Fig. 8(a), where the relative corrections of y t and λ in the m D -m T plane with fixed y 1 = 0.1, y 2 = 0.2 are presented. In Fig. 8(b), contours of the fine-tuning ∆ are demonstrated for evaluating finite naturalness.

V. TQFDM MODEL
In the TQFDM model [26,29,38], the dark sector involves one colorless left-handed Weyl triplet T and two colorless left-handed Weyl quadruplets Q 1 and Q 2 , obeying the following (SU(2) L , U(1) Y gauge transformations: We have the following gauge invariant Lagrangian: where The constants a ij , b ij , c ijk , and d ijk can be decoded from CG coefficients multiplied by a normalizing factor. Here we list their nonzero values: After the EWSB, the triplet will mix with the quadruplets due to the Yukawa coupling terms. More details can be found in our previous work [38]. The two-loop β-functions in the MS scheme are listed in Appendix A 3.
We investigate the influences of the TQFDM model on the EW vacuum, and find that the results are quite different from the previous two similar models. This is mainly due to the RGE of g 2 . For this kind of FEMDM models, the general one-loop β-function for g 2 can be written as where is the contribution from the SM sector, n f is the number of multiplets that transform under the r-dimensional representation of SU(2) L , and C(r) is the corresponding Dynkin index which defined as tr(t a r t b r ) = C(r)δ ab . An additional factor 1/2 is added because the multiplets in these models are Weyl spinors. With this formula we can write down the dark sector contributions in the three models: In the SDFDM and DTFDM models, β total (g 2 ) is negative, and g 2 becomes smaller as the energy scale µ goes up. In contrary, β total (g 2 ) is positive in the TQFDM model, and hence g 2 becomes larger and larger as µ increases and could hit a Landau pole at some high energy scale.
In Fig. 9(a) we show the evolution of y 1 , λ, and g 2 based on one-loop (dashed lines) and two-loop (solid lines) RGEs with fixed y 1 = y 2 = 0.5 and m T = m Q = 1000 GeV. According to our criteria for the status of the EW vacuum in Sec. III C, we get completely different conclusions for one-loop and two-loop RGE analyses: (1) based on one-loop RGEs, the EW vacuum is stable, because λ min > 0 and all parameters are perturbative; (2) based on two-loop RGEs, the theory breaks down, because the Landau pole of g 2 appears at around 10 17 GeV.
In Fig. 9(b) we demonstrate the status of the EW vacuum in the Y -M plane with Y = y 1 = y 2 and M = m T = m Q , based on two-loop RGEs in the TQFDM model. We find that most of the parameter space is excluded because of the non-perturbativity of λ. The solid lines indicate the energy scale where λ becomes non-perturbative. Among them, blue and green lines mean the non-perturbativity is caused by y 1 (y 2 ) and g 2 , respectively. The allowed regions are concentrated in a region with large masses and small Yukawa couplings. From Fig. 10(a), we can see that in this region g 2 gets a large negative corrections, which delays the appearance of the Landau pole. Nevertheless, in this region the issue of finite naturalness becomes prominent, as shown in Fig. 10(b).

VI. CONCLUSIONS AND DISCUSSIONS
In this paper, we have investigated the high energy behavior of the FEMDM models and their impacts on the stability of the EW vacuum. The calculations for the dark sector are carried out in the MS scheme, based on one-loop matching at the M t scale and two-loop RGEs. Differences between tree-level and one-loop matching, and between one-loop and two-loop RGEs are also compared. In addition, we have studied the effects of different matching scales on our results. We have found that the requirement of a stable (or metastable) EW vacuum and perturbativity would give strong constraints on the parameter space. Besides, the idea of finite naturalness is important for evaluating these FEMDM models.
For the SDFDM model, we have discussed the rationality and necessity of one-loop matching, and we have found that this had a significant effect on the stability of the EW vacuum, deserving careful calculations. Moreover, we have found that the effects of the matching scale Q varying from M t to 1000 GeV on the evolution of λ is small ( 2%) for µ < 10 5 GeV. Such a Q-dependence is expected to decrease if higher loop matching and RGEs are taken into consideration. The requirement of perturbativity gives a strong and almost mass-independent constraints on Yukawa couplings (|y 1 |, |y 2 | 0.5), and the constraints from metastability are even more stronger (|y 1 |, |y 2 | 0.3). Two phase diagrams with fixed Yukawa and mass parameters are exhibited to give a clear indication of the effects on the EW vacuum stability. In addition, the contours of the fine-tuning ∆ for evaluating finite naturalness are demonstrated, and if we let ∆ vary from −1 to −10, the corresponding upper bounds on mass parameters will vary from 1 − 2 TeV to 3 − 5 TeV. For the DTFDM model, general conclusions are similar to the SDFDM case. Nonetheless, because of more new states are introduced, the corrections to the initial values of running couplings are larger than those in the SDFDM model. As a result, the EW vacuum can be absolute stable for some parameter regions. Moreover, the fine-tuning ∆ receives more corrections and becomes worse.
For the TQFDM model, the situations are quite different from the other two models because of the opposite evolution behavior of g 2 . This could lead to a Landau pole of λ before the Planck scale. We have found that the requirement of perturbativity can exclude most of the parameter space, and the allowed regions are concentrated in a region with |y 1 |, |y 2 0.8 and m T , m Q 5 TeV, where the issue of finite naturalness becomes prominent. Furthermore, if EW fermionic multiplets in higher dimensional representations are introduced, the situation for the evolution of g 2 would become worse, rendering the breakdown of the theory at even lower energy scales. In addition, we have emphasized the importance of two-loop RGEs by comparing the results with those based on one-loop RGEs. We have found that the conclusions inferred from one-loop and two-loop RGE analyses can be completely different.

SDFDM Model
The contribution to the β-functions up to two-loop level in the SDFDM model are presented as follows.

DTFDM Model
The contribution to the β-functions up to two-loop level in the DTFDM model are listed as follows.

TQFDM Model
The contribution to the β-functions up to two-loop level in the TQFDM model are presented as follows.