Maximum value of the spin-independent cross section in the THDM+a

We investigate the maximum value of the spin-independent cross section ($\sigma_\text{SI}$) in a dark matter (DM) model called the two-Higgs doublet model + a (THDM+a). This model can explain the measured value of the DM energy density by the freeze-out mechanism. Also, $\sigma_\text{SI}$ is suppressed by the momentum transfer at the tree level, and loop diagrams give the leading contribution to it. The model prediction of $\sigma_\text{SI}$ highly depends on values of $c_1$ and $c_2$ that are the quartic couplings between the gauge singlet CP-odd state ($a_0$) and Higgs doublet fields ($H_1$ and $H_2$), $c_1 a_0^2 H_1^\dagger H_1$ and $c_2 a_0^2 H_2^\dagger H_2$. We discuss the upper and lower bounds on $c_1$ and $c_2$ by studying the stability of the electroweak vacuum, the condition for the potential bounded from the below, and the perturbative unitarity. We find that the condition for the stability of the electroweak vacuum gives upper bounds on $c_1$ and $c_2$. The condition for the potential to be bounded from below gives lower bounds on $c_1$ and $c_2$. It also constrains the mixing angle between the two CP-odd states. The perturbative unitarity bound gives the upper bound on the Yukawa coupling between the dark matter and $a_0$ and the quartic coupling of $a_0$. Under these theoretical constraints, we find that the maximum value of the $\sigma_\text{SI}$ is $\sim 5\times 10^{-47}$ cm$^2$ for $m_A = $ 600 GeV, and the LZ and XENONnT experiments can see the DM signal predicted in this model near future.


Introduction
One of the great achievements in Cosmology is the precise determination of the energy density of the dark matter (DM) by the Planck collaboration, Ωh 2 = 0.120 ± 0.001 [1]. The measured value is explained successfully by DM models that use the freeze-out mechanism [2], which have been widely studied for a long time. Those models generally predict non-zero DM-nucleon scattering cross section and have been searched by the direct detection experiments, such as the Xenon1T experiment [3]. However, no significant signals have been observed until now, and the null results set upper bound on the DM-nucleon scattering cross section. The latest result by the Xenon1T experiment gives a severe constraint on DM models.
If a DM particle is a gauge singlet fermion, χ, and couples to a scalar mediator, a 0 , with pseudo-scalar type interaction,χiγ 5 χa 0 , then it is possible to avoid this strong constraint from the Xenon1T experiment while keeping the success of the freeze-out mechanism [4,5]. The two-Higgs doublet model + a (2HDM+a) [6] is one of the models that realize this idea. 1 In addition to the introduction of the DM and the mediator, the Higgs sector is extended into the two-Higgs doublet model. The CP invariance is assumed in the dark sector and the scalar sector. Then, the dark sector and the visible sector can interact through the mixing between a 0 and the CP-odd scalar (A 0 ) in the two-Higgs doublet sector. The model predicts rich phenomenology [10][11][12][13][14][15][16][17], and it is summarized in Ref. [18].
The 2HDM+a predicts non-zero spin-independent DM-nucleon scattering cross section (σ SI ) at loop level [6,14,16,17,19]. In particular, it was shown that if c 2 , which is a quartic coupling between a 0 and a Higgs doublet field H 2 , is large enough, the model can be tested at the forthcoming direct detection experiments [19]. However, such a large coupling causes theoretical problems. If the coupling takes large negative value, the potential can be unbounded from the below. If the coupling is very large, it hits a Landau pole near the electroweak scale and the model loses predictability.
In this paper, we study the constraint on the scalar potential from the boundedness of the scalar potential, the stability of the electroweak vacuum, and perturbative unitarity. Using these constraints, we investigate the upper and the lower bounds on the scalar quartic couplings, and discuss the maximum value of σ SI . We show that the maximum value of σ SI is below the current constraint from the Xenon1T experiment and above the prospect of the LZ experiment [20] and the XENONnT experiment [21].
The rest of the paper is organized as follows. In Sec. 2, we briefly describe the 2HDM+a.
In Sec. 3, we investigate theoretical constraints on the model parameters. Conditions for the electroweak vacuum as the global minimum of the scalar potential, for the potential to be bounded from below, and for perturbative unitarity for the quartic couplings in the scalar potential are discussed. These conditions are used to find the upper and the lower bounds on c 1 and c 2 and the upper bound on the mixing angle between the two CP-odd states. In Sec. 4, we scan the model parameter space and find the maximum value of σ SI under the constraint discussed in Sec. 3.
Section 5 is devoted to our conclusion.

Model
The model contains a gauge singlet Majorana fermion χ as a DM candidate and a CP-odd gauge singlet scalar a 0 as a mediator. The standard model (SM) Higgs sector is extended into the two-Higgs doublet model. We assume CP invariance both in the dark sector and in the scalar sector.
This assumption guarantees that the Yukawa interaction between χ and a 0 is always pseudo-scalar interaction. The Lagrangian is given by terms with the SM fermions and gauge bosons), Since we assume the CP invariant scalar potential, all the couplings in the potential are real. In this paper, we assume that the thermal relic abundance of χ explains the measured value of the DM energy density [1], and g χ is fixed to realized it for a given parameter set by the freeze-out mechanism.
We impose the condition that the potential has the electroweak vacuum, This electroweak vacuum is realized if m 2 1 and m 2 2 satisfy the following relations.
where λ 345 = λ 3 + λ 4 + λ 5 . In the following, we assume m 2 1 and m 2 2 always satisfy these relations. It is also important that a 0 does not develop vacuum expectation value. Otherwise, the scalar-type Yukawa interaction is induced in the dark sector due to the scalar and pseudo-scalar mixing, and the model is strongly constrained from the direct detection experiments.
After the electroweak symmetry breaking, there are two CP-even scalars (h and H), two CP-odd scalars (a and A), a pair of charged scalars (H ± ), and three would-be Nambu-Goldstone bosons that are eaten by W ± and Z. The physical masses for h, H, a, A, and H ± are denoted to m h , m H , m a , m A , and m H ± , respectively. The two CP-even scalars are mixtures of the CP-even neutral components in H 1 and H 2 , and its mixing angle is denoted by α. Similarly, the two CP-odd scalars are mixtures of the CP-odd neutral components in H 1 and H 2 and also a 0 . Its mixing angle is denoted by θ.
We introduce the following notations for later convenience, Let us comment on the types of the Yukawa interaction. The model is classified into four types based on the Yukawa interaction between the two-Higgs doublet fields and the SM fermions, as in the two-Higgs doublet model with softly broken Z 2 symmetry [22][23][24]. In the following analysis, we choose the type-I Yukawa interaction where H 2 couples to the SM fermions but H 1 does not.
The following discussion is independent from the types of the Yukawa interaction because the type dependence is negligible for large σ SI region of the parameter space as we showed in our previous work [19], and the purpose of this paper is to find the maximum value of σ SI for a given parameter set.
We can express λ i (i = 1, 2, 3, 4, 5), κ, and m 2 a 0 by the mixing angles and mass eigenvalues. In the followings, we take the mixing angle in the CP-even scalars as α = β −π/2. This choice predicts the same hW W and the hZZ couplings as in the SM. We also take M = m H = m A = m H ± . This choice of the mass parameters enhances the custodial symmetry in the scalar potential, and thus the constraints from the electroweak precision measurements are automatically satisfied. With these parameter choices, the parameters of the scalar potential are given by As can be seen, |λ 1,2,3 | < 1. We can also show that |λ 4,5 | < 1 with a condition for the boundedness of the scalar potential (Eq. (3.17)) discussed in the next section.

Theoretical constraints on the scalar potential
In this section, we discuss the condition for the electroweak vacuum as the global minimum of the scalar potential, the conditions for the potential to be bounded from below, and the perturbative unitarity for the quartic couplings in the scalar potential. These constraints are used to find the upper and the lower bounds on c 1 , c 2 , and θ.

Vacuum structure
Vacua other than the electroweak vacuum can exist depending on the given parameter sets.
We study the vacuum structure at the tree level and impose the condition that the electroweak vacuum should be the global minimum. It is not necessary for the electroweak vacuum to be the global minimum if its lifetime is much longer than the age of our Universe. However, the lifetime is much shorter than the age of the Universe in most of the parameter space. 2 Therefore, we adopt the condition to be the global minimum in the current analysis.

3.1.1
In this case, the stationary condition for the scalar potential is given by At the vacuum with a 0 = 0, the potential energy is given by This should be larger than the potential energy at the electroweak vacuum, Therefore, we obtain the following condition for m 2 a 0 < 0, As a result, we obtain the upper bound on c 1 or c 2 for a given parameter sets.

One of H 1 and H 2 is zero
We investigate H 1 = 0 and H 2 = 0. Without lose of the generality, we can parametrize the vacuum as A stationary condition of this vacuum is given by This condition is obtained in both a 0 = 0 and a 0 = 0 cases. Since m 2 3 = 0, see Eq. (2.8), this condition implies σ 2 = 0. This is contradict to H 2 = 0. Therefore, we do not have vacua that satisfy H 1 = 0 and H 2 = 0.
In the same manner, we can show that we do not have vacua that satisfy H 1 = 0 and H 2 = 0.
We simplify the analysis as much as possible by using the gauge invariance in the potential.
Without lose of the generality, we can parametrize the Higgs fields as where σ 1 is positive and σ 2 , π 0 2 , and π + 2 are real numbers. In this case, since the analysis is complicated, we rely on numerical analysis.

Conditions for the potential to be bounded from below
The potential should be bounded from below. In other words, the potential should be positive for the region where the field values are extremely large. We find the following seven conditions for the boundedness of the scalar potential.
The derivation is given in Appendix A. 3 As can be seen, Eqs. The constraint on | sin θ| by using this result for m a = 100 GeV is shown in Fig. 1. We find that | sin θ| 0.21 for m A = 600 GeV, and | sin θ| 0.13 for m A = 1 TeV.
For s β−α = 1 and M = m H = m A = m H ± , we can also simplify other conditions with physical 3 The scalar potential discussed in Ref. [26] is the same as in this paper, but they find that the second condition in Eq. (3.15) should be applied for c1 or c2 < 0. The condition given in Ref. [28], which was derived from the result given in Ref. [29], is consistent with our result. The condition given in Ref. [30] is different from ours.

Perturbative unitarity
Constraints on scalar quartic couplings are often derived from the perturbative unitarity of two scalars to two scalars scattering processes. There are nine scalars in the model. Therefore, the two to two scattering matrix that only includes scalars is a 45 × 45 matrix. Since we consider the high energy limit and ignore the gauge couplings, the scattering processes are s-wave. In the following analysis, the Yukawa coupling, g χ , often takes O(1) value, and thus we also include DM two body initial and final states in the matrix. The DM particle takes two helicity states. In the high energy limit, we find that two DM particles into two DM particles processes are s-wave, and two DM particles into two scalars processes are p-wave. The former processes give stronger bound on g χ .

Spin-independent scattering cross section
We discuss the maximum value of σ SI under the constraints discussed in Sec. 3. We find upper bounds on c 1 and c 2 from the stability of the electroweak vacuum, and lower bounds from the boundedness of the scalar potential. Since g χ ∼ O(1) in the large σ SI regime [19], the perturbative unitarity also gives relevant constraint in the parameter space.
The left panel in Fig. 2 shows that the contours of σ SI for m DM = 800 GeV with the conditions discussed in Secs. 3.1 and 3.2. The other parameters except λ a are the same as one used in Fig. 8 in Ref. [19], namely m H = m H ± = m A = 600 GeV, m a = 100 GeV, t β = 10, θ = 0.1, and λ a = 1.
It is clearly shown that σ SI is larger in the larger |c 2 | region as discussed in Ref. [19]. It is also shown that there is an upper bound on σ SI from the condition discussed in Sec. 3.1. A large As discussed in Ref. [19], both σ SI and σv depend on the a-a-h coupling, g aah , that is given by for s β−α = 1 and m H = m H ± = m A . Combining these two equations, we find This condition is not satisfied with the large g aah , and thus the large g aah induces the a 0 = 0 vacuum. On the other hand, the large g aah is necessary to obtain the larger σ SI . Therefore, there is a correlation between σ SI and the condition of the stability of the electroweak vacuum.
We can also see from and t β = 10 in the following analysis, but the following results are insensitive to the choice of t β .
The larger λ a allows us to take larger g aah while keeping a 0 = 0, which can be seen from Eq. (4.4). On the other hand, the larger λ a implies the breakdown of perturbative calculation at a higher energy scale. In our analysis, g χ is typically O(1) to obtain the measured value of the DM energy density, and it also implies the breakdown of perturbative calculation at a higher energy scale. We calculate the running of the couplings at the 1-loop level and estimate the cutoff scale Λ as the highest scale that satisfies Eqs. (3.11), (3.28), and (3.32). In the calculation, we assume that the input parameters are given at µ = m A . The beta-functions of the couplings we used are given in Appendix B. The smaller λ a at µ = m A becomes negative at higher scale because a 0 couple to the fermionic DM that gives a negative contribution to the beta function of λ a . On the other hand, the beta function is proportional to λ 2 a and positive for the larger λ a . The cutoff scale gives the upper and the lower bounds on λ a at µ = m A . Figure 3 shows the contours of σ SI in λ a -c 2 planes. It is shown that the larger λ a at µ = m A keeps its value positive at any higher scale. We find that it is easy to make the cutoff scale higher than O(100) TeV by choosing λ a 1.5. Thus we can expect that unknown UV physics does not modify our results for λ a 1.5. We also find that σ SI is maximized along the boundary of the orange shaded region where the electroweak symmetry is not broken. For c 1 = c 2 , Eq. (3.5) is  (10) TeV.
The other color notation is the same as in Fig. 2. simplified as In the following analysis, we choose c 1 = c 2 = 0.99c * for given parameter sets. This choice of c 1 and c 2 maximizes σ SI . Figure 4 shows the contours of σ SI in λ a -θ plane. We find that σ SI is larger in the smaller θ regime. This is because the smaller θ requires larger g χ to obtain the measured value of the DM energy density. The left and right panels are for m A = 600 GeV and 1 TeV, respectively. We find that σ SI is almost independent from m H,H ± ,A for θ < 0.01. This is because the heavier scalars almost decouple both from the DM annihilation processes and from the loop contributions to σ SI .
The cutoff scales are the only difference if we change m A ; a larger m A predicts higher cutoff scales.
In the following analysis, we take θ = 0.001. With this choice, σ SI is maximized and is independent from m A . We also take m A = 600 GeV in the following, which gives us a conservative bound from the RGE analysis.

Conclusion
The 2HDM+a is a DM model that can explain the measured value of the DM energy density by the freeze-out mechanism and can avoid the constraint from the Xenon1T experiment. The leading order contribution to σ SI is given at the loop level, and σ SI can be large enough for the model to be tested by the forthcoming direct detection experiments.
In this paper, we have investigated the maximum value of σ SI under theoretical constraints.
We take into account the stability of the electroweak vacuum, the condition for the potential bounded from the below, and the perturbative unitarity of two to two scattering processes. As shown in Fig. 2, large |c 1 | and |c 2 | make σ SI larger. However, the condition for the stability of the electroweak vacuum gives upper bounds on c 1 and c 2 , and the potential boundedness condition gives lower bounds on them. As a result, there exists the maximum value of σ SI for a given parameter set. It is also shown that σ SI is maximized for c 1 = c 2 = c * , where c * is the maximum values of c 1 and c 2 to keep the electroweak vacuum as the global minimum of the scalar potential.
With this choice, the result is insensitive to t β . We found that a smaller θ makes σ SI larger, as shown in Fig. 4. For the small θ regime, σ SI is almost independent of m H,H ± ,A . Finally, in Fig. 5, we found that σ SI can be larger than the prospects of the LZ and XENONnT experiments for m DM 600 GeV. We also found that the perturbative unitarity gives an upper bound on m χ . The maximum value of the σ SI is ∼ 5 × 10 −47 cm 2 for m A = 600 GeV where the cutoff scale of this model is estimated as 100 TeV. Therefore, if the LZ and XENONnT experiments observe the DM signal in future, then this model predicts 600 GeV m DM 1TeV.

A Condition for the potential to be bounded below
The potential should be bounded below, namely the potential should be positive for the region where the field values are extremely larger. In this section, we derive the condition for the bounded below.
By imposingṼ > 0, we find constraints on the parameters.
The sign of the left-hand side is determined by the sign of the terms in the big parenthesis in the right-hand side. It takes minimum if the terms depending on θ vanish, namely, √ a cos θ− √ b sin θ = 0. Its minimum value is c + √ ab. Since sin θ cos θ > 0, if c + √ ab > 0 then the right-hand side is always positive.
Since λ 2 > 0 and λ a > 0 are already guaranteed, this is positive if This is Eqs. (3.14).