On the Viability of Minimal Neutrinophilic Two-Higgs-Doublet Models

We study the constraints that electroweak precision data can impose, after the discovery of the Higgs boson by the LHC, on neutrinophilic two-Higgs-doublet models which comprise one extra $SU(2)\times U(1)$ doublet and a new symmetry, namely a spontaneously broken $\mathbb{Z}_2$ or a softly broken global $U(1)$. In these models the extra Higgs doublet, via its very small vacuum expectation value, is the sole responsible for neutrino masses. We find that the model with a $\mathbb{Z}_2$ symmetry is basically ruled out by electroweak precision data, even if the model is slightly extended to include extra right-handed neutrinos, due to the presence of a very light scalar. While the other model is still perfectly viable, the parameter space is considerably constrained by current data, specially by the $T$ parameter. In particular, the new charged and neutral scalars must have very similar masses.


Introduction
The smallness of neutrino masses suggests a mass generating mechanism distinct from the usual Higgs mechanism, which resides in a scale different from the electroweak one. From neutrino oscillation experiments, we know that neutrinos are massive and that mass and flavor eigenstates do not coincide. Besides, other terrestrial [1,2] and cosmological [3,4] experiments indicate that neutrino masses should be below the eV scale. Therefore, if the same Higgs mechanism is responsible for the top and neutrino masses, then the Yukawa couplings would span twelve orders of magnitude, evincing an unpleasant and inexplicable hierarchy.
A well known alternative is the seesaw mechanism [5][6][7]. In this scenario, the light neutrino masses are suppressed by some heavy physics, for instance, right-handed Majorana neutrino masses [5,6,8,9]. What typically happens is that the scale at which new physics can be found is extremely high, much above the TeV scale, rendering the model intangible, except for the possible presence of neutrinoless double beta decay. 1 The latter could also originate from some physics that do not comprise the main contribution to neutrino masses [13,14], and hence it does not consist of a test of the seesaw mechanism by itself.
Another possibility is to generate neutrino masses by a copy of the Higgs mechanism, having a second Higgs doublet, but with a much smaller vacuum expectation value (vev).

JHEP12(2015)160
This can be achieved in a two-Higgs-doublet model (2HDM) where one of the scalars gives mass to the charged fermions, while the other one acquires a very small vev and generates neutrino masses with O(1) Yukawa couplings, a neutrinophilic 2HDM. As a consequence, neutrino masses would generically require new physics at the TeV scale (or even lower). For instance, by imposing a lepton number symmetry and adding three right-handed neutrinos which carry no lepton number, a type I seesaw mechanism can be realized below the TeV scale [15]. Moreover, lepton number could be conserved and a Z 2 symmetry [16,17] or a global U(1) [18] could be used to prevent the SM Higgs boson to couple to neutrinos, yielding Dirac neutrinos. Also, the 2HDM could be augmented by a type III seesaw and a µ − τ symmetry, giving rise to interesting LHC phenomenology [19]; or by a singlet scalar and a Z 3 symmetry, possibly generating lepton flavor violating signals [20]. It is important to note that such models are stable against radiative corrections [21,22].
On general grounds, a new symmetry is typically invoked to prevent the first scalar doublet from coupling to neutrinos as well as to enforce the second one to interact only with them. These models introduce a minimal new field content which should materialize as particles below the TeV scale. The presence of such a low scale in the theory might have important phenomenological consequences, like the presence of light scalar particles (for instance, supernova energy loss strongly constrains such scenarios [23]). After the discovery of a 125 GeV scalar by the LHC experiments, new limits from electroweak precision data can be derived on the allowed parameter space of such models. The purpose of this manuscript is to investigate to what extent these minimal neutrinophilic 2HDMs can survive electroweak precision data scrutiny.
In section 2 we briefly review the neutrinophilic 2HDMs which we will study in this work. In section 3 we describe the theoretical and experimental constraints that will be imposed on these models in section 4. Finally, in section 5 we present our conclusions.

Neutrinophilic two-Higgs-doublet models
We first start by making general considerations on the 2HDM and the link to neutrino masses. The most general scalar potential for a 2HDM is where Φ 1 and Φ 2 are two scalar doublets with hypercharge Y = +1. For the vacuum expectation values of the two scalars, we adopt the notation Φ 1 = v 1 / √ 2, Φ 2 = v 2 / √ 2, and we pick Φ 2 to be the one responsible for neutrino masses. In order to have sizable Yukawa coupling for neutrinos, it is required that In principle, the parameters m 2 12 , λ 5 , λ 6 , and λ 7 can be complex. Nevertheless, in all models we analyze, the symmetries will forbid both λ 6 and λ 7 , and only m 2 12 or λ 5 will be allowed to be non-zero. A single phase of the aforementioned parameters can always be JHEP12(2015)160 absorbed in a redefinition of the scalar fields, and therefore we can take all scalar potential parameters to be real without loss of generality.
To forbid the coupling between neutrinos and Φ 1 , a symmetry is called for. In this minimal setup, there are two straightforward examples. The first possibility is a Z 2 symmetry under which only Φ 2 and the right-handed neutrinos are charged, forcing m 12 = λ 6 = λ 7 = 0. An alternative is to trade the Z 2 by a global U(1), yielding, in principle, m 12 = λ 5 = λ 6 = λ 7 = 0. In this case, to avoid the presence of a massless Goldstone boson, a soft breaking is introduced by having a non-zero but small m 12 . On the other hand, if a softly broken Z N , N > 2, symmetry is postulated, λ 5 might be forcefully zero as well, making this case identical to the U(1) scenario. Therefore, the phenomenology of a softly broken Z N >2 model is identical to the softly broken U(1) case. Anyhow, in all realizations we will study here λ 6 = λ 7 = 0, so these couplings will be disregarded henceforth.
One last option that one could consider would be to gauge the U(1) symmetry, avoiding the massless Goldstone boson. Nevertheless, in such a scenario, the corresponding gauge boson as well as one of the neutral scalars would be extremely light, with mass around the v 2 scale. This seems, at first glance, phenomenologically quite problematic. We do not investigate this possibility here as it would require a completely different study compared to the other two cases.
The two complex scalar SU(2) doublets can be written as After electroweak symmetry breaking, three Goldstone bosons become the longitudinal modes of the W and Z bosons. Then, the remaining scalar spectrum is composed of two charged particles, H ± , two CP-even neutral bosons, h and H, and one CP-odd neutral boson, A. The physical fields are given by where the angles α and β are associated with the rotations that diagonalize the mass matrices where λ 345 ≡ λ 3 + λ 4 + λ 5 . We will see below that both α and β are expected to be very small. Hence, h behaves very similarly to the SM Higgs, while the scalars H, A, H ± develop neutrinophilic interactions in the Yukawa sector, as described below where m ν i are neutrino masses and U i is the PMNS matrix. As we will see in section 4, the tree-level stationary conditions on the potential, ∂V /∂Φ i = 0, can be used to write the diagonal mass parameters m ii as functions of m 2 12 , the quartic couplings and the vevs. With that in mind, we can consider the quartic couplings as free parameters and express them in terms of the physical masses, vevs and mixing angles [24] (see appendix A).
Next we describe the two specific realizations of the neutrinophilic scenarios that will be studied in this paper.

Neutrinophilic 2HDM: Z 2 symmetry
The model to be studied was proposed by Gabriel and Nandi [16]. 2 It consists of a 2HDM where both the right-handed neutrinos and one of the scalar doublets, Φ 2 , are charged under a Z 2 symmetry. The consequence is that the masses of the charged fermions come solely from the Φ 1 vev, and neutrinos, which are Dirac fermions in this scenario as the authors impose lepton number conservation, couple exclusively to Φ 2 . This extra symmetry can, in principle, be dropped allowing for Majorana neutrinos with a low scale realization of the seesaw mechanism. We will also investigate this possibility in our analysis.
In the scalar potential (2.1) of this model, the parameters m 2 12 and λ 6,7 will vanish due to the Z 2 symmetry. The smallness of neutrino masses is explained by the very low scale at which Z 2 is broken, preferably v 2 O(eV). 3 A tiny v 2 /v 1 ratio and the absence of an explicit breaking m 2 12 term leads to almost no mixing between the doublets. The smallness of tan β and tan α can be seen from eqs. (2.5) and (2.6) after imposing v 2 /v 1 → 0. Therefore, apart from its couplings to neutrinos, Φ 1 behaves almost identically to the SM Higgs doublet, so we do not expect any observable deviation from the Higgs couplings to the SM particles, except possibly the loop induced couplings, e.g. hγγ.
The second doublet displays some interesting features. Through the Yukawa coupling, the neutral components couple almost only to neutrinos, while the charged scalars mediate interactions between neutrinos and charged leptons (see eq. (2.7)). The Yukawas are ideally expected to be of O(1). The neutral scalars couple to the W and Z bosons, but notice that triple gauge couplings (TGCs) involving only one scalar are highly suppressed by the small vev, v 2 . Obviously, TGCs with two scalars and one gauge boson are present and may provide a sizeable pair production cross section at colliders, for instance pp → A * → H + H − at the LHC.
The scalar spectrum of this model is quite constrained. By setting m 2 12 = 0 in eqs. (A.6)-(A.8), as well as sin 2 α, sin 2 β 1, we notice that: (i) h is identified as the 125 GeV Higgs particle found at the LHC, and this essentially fixes are bounded to be about or below the TeV scale. When we analyse the viability of this model in section 4, it will turn out that oblique parameters will play a decisive role in constraining it, due to the peculiar structure of the scalar spectrum. The sensitivity of the S parameter to the presence of a very light neutral scalar, m H ∼ O(v 2 ), will essentially rule out the model.

Neutrinophilic 2HDM: softly broken global U(1) symmetry
The second model we study was proposed by Davidson and Logan [18]. Analogously to the other scenario, both Φ 2 and right-handed neutrinos are charged under a new global U(1). The model spans λ 5,6,7 = 0 and a small m 2 12 which breaks the symmetry softly and generates neutrino masses. The presence of the soft breaking mass term, is required in order to avoid a massless Goldstone boson which might create problems with cosmology and electroweak precision data. Neutrinos are Dirac particles, as the Majorana mass term is strictly forbidden by the new U(1). From eq. (A.8), we write and we observe that to obtain simultaneously v 2 ∼ eV and m A ∼ O(100 GeV) one would need m 2 12 ∼ (200 keV) 2 . As said before, to avoid the issues of having a massless Goldstone, instead of softly breaking the new U(1) symmetry, one could also envisage to gauge it. Nonetheless, the theory would contain a very light vector resonance as a consequence of the small vev, and it is not clear if such a model can satisfy all neutrino data and astrophysical constraints. We do not explore this possibility here.
The presence of a non-zero m 2 12 term makes this case fairly different from the last one. From eq. (A.7), we notice that the mass of the neutrinophilic scalar, m H , increases with M , and therefore the H mass in this scenario is not bounded by v 2 as in the previous case. As we will see later, this will ease the constraints from the oblique parameters. Combining eq. (A.2) with the definition M 2 = m 2 12 /(sin β cos β), and imposing tan β = v 2 /v 1, we obtain which indicates that To grasp the impact of this conclusion, assume that m 2 12 = m 2 H v 2 /v. Hence, from eq. (A.8) we see that m A ≈ m H , so the neutrinophilic CP-odd and CP-even scalars, A and H, are degenerate in mass. We emphasize that this degeneracy by itself is not a fine tuning of the model: the degenerate spectrum arises naturally given the symmetries of the scalar potential and the hierarchy between the vevs. As a last comment, we emphasize that since m 2 12 is the only source of U(1) breaking, it is natural in the t'Hooft sense -m 2 12 only receives radiative corrections proportional to itself [21,22].

Theoretical constraints
There are a number of conditions to be fulfilled by the scalar potential. These will be used to constrain the parameter space, ultimately restricting the range of physical scalar masses, having an important impact on the phenomenology of the models. To have stability at tree-level, the following constraints should be fulfilled [28] In addition, the stationary conditions ∂V /∂Φ i = 0 read which allow us to write m 2 ii as functions of m 2 12 , v 1 and v 2 . If m 2 12 = 0, it is easy to see that there are at least two equivalent stable solutions, (v, 0) or (0, v) (although they may not be the global minima). In this case, the vev is precisely the electroweak scale, one of the scalars is exactly the Higgs and the other one is inert. For m 2 12 = 0, these equations cannot be solved analytically. Nevertheless, if m 2 12 v 2 a perturbative approach yields and a symmetric solution interchanging the indices 1 ↔ 2, which reveals that the small vev necessary to satisfactorily explain small neutrino masses might require a correspondingly small m 2 12 parameter. This can be understood intuitively, as the breaking of the U(1) happens only through the soft breaking term m 2 12 . In general, there can be more than one solution satisfying the stationary conditions (3.2), and hence different non-trivial and nondegenerate minima (v 1 , v 2 ) and (v 1 , v 2 ) might coexist. It is possible to check analytically if the chosen vacuum is the deepest one in the potential for a 2HDM with λ 6,7 = 0 [29]. In this case, the potential describes a Z 2 symmetry softly broken by m 2 12 . Both models we deal with here are special cases of such scenario. In the absence of an explicit breaking, that is m 2 12 = 0, there can be multiple minima, but they are degenerate and hence stability is not threatened. This is the case of the Z 2 model we analyze. For the softly broken U(1) model, it can be shown that the chosen vacuum is the deepest one (at tree-level) if and only if the following condition is satisfied [29]: Although for a general 2HDM scenario this bound may be important, for the neutrinophilic case we have checked that it does not lead to any significant effect on the parameter space, after the other constraints are taken into account, but we include it in the analysis of the softly broken U(1) model for completeness.

JHEP12(2015)160
Another theoretical requirement is to satisfy the tree-level pertubative unitarity condition [30][31][32]. If the quartic couplings are too large, the lowest order amplitudes for scalarscalar scattering may violate unitarity at high enough scales, requiring additional physics to mitigate this issue. To obtain the constraint, the scalar-scalar S matrix is computed and the following conditions are imposed on its eigenvalues To have an idea of the impact of these bounds, one can conservatively assume that all |λ i | should be smaller than 8π (the actual bound is always more stringent than that). Some authors prefer to use a stronger limit of 4π. We checked that this does not change very much the allowed regions. Evidently, even if tree-level unitarity is satisfied, loop corrections could still play an important role leading to violation of unitarity at some scale and thus demanding the presence of new physics below such energies. This could be particularly relevant when some of the tree-level constraints are just barely satisfied, as the size of the quartic couplings could enhance the loop contributions. Nevertheless, we only take into account unitarity constraints at tree-level, as a full one loop evaluation of the parameter space is beyond the scope of this manuscript.

Electroweak data constraints
Oblique parameters. The impact of a second Higgs doublet in the so-called electroweak precision tests (EWPT), encoded in the Peskin-Takeuchi parameters S, T , and U [33], has been studied in the literature to a great extent (see for instance refs. [34][35][36]). These are radiative corrections to the gauge boson two point functions, known as oblique corrections. For the precise expressions of S, T , U , we point the reader to the aforementioned references.
The S parameter encodes the running of the neutral gauge bosons two point functions (ZZ, Zγ and γγ) between zero momentum and the Z pole. Therefore, it should be specially JHEP12(2015)160 sensitive to new physics at low scales, particularly below the Z mass. Thus, we expect it to be important in the presence of very light neutral scalars, as is the case for the Z 2 model. The T parameter measures the breaking of custodial symmetry at zero momentum, that is, the difference between the W W and the ZZ two point functions at q 2 = 0. It usually plays a significant role in constraining the parameter space of particles charged under SU(2) L . Splitting the masses of particles in a doublet breaks custodial symmetry and affects T . As we will see later, in the softly broken U(1) scenario, the T parameter will provide the major constraint on the mass splitting m H ± − m A , forcing the scalar spectrum of this model to be somewhat degenerate. Last, and this time least, the U parameter (or better, the combination S + U ) is somewhat similar to S but for the W bosons, being sensitive to light charged particles in the loops. Given the fact that light charged particles are excluded by LEP data [37,38], usually U is the least important of these three precision parameters, having a minor impact on the model phenomenology, we have checked that this is indeed the case for all scenarios analyzed here.
To evaluate the impact of the EWPT on the neutrinophilic 2HDM scenarios, we calculate S, T , and U using the results available in ref. [36], and we use the latest GFITTER values for the best fit, uncertainties and covariance matrix [39], composing the χ 2 function as with X i = ∆S, ∆T, ∆U and the covariance matrix σ 2 ij ≡ σ i V ij σ j , in which (σ 1 , σ 2 , σ 3 ) = (0.11, 0.13, 0.11). As we are interested in the goodness of fit of the model to the EWPT data, the 1, 2, and 3σ regions are calculated using χ 2 = 3.5, 8.0, 14.2, respectively.
Higgs invisible width. When the first doublet acquires a vev, triple scalar vertices like hSS (S = H, A) are induced. Therefore, light neutral scalars with 2m S < m h could contribute to the Higgs invisible width h → SS, and sequentially S →νν. Because of the small tan β of the model, the Higgs boson couplings to the Standard Model particles is basically unchanged. Hence, the contribution to the Higgs total width due to the invisible decay will suppress all Standard Model branching fractions by the ratio Γ SM h /Γ new h . In this scenario, as the only modification to the Higgs branching fractions is the addition of an invisible channel, the LHC 8 TeV data bound is BR(h → invisible) < 0.13 at 95% CL [40].
In our framework, the decay rate of such a process is given by [41] Γ While the couplings between the SM Higgs and the SM fermions, g hf f = m f /v, are well below one due to the suppression by the EW scale (except for the top, to which the Higgs cannot decay), the trilinear scalar couplings are typically much larger, g hSS ∼ m 2 h /v ∼ 60 GeV, unless there is some sort of cancellation happening [41]. Therefore, SM Higgs decays to lighter scalars may have an important phenomenological impact, see e.g. ref. [42], specially because the total Higgs width in the Standard Model is predicted to be very small, around 4.07 MeV [43].
Higgs to diphoton. The charged scalars will contribute to the h → γγ width, and thus we also analyse the impact on this observable. 4 The h diphoton width is a destructive interference effect mainly between W and top loops, where the latter dominate. Charged scalars contribute with the same sign as the W , and their contribution usually do not overcome the top one. Therefore, we expect h → γγ to be somewhat suppressed in most cases. The expression for the h → γγ width at one loop can be found in many papers, see, for instance ref. [44]. For reference, the current ATLAS+CMS combination value of the Higgs to diphoton signal strength is µ γγ = 1.16 +0. 20 −0.18 [45].
Z invisible width. We also have to consider possible extra contributions to the Z invisible width coming from the decays Z → Sνν with S = A, H and m S < m Z . In the model with a softly broken U(1) symmetry, the expression for Γ(Z → Sνν) = Γ(Z → Aνν) + Γ(Z → Hνν) can be easily calculated and reads where m S is the mass of the neutrinophilic scalars H and A, which are degenerate in mass, and the total width is given by (3.13)

JHEP12(2015)160
We also define λ(a 2 , b 2 , The ratio between m 2 ν,tot ≡ m 2 ν i and v 2 2 arrives from the neutrino Yukawas. Clearly, if the Yukawas are small, both widths vanish, so we expect this bound to be more significant for lower v 2 and larger neutrino masses. To constrain extra contributions from new physics to the Z invisible width, we use LEP result Γ exp (Z → invisible) = 499.0(15) MeV and the Standard Model prediction Γ SM (Z → invisible) = 501.69(6) MeV [38], which yields Γ NP (Z → invisible) < 1.8 MeV at 3σ (notice that there is a mild 2σ discrepancy between the data and the SM predicted value). In the case of the Z 2 symmetry model, one must take care while doing the computation, since m H m Z , as the expression for the width has an infrared divergence, which cancels out with radiative conditions. As we will see in section 4, the other constraints will exclude most of the parameter space of this model. For this reason we will not discuss the constraints from the Z invisible width in this scenario.
Collider bounds on charged scalars. The charged scalars can be pair produced directly at colliders via s-channel off shell photon or Z exchange. Due to the neutrinophilic character of the second Higgs doublet and small admixture with the SM degrees of freedom, the charged scalars decay almost only to ν. Therefore, we use the corresponding LEP bound, i.e. m H ± > 80 GeV [37,38].
It is not clear how LHC data improves the situation. There has been some studies on the LHC sensitivity to such charged scalars, mainly focused on 14 TeV center of mass energy [46][47][48][49][50], but to the best of our knowledge, there has been no dedicated experimental search for charged scalars in neutrinophilic 2HDMs. As v 2 is very small, the main production modes of H ± would be pair production through vector boson fusion or off-shell s-channel photon and Z exchange, and the tipical t → H + b would be absent due to small tan β. The LHC sensitivity then would come mainly from opposite sign dilepton plus missing energy, which has SM W pair production as an irreducible background. Moreover, the branching ratios of the charged scalar depend on the neutrino masses and the mass ordering. If the τ ν branching ratio is dominant, the sensitivity is expected to be smaller. Therefore, to be conservative, we will scan the parameter space considering only the LEP bound.
Anomalous magnetic moments and other constraints. In principle, the charged scalars could also contribute to charged lepton g−2 values, but the corresponding amplitude at one loop is suppressed by m 4 /m 4 H ± (see ref. [51] for a recent analysis on the impact of a second Higgs doublet on the muon g − 2). We have checked that the 1-loop contribution to both muon and electron g −2 is negligible due to that suppression, while the tau g −2 is not measured with enough precision to pose a bound. For a general 2HDM, it has been noticed that two loop Barr-Zee diagrams [52] can be more important than 1-loop contributions, JHEP12(2015)160 but this is not the case in the neutrinophilic 2HDM, as the charged lepton couplings to H and A are suppressed by tan β. 5 Therefore we conclude that the electron, muon and tau g − 2 measurements do not pose any bound on this scenario.
Flavor physics constraints have also been studied in the literature. The charged scalars will mediate lepton flavor violating decays. In µ → eγ, for instance, the additional branching ratio is proportional to (m H ± v 2 ) −4 . Because of that, it is always possible to evade this bound for large enough values of v 2 (or m H ± ). The limits we derive here are independent of a particular choice of v 2 as they concern directly the spectrum.

Analysis of the models
For each model we generate ≈ 10 7 points. For each of those points we calculate the corresponding scalar potential parameters and verify if they fulfill the constraints described in section 3. We only show on our plots the allowed points, which are about 10% of the generated ones. Unless stated otherwise, the points are color coded accordingly to the fit to EWPT data: blue, green, and red correspond to the 1, 2, and 3σ allowed regions, while gray points are excluded at 3σ or more.
2HDM with a Z 2 symmetry. Let us first discuss the results for the 2HDM with a Z 2 symmetry. As discussed in section 2.1, the model has a very light neutral scalar. In fact, we verified that eq. (A.7) and the perturbative unitarity conditions (3.5) require m H 10×v 2 . Moreover, as the scalar potential parameters λ i and m 2 ij can be written in terms of the physical masses and the vevs, we perform a scan in the physical parameter space, imposing the following conditions 0.01 eV < m H < 1 GeV, 124.85 GeV < m h < 125.33 GeV, 70 GeV < m H ± < 1 TeV, 1 GeV < m A < 1 TeV, −π/2 < α < π/2, 0.01 eV < v 2 < 1 MeV.
The Higgs mass range is taken from the ATLAS+CMS measurements combination in ref. [53]. Note that although α has to be small we did a scan over the whole physical range of this parameter, since we wanted to be as general as possible. By using a logarithmic prior, for convenience, we found that indeed α is small.
Using the power of perturbative unitarity constraints, we found that the CP-odd and charged scalars are restricted to be below ∼ 600-700 GeV. This can be easily understood from eqs. (A.8) and (A.9). Since M 2 ∝ m 2 12 = 0 and λ 4,5 cannot be too large, the masses cannot go arbitrarily above the electroweak vev.

JHEP12(2015)160
Moreover, the presence of a very light scalar in the spectrum, below the GeV scale, yields a substantial negative contribution to the S parameter. The impact of the EWPT can be seen in figure 1, where all points scanned were projected in the S × T plane and the allowed region by EWPT was drawn. Remarkably, only very few points (in red) were found which provide a viable model, within the 3σ allowed region for the EWPT. From our scan, it can be concluded that: the T parameter strongly prefers m A ≈ m H ± or a lighter H ± with m H ± ∼ 150 GeV together with a m A > 300 GeV; while the S parameter, although it depends very mildly on the charged and pseudoscalar masses, exhibits a slight preference to this latter region. All in all, the values of S are always below ∼ −0.25, revealing a tension with EWPT always above the 2.97σ level. 6 As an example, we obtained the following scalar spectrum, which is allowed at 2.99σ: From this analysis, we can conclude that the 2HDM with a Z 2 symmetry is definitely very disfavored by data. It is not even clear that the region found which is in the 3σ border of EWPT is really viable. A closer look into this region of the parameter space reveals that these points suffer from at least one of the following worrisome situations: (i) the e 1 scattering amplitude, in eq. (3.5), is on the verge of violating unitarity, with at least about ∼ 98% of the bound saturated; (ii) the same for a + scattering amplitude, with at least ∼ 98% of the bound saturated; (iii) the stability condition is very fragile, with the third condition of eq. (3.1) satisfied with a relative difference of less than ∼ 4 × 10 −4 ; and (iv) the same but for the second condition of eq. (3.1), satisfied with a relative difference of less than ∼ 0.05. Therefore, given this delicate region of the parameter space, it would be important to include radiative corrections to see if the stability and unitarity of the model still holds at one loop. Notice that by using 4π as the perturbative unitarity limit this small region disappears. A possible way to evade these problems could be to have a larger v 2 so that the mass spectrum, specially m H , becomes more flexible. Nevertheless, unless v 2 O(GeV), the problem does not disappear, strongly disfavoring this minimal model as an explanation for neutrino masses.
One could now be tempted to include a right-handed neutrino contribution, dropping the lepton number conservation symmetry of the model. In fact, as v 2 is small, it may be possible to have a low-energy realization of the type I seesaw scenario which leads to observable sterile neutrino phenomenology, and hopefully could increase a bit the value of the S parameter to make the model viable. As the effect on S grows with the mass of the fermions in the loop, we make a distinction between two regimes: the right-handed neutrinos can be below or above the GeV scale. In the first, what happens is that the contribution to the S parameter is suppressed by the ratio between these small masses and JHEP12(2015)160  the Z mass and can be neglected (for instance, the active neutrino contribution to S is virtually zero). In the second case, although the sterile neutrino masses might be large, the coupling to the Z is suppressed by the active-sterile mixing which generically goes as the ratio between the active to sterile neutrino masses, m ν /m N . Therefore the impact of right-handed neutrinos is never large enough to substantially change the S parameter. 7 For completeness, we also show in the right panel of figure 2 the impact of the Higgs invisible width measurement in the tan α × tan β plane. Given the preference for heavier S, we will consider the case where only h → HH is present. From eq. (3.11), since m 2 12 = 0 in this model, the g hHH coupling can be rewritten in the limit of small β and α as which can be sizable only if α β, explaining the behavior of the excluded region (orange) in figure 2. Since the ratio α/β is already constrained by the theoretical limits (see eq. (2.9)), this constraint turns out to be less stringent than the others. As a last comment,

JHEP12(2015)160
the charged scalars could also have an impact on h → γγ. In the small 3σ allowed region, the modifications to the diphoton width are generically between ±10%, depending on the precise values of λ 3 . This quartic coupling only affects m H , so it is only weakly bounded by perturbative unitarity.
Finally, since the smallness of m H causes the tension with EWPT, one may wonder what is the impact of loop corrections on the scalar spectrum of this model. Generically, in a 2HDM with Z 2 symmetry, the charged and CP-odd mass matrices are not modified by one-loop corrections. The CP-even matrix receives radiative corrections of the form [55] where the second term comes from the one-loop effective potential. Since ∆m 2 ij are solely functions of masses and quartic couplings, the dependence of the CP-even mass matrix on v 1,2 is preserved at one-loop level, implying a small value for m H if v 2 is small. We have checked by explicit calculations that the corrections to m H are at the most a factor 100, which is still insufficient to solve the problem with the S parameter.
2HDM with a global U(1) symmetry. We now focus on the phenomenology of the softly broken U(1) model. A non-zero m 2 12 term allows for heavier H, presenting a major change in the phenomenology with respect to the previous model. Without the requirement of a light scalar, we enlarge the scanned region accordingly. The absence of the λ 5 quartic coupling makes the pseudoscalar degenerate in mass with H (to first order in v 2 ). Therefore we perform an initial scan of the spectrum parameter space, this time in the region 10 GeV < m H < 1 TeV, 124.85 GeV < m h < 125.33 GeV, 70 GeV < m H ± < 1 TeV, m A = m H , −π/2 < α < π/2, 0.01 eV < v 2 < 1 MeV, as well as a second scan with m H ± and m A heavier then 1 TeV and almost degenerate. We follow the same procedure as before, showing only the points allowed by perturbative unitarity and stability constraints. The results are presented in figure 3.
In contrast to the previous case, due to the possibility of obtaining a heavier H in the mass spectrum, there is a region of the parameter space of this model which passes the electroweak precision tests and theoretical constraints. The behavior of the T parameter is similar to the previous scenario: either the mass splitting between A and H ± is at most ∼ 80 GeV, or the charged scalar is around 100 GeV while m H = m A > 150 GeV, with negative values of T for larger m H ± . This explains the strong correlation on the allowed region in the upper left panel of figure 3. We also present the projection of these points in the S × T plane in the upper right panel of figure 3. For the α and β parameters we find that the allowed region is tan β 10 −6 and α 5β. As discussed in the previous sections, this model can also accommodate a pair of neutral scalars (S = H, A) satisfying m S < m h /2 if m 2 12 is small enough. In this case, the constraints coming from the Higgs invisible decays are similar to those described for the model with a Z 2 symmetry and turn out to be relatively weak. On other hand, the Z invisible width can provide valuable constraints when the channel Z → Sνν is open. To perform this analysis we scan over the oscillation parameters, imposing the perturbativity condition Γ S < m S /2. We show on the bottom right panel of figure 3 the excluded region (orange points) under these assumptions in the m H × v 2 plane. The region m S < m Z /2 is completely excluded, because in this case we integrate over the poles of the off-shell scalars in Z → H(A * → νν) and Z → A(H * → νν), enhancing the decay rate by orders of magnitude.

JHEP12(2015)160
For a heavy enough H ± , as can be seen in the lower left panel of figure 3, the h → γγ signal strength is diminished by about ∼ 5%. We can understand this non decoupling feature by noticing that the h H + H − coupling is −iλ 3 v, which in turn has a correlation with m H ± , specially in the larger mass region. This can be understood by noticing that, in eq. (A.3), for large m H ± , we have Typically, α 5β, which corresponds to a strong correlation between λ 3 and m H ± , and

JHEP12(2015)160
this is the denser region around µ γγ = 0.95. However this is not always the case, and the correlation is lost when the ratio of sines is closer to 1, now corresponding to the sparser points with a much weaker correlation. Nevertheless, we see that for a heavy enough charged scalar, the contribution to the Higgs diphoton width is always negative. One could ask if it is also possible to have a Majorana mass term, since the U(1) symmetry is softly broken. First, as pointed out in refs. [54,56], the impact of heavy right-handed neutrinos via loop effects on electroweak precision observables is very small. Therefore, there is no significant interplay between this and the scalar sector of the model, and thus the phenomenology studied here would be essentially unchanged. On the other hand, if we consider an UV completion that simultaneously breaks the symmetry and originates a Majorana mass term, we find that such scenario is non-minimal, i.e., at least two new fields have to be included.
Comments on non-minimal models. Due to the large number of possible variants, performing exhaustive analises of non-minimial models is unpractical and well beyond the purpose of this paper. Nonetheless, we may glimpse the phenomenology of some representative cases.
A neutrinophilic 2HDM with a softly broken Z 2 symmetry would surely be allowed by data, in contrast to the spontaneously broken Z 2 scenario. Such model would be more general than the two models considered here, as it would span a non-zero value of both m 12 and λ 5 . As can be seen from eqs. (A.7) and (A.8), the simultaneous presence of these terms in the scalar potential lifts the degeneracy between the neutral scalars A and H. In fact, we have checked that the allowed region in the plane m A × m H ± is very similar to the one exhibited in figure 3 (top left panel), except for the fact that the T parameter now implies a correlation only between m H ± and m A . If the neutral scalar H decays dominantly to neutrinos, as it is likely to happen, it would be very difficult to probe it by resonant production at colliders.
Another way of evading our limits would be to enlarge the particle spectrum of the spontaneously broken Z 2 model (or generically any Z N ), for instance, by adding a scalar singlet S, doublet Φ 3 or triplet ∆, all charged under the new symmetry. In the singlet case, a triple or quartic term SΦ † 1 Φ 2 or S 2 Φ † 1 Φ 2 could be present in the potential for a judicious choice of charges. After the singlet acquires a vev, this term would play a role similar to the m 12 soft breaking term allowing for larger values of m H . However, the quartic S † SΦ † 1 Φ 1 , always present, would induce a Higgs-singlet mixing. This would diminish all Higgs couplings to fermions and gauge bosons by a factor sin θ where θ is the corresponding mixing angle. The mixing is constrained by Higgs production cross section measurements to be sin 2 θ 0.2 [45]. In the case of adding a third doublet, the triple coupling is impossible, but a quartic one could be present. Last, in the case of a scalar triplet, although a triple coupling would be possible, a large triplet vev could irrevocably disturb electroweak precision tests, especially the T parameter.

Conclusion
We performed an analysis of the minimal neutrinophilic two-Higgs-doublet models which can accommodate neutrino masses by means of the tiny vev of the additional Higgs dou-JHEP12(2015)160 blet. The models studied here differ among themselves by the symmetry that forbids the couplings between neutrinos and the scalar which gets the electroweak scale vev. The cases studied here span a discrete Z 2 and a softly broken global U(1) symmetry.
The bounds considered come both from theory and experiment. The unitarity perturbative requirement at tree-level strongly constrains the scalar mass spectrum of these models, either by the presence of a very light neutral scalar (m H ∼ v 2 ), in the Z 2 model, or with a degeneracy between the scalar and pseudoscalar particle masses (m H = m A ), in the global U(1) scenario.
If there is no additional particle content, the Z 2 symmetry model was found to be in severe tension with the electroweak precision tests, due to the very light neutral scalar, which generates a large negative contribution to the S parameter. The inclusion of a Majorana mass term for the right-handed neutrinos, providing a low scale realization of the seesaw type I mechanism, does not save the model, as the right-handed neutrino contribution to the S parameter is always negligible. Therefore, we conclude that the neutrinophilic 2HDM with a spontaneously broken Z 2 symmetry is strongly disfavored by data.
The analysis of the model with an explicit broken global U(1) symmetry reveals a region of the parameter space which is allowed by all bounds considered. Due to the set of constraints and the symmetries of the model itself, the spectrum is quite limited. The U(1) symmetry predicts that the neutrinophilic scalar is degenerate in mass with the pseudoscalar, m H = m A . Besides, the electroweak precision tests play a very important role, specially the T parameter which is sensitive to the absolute mass splitting of the pseudoscalar and the charged scalars, limiting it to be at most ∼ 80 GeV. Therefore, an important consequence of the theoretical and experimental constraints is that, if the new scalars are above ∼ 400 GeV, all these particles should have very similar masses. Moreover, the Z invisible width excludes the region m H = m A < m Z /2. Besides, the h → γγ branching fraction might be modified by about ±30% for m H ± < 200 GeV, while for heavier H ± , above 500 GeV, this ratio can be atmost 1 or lower by 5%. Finally, we stress that this model can be well within the reach of LHC 13 TeV, by probing the h → γγ branching fraction of by direct pair production of the charged scalars, if they are below O(300 GeV).