Multi-Higgs doublet models: physical parametrization, sum rules and unitarity bounds

If the scalar sector of the Standard Model is non-minimal, one might expect multiple generations of the hypercharge-1/2 scalar doublet analogous to the generational structure of the fermions. In this work, we examine the structure of a Higgs sector consisting of N Higgs doublets (where N \geq 2). It is particularly convenient to work in the so-called charged Higgs basis, in which the neutral Higgs vacuum expectation value resides entirely in the first Higgs doublet, and the charged components of remaining N-1 Higgs doublets are mass-eigenstate fields. We elucidate the interactions of the gauge bosons with the physical Higgs scalars and the Goldstone bosons and show that they are determined by an Nx2N matrix. This matrix depends on (N-1)(2N-1) real parameters that are associated with the mixing of the neutral Higgs fields in the charged Higgs basis. Among these parameters, N-1 are unphysical (and can be removed by rephasing the physical charged Higgs fields), and the remaining 2(N-1)^2 parameters are physical. We also demonstrate a particularly simple form for the cubic interaction and some of the quartic interactions of the Goldstone bosons with the physical Higgs scalars. These results are applied in the derivation of Higgs coupling sum rules and tree-level unitarity bounds that restrict the size of the quartic scalar couplings. In particular, new applications to three Higgs doublet models with an order-4 CP symmetry and with a Z_3 symmetry, respectively, are presented.


Introduction
The discovery of the Higgs boson at the Large Hadron Collider (LHC) appears to complete the story of the Standard Model (SM) of particle physics [1,2]. In particular, the subsequent experimental measurements of the properties of the observed scalar with mass 125 GeV are so far consistent with those of the SM Higgs boson [3]. The current Higgs data set is still statistically limited, so a more precise statement is that the observed scalar behaves as a SM-like Higgs boson, to within an accuracy of about 20%. Future experimental studies of the Higgs boson at the LHC will continue to search for deviations from SM behavior, as well as evidence for additional scalar states that might comprise an extended Higgs sector.
Nevertheless, there are numerous reasons to suspect that there must exist new physical phenomena beyond the SM. For example, the SM cannot accommodate dark matter [4], massive neutrinos [5], baryogenesis [6] and the gravitational interaction. Indeed, there is no fundamental understanding of how the electroweak scale arises, and why this scale is many orders of magnitude smaller than the Planck scale [7]. However, the casual observer examining the structure of the SM in its present form might also be puzzled that the set of fundamental scalar fields consists of a single neutral CP-even Higgs boson. After all, the SM employs a direct product of three separate gauge groups and the elementary fermionic matter comes in three generations. Why should one expect a scalar sector to consist only of a single physical state?
In light of the non-minimal structure of the fermionic and gauge bosonic sectors of the SM, one is tempted to suppose that the scalar sector is likewise non-minimal. It is a simple matter to construct an extension of the SM that incorporates an enlarged Higgs sector. In order to preserve the tree-level relation, ρ 0 ≡ m W /m Z cos θ W = 1 (which is confirmed by the electroweak data after accounting for electroweak radiative corrections [8]), the electroweak quantum numbers of the Higgs scalar multiplets are constrained [9,10]. For example, a Higgs sector employing hypercharge-1/2 scalar doublets and hypercharge-zero scalar singlets yields ρ 0 = 1, independently of the vacuum expectation values (vevs) of the neutral scalar fields. Following the generational pattern of the fermions, we shall simply replicate the SM Higgs doublet and consider an extended Higgs sector consisting of N hypercharge-1/2 scalar doublets. For N ≥ 2, one must further require that at the minimum of the scalar potential, only the neutral components of the N scalar fields acquire vevs [11,12]. 1 In this paper, we wish to explore various relations among Higgs couplings and bounds on scalar masses that arise in an N Higgs doublet extension of the SM. Electroweak gauge invariance plays a central role in determining the structure of the Higgs couplings. Our primary focus here will be the couplings of Higgs bosons to gauge bosons and the cubic and quartic scalar self-couplings. The couplings of the Higgs bosons to fermions is governed by the Yukawa couplings, which must be highly constrained in order to avoid tree-level flavor-changing neutral currents mediated by neutral scalars [13,14]. In this paper, we shall simply postpone the consideration of the scalar-fermion interactions [15].
Our main tool for obtaining relations among Higgs couplings and constraints on Higgs masses is tree-level unitarity. If one computes the scattering amplitudes for 2 → 2 scattering of gauge and Higgs bosons, under the assumption that all Higgs couplings are independent of one another, then one finds that some of the scattering amplitudes grow with the center of mass energy. Such behavior is not consistent with unitarity. Of course, there is no paradox here since the assumption of independent Higgs couplings is incorrect. Electroweak gauge invariance imposes relations among the couplings that guarantee that the bad high energy behavior of any scattering amplitude must exactly cancel. One can turn this argument around and derive relations among the Higgs couplings that are required to cancel the bad high energy behavior of all scattering amplitudes [16,17]. This procedure allows one to deduce a variety of sum rules that relate various Higgs couplings [18][19][20].
Having canceled the bad high energy behavior, one finds that scattering amplitudes in the high energy limit either approach a constant value or vanish in the limit of large center of mass energy. In the former case, the condition of tree-level unitarity imposes an upper limit on the value of this constant. Ultimately, one can show that this constant is a function of dimensionless quartic couplings that appear in the scalar potential. Thus, the imposition of tree-level unitarity yields an upper bound on the values of various combinations of quartic scalar couplings. This in turn can provide upper bounds on some combinations of scalar masses [21,22].
In section 2 we consider the most general N Higgs doublet model (NHDM). We explicitly write out the couplings of the Higgs bosons to the gauge bosons, which arise from the scalar field kinetic energy terms after replacing the ordinary derivative with the gauge covariant derivatives of the electroweak theory. Remarkably, these couplings are controlled by an N × 2N matrix B, whose physical significance is explained below. We also note the appearance of the 2N × 2N matrix A = Im(B † B). The matrix A is an orthogonal antisymmetric matrix that is governed by N (N − 1) parameters. These parameters are independent of the basis of scalar fields used to define the model, and thus are related to physical observables. 2 The matrix B is governed by 2(N − 1) 2 physical parameters, which include the N (N − 1) parameters that already appear in the matrix A. In addition, the matrix B depends on N − 1 (unphysical) phases that can be eliminated by appropriately rephasing the physical charged Higgs fields of the model. We also examine some aspects of the scalar self-couplings. We find that some specific cubic and quartic couplings involving Goldstone boson fields also depend exclusively on the matrices A and B. This is a consequence of the fact that these interactions terms are related by gauge-fixing to terms appearing in the gauge-covariant scalar kinetic energy terms. The structure of the other cubic and quartic couplings are not as simple, and involve more complicated invariant expressions involving the coefficients of the scalar potential.
In section 3, we derive a variety a sum rules involving the couplings of gauge bosons and Higgs bosons and the couplings of Goldstone bosons and Higgs bosons. In section 4, we present an efficient technique to impose the condition of tree-level unitarity, leading to upper bounds on various combinations of couplings and scalar masses. We apply this to known results and present new results for the three Higgs doublet models (3HDMs) with Z 3 symmetry and with order-4 CP symmetry, respectively. Conclusions are given in section 5. In our analysis of the NHDM, one can define a new basis of scalar fields such that the neutral scalar vev resides entirely in one of the Higgs doublet fields, denoted by Φ H 1 . This is the well-known Higgs basis [23][24][25][26], in which the charged component of Φ H 1 is identified as the charged Goldstone boson field and the imaginary part of the neutral component of Φ H 1 is the neutral Goldstone boson field. The Higgs basis is not unique, since one is free to make an arbitrary U(N − 1) transformation among the remaining N − 1 doublet fields. In particular, one can employ this transformation to diagonalize the physical charged Higgs squared-mass matrix. This procedure yields the charged Higgs basis, as discussed in appendix A. Note that the resulting basis is unique up to an arbitrary separate rephasing of the N − 1 scalar doublets that contain the physical charged Higgs boson fields.
In appendix B, we apply the analysis of the NHDM given in section 2 to the two-Higgs doublet model (2HDM). We first discuss the complex two-Higgs doublet model (C2HDM), in which a Z 2 -symmetric scalar potential is softly broken by a complex squared-mass parameter. We then generalize to the most general 2HDM, which is treated using the basisindependent formalism of ref. [27]. In both cases, we display the explicit 2HDM forms for the matrices A and B and exhibit the unphysical parameters in the matrix B that can be eliminated by an appropriate rephasing of the physical charged Higgs fields. In appendix C, we demonstrate how to count the number of parameters that govern the matrices A and B of the NHDM and identify which of these parameters are physical.
In our derivation of coupling sum rules and unitarity relations, the coupling of the Goldstone bosons fields to the physical Higgs fields play an important role. Although the initial forms of these expressions are quite complicated, the corresponding cubic couplings end up reducing to remarkably simple forms. As an example, we provide in appendix D the details of the derivation and simplification of the coupling of two neutral Goldstone fields and a physical Higgs scalar. Finally, in appendix E, we rederive the sum rules obtained directly from the NHDM interaction Lagrangian using an alternative method, which imposes the cancellation of bad high energy behavior in the 2 → 2 scattering amplitudes of processes involving the gauge and Higgs bosons. The relation between sum rules involving the neutral Higgs boson couplings to W + W − and ZZ is clarified in appendix F.

N Higgs doublet models
In this section we discuss the bosonic Lagrangian of the most general NHDM. The field content consists of the SU (2) L ×U (1) Y gauge bosons and N hypercharge-1/2 Higgs doublet fields, parameterized as, , for k = 1, . . . , N . (2.1)

The scalar potential
For the scalar potential, we follow the notation of [25,28]: where, by hermiticity, µ ij = µ * ji , λ ij,kl ≡ λ kl,ij = λ * ji,lk . (2.3) Using eq. (2.1), the Higgs potential becomes where, 3 and is the mass matrix for the charged scalar fields. Requiring the absence of linear terms yields the stationarity condition, The vacuum of these models has been studied in ref. [11]; here we assume only that the electromagnetic U (1) em remains unbroken. 4 Expanding the neutral fields in terms of their real and imaginary components, the second and third terms of V 2 may be written as

12)
3 It is useful to note that, because (M 2 ± )ij and λ ik,lj v k v * l are hermitian in (i, j), the second term in eq. (2.7) may be written as where [cf. footnote 3], Using eq. (2.3), we conclude that M 2 ± is hermitian, the matrices M 2 R and M 2 I are real and symmetric, and M 2 RI is a general real matrix. Thus, the mass matrix in eq. (2.12) is real and symmetric, as expected. Our mass matrices agree with those in ref. [11], after noting that their ν d k = v k / √ 2. Our eq. (2.15) differs by a minus sign in the last term with respect to a similar eq. (A17) of ref. [32]. However, this sign error is the result of a misprint, in light of the agreement in signs between our results and eqs. (A18)-(A22) of ref. [32].
For later use, we note that Eqs. (2.5)-(2.9) are written in with respect to a generic basis of the scalar doublet fields. One can define a new set of charged and neutral scalar fields denoted, respectively, by S ± a (a = 1, . . . , N ) and S 0 β (β = 1 . . . , 2N ), via Note that the transformation given in eq. (2.20) results in a real orthogonal similarity transformation of the 2N × 2N symmetric squared-mass matrix given in eq. (2.12). That isṼ is a 2N × 2N real orthogonal matrix. As a result, we find Hence, it follows that, Note that we have used the convenient notation of refs. [33,34], which in turn was inspired by refs. [32,35]. In addition, Finally, one can show that the matrix Im V † V is antisymmetric. Moreover, using eqs. (2.23), this matrix satisfies, (2.26) Thus, Im V † V is the only nontrivial piece of V V † and V † V . As we shall see below, it has the crucial role of controlling scalar couplings involving the Z boson or its corresponding Goldstone boson. Combining the antisymmetry with eq. (2.26), we conclude that Im V † V is orthogonal In particular, all of its entries satisfy This will be of interest later. Let us consider matrices such that where we have definedv With these choices, the unitarity of U implies where the last equality holds because Im(V † V ) αβ is antisymmetric in αβ. We wish to study the squared-mass matrix of the charged scalars with respect to the charged scalar fields S + a , This means that with respect to the charged scalar fields S + a [reached by transformations with (2.29)], the first row and first column of the transformed squared-mass matrix of the charged scalars vanishes. This identifies S ± 1 with the charged would-be Goldstone boson G ± , Next, we turn to the squared-mass matrix of the neutral scalar fields.
Using eq. (2.30), we start by looking at For the first equality, we have used eqs. (2.13) and (2.15). To reach the last line of eq. (2.40), we have used eq. (2.31) and the stationarity condition given in eq. (2.11). Similarly, Multiplying eq. (2.40) by (ReV T ) αi , multiplying eq. (2.41) by (ImV T ) αi , and summing over i, we conclude from eq. (2.39) that where the first equality holds since M 2 N is a real symmetric matrix. This means that with respect to the scalar fields S 0 β , [reached by transformations with (2.30)], the first row and first column of the transformed squared-mass matrix of the neutral scalars vanishes. This identifies S 0 1 with the neutral would-be Goldstone boson G 0 , One can choose matrices U and V in such a way that the transformations in eqs. (2.19) and (2.20) yield the charged and neutral scalar mass eigenstate fields, respectively, Since we have identified S ± 1 = G ± and S 0 1 = G 0 , it follows from our above analysis that the matrices U and V must satisfy eqs. (2.29) and (2.30), respectively. In this case, S ± a (a = 2, . . . , N ) and S 0 β (β = 2, . . . , 2N ) denote the fields of the physical charged and neutral scalar particles, respectively. Their corresponding masses are m 2 ±,k (k = 1, 2, . . . , N ) and m 2 β (β = 1, 2, . . . , 2N ).
These are the only combinations of quartic couplings (and vevs) that one can obtain from the diagonalization of the scalar squared-mass matrices. Thus, only those cubic and quartic terms of the scalar potential involving these combinations will be related to scalar masses. Eqs. (2.46)-(2.47) constitute a crucial result of our paper, since, they will enable us to relate the gauge-Higgs couplings with the scalar-scalar couplings. 6

Gauge-Higgs couplings
When expressed in terms of the physical gauge fields, the gauge covariant derivative may be written as where g is the SU (2) coupling constant, c W = cos θ W , s W = sin θ W , e is the electric charge of the positron, Q is the charge operator, and 7 when acting on SU(2) doublets. Note that the signs of the coupling constants and gauge fields above correspond to choosing all the η k equal to +1 in the notation of ref. [36]. This coincides with the conventions of ref. [9], but differs in the signs in g from refs. [37,38]. 8 This also has an impact on any Feynman rules proportional to M W or M Z . The kinetic term for the scalar fields is with c 2W = cos (2θ W ) and

55)
are matrices of dimension N × 2N and 2N × 2N , respectively. 8 The signs in refs. [25,33] correspond to yet another choice, which yields an unexpected sign in the relation tan θW = −g ′ /g. 9 We employ the notation where S + The matrix A is basis-independent and hence physical, whereas the matrix B is basisindependent up to unphysical phases that can be absorbed into the definition of the physical charged Higgs fields, S ± a (for a = 2, 3, . . . , N ). Further details will be provided in the next section. Eqs. (2.51)-(2.54) agree with eqs. (29a)-(29p) of ref. [33], if we notice that, because of the different sign in the coupling g, we have g = −g GLOO , M W = −M GLOO W , and M Z = −M GLOO Z , where the subscript "GLOO" stands for ref. [33]. The only exception is in the ZS γ S β coupling given in eq. (2.53) which, after the difference in notation is properly accounted for, disagrees with the sign of eq. (29h) of ref. [33]. We have checked in both notations that their incorrect sign can be attributed to a misprint.

The A and B matrices
In section 2.2 we showed that the couplings arising from the kinetic Lagrangian depend exclusively on the matrix B in the charged scalar sector and on A = Im(B † B) in the neutral scalar sector. A and B are defined in terms of the matrices U and V , which relate the charged and neutral fields in a generic basis of scalar fields, {Φ k }, to the corresponding mass-eigenstate scalar fields, respectively. The choice of basis is of course arbitrary. For example, another set of scalar fields where X is some N × N unitary matrix, could have been employed. The scalar field kinetic energy terms are invariant with respect to eq. (2.57), since which follows from X † X = 1 N ×N . Consequently, the interactions of the scalars with the gauge bosons given by eqs. (2.51)-(2.54) are basis-independent. Indeed, any physical observable cannot depend on the choice of basis. We would like to use these observations to address the behavior of the matrices A and B under an arbitrary change of basis. The interaction Lagrangian given by eqs. (2.51)-(2.54) is written in terms of the scalar mass eigenstates S ± a and S 0 β , which are related via eqs. (2.19) and (2.20) to the scalar fields in a generic basis. The diagonalization of the neutral scalar squared-mass matrix is given by eq. (2.39) and yields real neutral scalar mass-eigenstate fields. The overall sign of the neutral scalar mass-eigenstate fields are not physical. However, the standard practice is to fix this sign by appropriately restricting the range of the angles that parameterize the diagonalization matrix. Having adopted this convention where the sign of the neutral scalar mass-eigenstate fields are fixed, it follows that the matrix A that appears in eqs. (2.52) and (2.53) is basis-independent and hence physical.
In contrast, the diagonalization of the charged scalar squared mass matrix yields complex mass-eigenstate charged scalar fields. By convention, the phase of the charged Goldstone field is fixed. In particular, it is convenient to choose X = U in eq. (2.57), which yields the scalar field basis, where S + 2 , . . . , S + N are the physical (mass-eigenstate) charged Higgs fields with corresponding masses m 2 ±,i . This is called the charged Higgs basis and has two defining properties: 1. S + 2 , S + 3 , . . ., S + N , are the charged scalar mass-eigenstate fields; 2. the first doublet field, Φ C 1 , has the massless would-be Goldstone boson G + as its charged (upper) component, and its phase has been chosen such that the (real and positive) vev, v ≃ 246 GeV, is contained in the real part of its neutral (lower) component.
The Higgs basis (which is defined by property 2 alone) and the charged Higgs basis are discussed in detail in appendix A. Notice that the fields H 0 and ϕ C0 2 , . . . , ϕ C0 N are not the neutral scalar mass-eigenstate fields. B is the matrix that transforms these fields, written in the charged Higgs basis, into the neutral scalar mass-eigenstates.
The charged Higgs basis is unique up to a possible rephasing of the charged Higgs fields, S + a → e iχa S + a [for a = 2, 3, . . . , N ]. That is, the charged Higgs basis is a family of scalar bases that is characterized by N − 1 phases, χ a , as discussed at the end of appendix A. In eqs. (2.53) and (2.54), the invariant combinations B aβ S − a and its charged conjugate appear. Hence, it follows that the matrix B is not quite basis-independent, since B aβ → e iχa B aβ [for a = 2, 3, . . . , N ] under the rephasing of the charged Higgs fields to preserve the invariance of eqs. (2.53) and (2.54).
In contrast, the mixing matrices U and V that relate the charged and neutral fields in a generic basis of scalar fields, {Φ k }, to the corresponding mass-eigenstate scalar fields, respectively, are basis dependent. It is instructive to examine the question of basis dependence and work out the explicit forms of the matrices A and B in the more familiar 2HDM. In this case, the matrix U contains the angle β = tan −1 (v 2 /v 1 ). As discussed at length in ref. [27], this means that the angle β does not in general have a physical meaning, since the vevs v 1 and v 2 (and hence tan β) transform under the basis transformation given by eq. (2.57). Similarly, the mixing angle parameters appearing in the matrix V that transforms the neutral scalar fields of the generic basis into the neutral scalar masseigenstate fields are also basis dependent. To find basis-independent quantities, one must consider the neutral mixing angle parameters relative to the angle β. That is, under a basis transformation, both the neutral mixing angle parameters and the angle β shift by the same amount so that their difference is invariant. Thus, one way to determine the invariant neutral mixing angle parameters is to work in the Higgs basis, in which β = 0. Further details can be found in appendix B. The matrix B is the NHDM generalization of the neutral mixing angle parameters relative to the angle β. As noted above, the matrix B is almost basis-invariant, since unphysical phases remain that reflect the possible rephasing of the charged Higgs fields.
From section 2.1, we can deduce the following properties of the matrices A and B. First, it is convenient to re-express A in terms of the matrixṼ defined in eq. (2.21). We introduce the 2N × 2N orthogonal antisymmetric matrixJ, which in block form is defined byJ Then A = Im(V † V ) can be rewritten as, It immediately follows that A is a real orthogonal and antisymmetric matrix; i.e., In particular, the orthogonality of A implies that |A αβ | ≤ 1. Given an N × N unitary matrix U , it is always possible to represent this matrix by the following 2N × 2N real orthogonal matrix, Henceforth, we shall always identify the real orthogonal representation of a unitary matrix with a subscript R and a tilde (to indicate that the dimensionality of the matrix has been doubled). Using this notation, it is convenient to construct the 2N × 2N matrixB, SinceŨ R andṼ are real orthogonal 2N × 2N matrices, it follows thatB is also a real orthogonal 2N × 2N matrix. Moreover, Noting that B = U † V and A = Im(B † B), one can write Finally, we note that The central point of this section is the following. Unless there is an (imposed symmetry) reason to single out some specific basis, the best way to count parameters and to set up a numerical simulation is to write the potential in the charged Higgs basis ab initio.
As seen from appendix A, the charged Higgs basis is (almost) unique, up to the separate rephasing of the N − 1 doublets with zero vev. As a result, all the parameters of the scalar potential, when written in the charged Higgs basis, are either invariant or pseudo-invariant quantities with respect to arbitrary scalar basis transformations. Here, pseudo-invariant means invariant up to an overall phase that arises from the rephasing of the N − 1 doublets that contain the physical charged Higgs fields. It is straightforward to construct invariants from appropriate products of pseudo-invariants in which the overall phase ambiguity cancels. All such observable quantities are potential experimental observables. This provides the generalization of the parameters Y 1 , Y 2 , Y 3 and Z 1 , Z 2 , . . . , Z 7 championed for the 2HDM in ref. [39].

Parameter counting
We begin by asking the following question. How many parameters govern the matrices A and B and how many of these parameters are physical? To address this question, we first examine the matrices V and U , which transform the neutral and charged scalar fields of a generic basis to the corresponding scalar-mass eigenstates, respectively.
Starting from a scalar basis, {Φ k }, one can compute the neutral scalar squared-mass matrix as shown in eq. (2.45). The real orthogonal 2N × 2N diagonalizing matrixṼ is related to the matrix V defined in eq. (2.20) via eq. (2.21). With respect to a new scalar basis {Φ ′ ℓ }, one obtains a new matrix V given by Note that the same result can be obtained by employing eqs. (2.61) and (2.69), since after noting thatX RJX T R =J (the latter makes use of the fact that X is unitary). That is, A is basis-independent. In particular, if we choose X = W in eq. (2.71), thenṼ ′ =R c , which can be expressed in terms of N (N −1) parameters. 10 Consequently, eq. (2.73) implies that A =R T cJRc , which depends on N (N − 1) physical parameters. 11 Likewise, starting from a scalar basis, {Φ k }, one can compute the charged scalar squared-mass matrix as shown in eq. (2.44). The N × N unitary diagonalizing matrix U defined in eq. (2.19) depends on the basis. One also must take into account that the charged Higgs basis is not uniquely defined, as discussed at the end of appendix A. Hence, the physical charged Higgs fields may acquire phases under the basis change. If we perform a basis transformation given by eq. (2.57), we must allow for the possibility that S + a → e iχa S + a [for a = 2, 3, . . . , N ]. We can write the latter as where the primes indicate quantities associated with the transformed basis and That is, one obtains a new unitary diagonalizing matrix U ′ given by where there is an implicit sum over the repeated index k.
We can now compute the transformation of B under a change of scalar basis. With respect to the basis {Φ ′ ℓ }, we have after making use of eqs. (2.68) and (2.78). If we consider the charged Higgs basis where X = U , then eq. (2.77) yields Table 1. Number of parameters in the Y and Z coefficients of the Higgs potential.
after making use of eqs. (2.64) and (2.71) with X = U . The N − 1 phases χ a (for a = 2, 3, . . . , N ) that appear U ′ = K (and similarly are contained inŨ ′ T R ) are unphysical and can be absorbed into the definition of the charged Higgs fields S + a . Since W = U (unless the neutral Higgs fields are mass eigenstates in the charged Higgs basis), it follows thatB contains additional parameters beyond the N (N − 1) physical parameters that determine the matrixR c . As shown in appendix C.2, the matrix B (andB) depends on an additional (N − 1)(N − 2) physical parameters. That is, after absorbing the N − 1 phases into a redefinition of the charged Higgs fields, there are 2(N − 1) 2 physical parameters remaining in the matrix B (andB).
The C2HDM provides an interesting example, discussed in detail in appendix B. One sees in eq. (B.49) thatB (or equivalently B) depends on 3 angles, where one of the angles is unphysical and corresponds to the freedom to rephase the second scalar doublet with zero vev. In contrast, A in eq. (B.56) depends only on the two invariant angles contained in B. The case of N = 2 is special in that the parameters that define the matrix A corresponding precisely to the physical parameters that appear in the matrix B.
The 2(N − 1) 2 physical parameters contained in the matrix B are sufficient to parameterize all the gauge boson-Higgs boson interactions. But, the three-scalar and four-scalar interactions derived from the scalar potential necessarily involve additional parameters. We have already emphasized that the charged Higgs basis is especially useful in identifying the invariant (and pseudo-invariant) scalar self-coupling coefficients. In particular, in the charged Higgs basis, the scalar potential is given by The number of relevant parameters is shown in table 1.
The number of parameters gets reduced for two reasons. First, our definition of the charged Higgs basis allows for a rephasing of the N − 1 scalar doublets with zero vevs. This reduces the number of phases by N − 1. In addition, the stationarity conditions written in the charged Higgs basis relate some Y and Z parameters. More generally, eq. (A.16) can be used to obtain a relation of the Y and Z parameters with the charged scalar masses, which can be used to trade in the Y ij for the Z parameters and the charged scalar masses. Hence, using the charged Higgs basis as parameters, we need only the 1 4 N 2 (N 2 + 3) magnitudes and the 1 4 For example, in the 2HDM and in the notation of ref. [39], we find The seven magnitudes correspond to |Z 1 |, |Z 2 |, . . . , |Z 7 | and the two independent phases are Im(Z * 5 Z 2 6 ) and Im(Z * 5 Z 2 7 ), first identified in ref. [24] as basis invariant measures of CP violation. Although not independent in the case Z 5 , Z 6 , Z 7 = 0, the possibility of Z 5 = 0 is only covered by considering in addition the third invariant phase, Im(Z 6 Z * 7 ) [24,26]. Many different parameter choices exist in the literature. For example, one can employ: (i) m ±a and the quartic coefficients Z of the scalar potential in the charged Higgs basis; or (ii) the quadratic coefficients Y not fixed by the scalar potential minimum conditions and all of the quartic coefficients Z; or (iii) the physical parameters in B, the physical charged and neutral scalar masses, and (if needed) some invariant combination of scalar selfcouplings. In extended Higgs sectors with additional symmetries (where the basis in which the symmetries are manifest becomes physical), one can also employ the various ratios of scalar vevs in the list of parameters. For example, in the CP-conserving 2HDM with a softly-broken Z 2 -symmetric scalar potential, some authors choose parameters consisting of β = tan −1 (v 2 /v 1 ), β − α [which appears in the matrix B in our notation], the charged scalar mass m ± , the neutral scalar masses (m h , m H , and m A ), and the soft Z 2 breaking squared-mass term m 2 12 .

Scalar self-couplings
In this section, we demonstrate that some of the scalar self-couplings are related to the kinetic terms obtained in section 2.2. In particular, it is convenient to enquire which scalar couplings may be written exclusively in terms of the matrix B (and A), that appear in the kinetic terms, and the scalar masses m ±a and m β . We first consider the cubic scalar self-couplings. The scalar potential in eqs. (2.4)-(2.9) can be expressed in terms of the physical fields using the mixing matrices U and V in eqs. (2.19)-(2.20). The cubic vertex interactions can therefore be written as In contrast to the terms of the interaction Lagrangian that couple the scalar and vector bosons given in eqs. (2.51)-(2.54), additional structures appear beyond those combinations of U and V that define the matrices A and B. However, if we focus on the cubic couplings that involve at least one Goldstone field, the form of the cubic interaction terms simplify significantly and can be expressed in terms of the A and B matrices and the squared masses of the physical scalars. 12 Indeed, this is to be expected, as these interaction terms are related by gauge-fixing to the pure gauge boson terms arising from the kinetic Lagrangian, which can be expressed in terms of the matrices A and B, as shown in eqs. (2.51)-(2.54).
In order to simplify the non-vanishing cubic potential terms we employ eqs. (2.22)-(2.35), from which we obtain the further useful relations and Moreover, we must use the crucial relations given in eqs. (2.46)-(2.47), relating some combinations of quartic couplings with the scalar masses. This simplifies considerably the final expressions. The end result is, 13 where there is an implicit sum over repeated indices. Achieving the simplified forms of the couplings above is rather laborious. For example, the explicit derivation of eq. (2.93) is given in appendix D.
The cubic couplings that are present in the kinetic lagrangian are also present in parallel with the scalar potential. A comparison is shown in table 2, where we have collected the relevant Feynman rules. We notice that both sets of couplings depend on the same parameters: (2.94) Kinetic lagrangian

Scalar potential Coupling
Feynman rule Coupling Feynman rule This is not surprising, since gauge boson couplings and Goldstone boson couplings are related by the gauge-fixing, as previously noted. Although the equivalence theorem [41] is not a requirement on the couplings but rather on full processes [42], the relation between these couplings insures that the equivalence theorem is satisfied. For the calculation of the quartic couplings we use the general expression It is instructive to focus on the quartic couplings that involve the Goldstone boson fields. After some very long simplifications, we find the non-vanishing Goldstone-scalar quartic couplings that involve an even number of Goldstone boson fields listed below. : The corresponding expressions involving an odd number of Goldstone fields are more complicated and will not be given here.
In contrast to the cubic couplings given in eqs. (2.90)-(2.93), additional structures beyond the matrices A, B and the squared masses of the physical scalars appear in the expressions above. For example, the quadratic coefficient of the scalar potential in the charged Higgs basis, Y , appears on the right hand side of eq. (2.97). Note that there is no simple relation between quartic couplings in the kinetic Lagrangian and quartic terms from the scalar potential. For example, as in the SM, there is a (G 0 ) 4 coupling, but no Z 4 coupling. This does not violate the equivalence theorem since the processes involving quartic couplings typically involve other diagrams. In particular, only the sum of all contributing diagrams must obey the equivalence principle.

Coupling relations and sum rules from the Lagrangian
One can obtain a large number of sum rules from the kinetic part of the Lagrangian and the scalar potential of the NHDM. In the case of the 2HDM, two of the sum rules that are usually exhibited (see, e.g., refs. [19,20,43,44]) are where k identifies some specific neutral scalar physical field. It is easy to find the corresponding sum rules in the general NHDM: where the indices follow the same notation as above. In this section we use a simplified notation for the couplings, strictly related to the matrices A and B, in which some coupling [X a Y b Z c ] is identified as the term in the Lagrangian that depends explicitly on family type indices. For example, in the Lagrangian term . In cases where the corresponding coupling also exists in the SM, this procedure means simply that we have divided out by the SM coupling C 1 .
In addition to these, we have found many sum rules for an arbitrarily extended scalar doublet sector. For example, Using this result in eq. (3.3), and recalling that Further equations arising from the kinetic Lagrangian are The relations among different couplings in the kinetic Lagrangian imply that one can write eq. (3.9) in terms of massive fields (i.e., without the Goldstone bosons) as We have checked that all our relations are verified when using the couplings in the special case of the C2HDM. Notice that eqs. (3.5) and (3.8) are not proper sum rules, but rather relations among couplings valid in our particular model. These sum rules have been found by employing the NHDM interaction Lagrangian. But we know that sum rules can also be found for generic Lagrangians by using unitarity arguments, as in refs. [19,20]. The sum rules we have derived in eqs. (3.2) and (3.6)-(3.8) can also be found in that fashion. We will revisit this question in section 3.2 and in appendix E.
In the quartic couplings sector, we have found further interesting sum rules. For example, using the couplings in eq. (2.97), we find 14 which has the peculiar feature that a sum of couplings involving solely charged particles yields a contribution that depends on the masses of the neutral fields. Similarly, using the couplings in eq. (2.102), we find We can recombine a number of these sum rules as 14) or N a=1 because the following equation is satisfied, The sum rules involving quartic couplings presented in eqs. 14 Note that tr[µ] = tr[Y ] is a basis-invariant quantity.

Sum rules from unitarity
The constraints from unitarity have had an historical impact on the development of particle physics. The idea is that the scattering among vector bosons and/or scalars cannot grow with energy and must obey the optical theorem. Imagine that one writes the most general effective couplings between these states, up to dimension four. Forcing the sum of all terms growing like the fourth power of the center of mass energy (E) to vanish immediately restricts the vector bosons to originate from a gauge theory [16,17]. Requiring that all E 2 terms vanish forces the scalar-gauge couplings to originate from a gauge theory and further constrains the couplings [18,19]. But unitarity also limits the value of the constant in the E 0 terms. This has been used in the past to place limits on (combinations) of the scalar masses and couplings [18,[20][21][22].
In the previous section we derived the coupling relations (3.5) and (3.8) and many sum rules directly from the Lagrangian of the NHDM. We have checked that most sum rules involving triple couplings can be obtained from those presented in ref. [19], by applying them to the NHDM and cycling through all possible indices. The exception arises in applying eqs. (3.2) and (3.3) for the case of V V = ZZ, which are presented in ref. [19] under the assumption of CP conservation. For completeness, we include in appendix E a careful and detailed derivation of the sum rules given in eqs. (2.4)-(2.6) of ref. [19]. The derivation of these sum rules does not make use of CP conservation. Nevertheless, when applied to models with scalars in their section IV, the authors of ref. [19] focus on models that conserve CP. In particular, they point out that requiring CP conservation, in the case of a CP conserving Higgs sector. In contrast, the general NHDM may (or not) violate CP; nevertheless, the same sum rules apply. We include in appendix F a proof that , which makes no assumption concerning CP conservation, and yet is valid for scalars in any representation of SU (2) L , provided that the theory satisfies m 2 W = m 2 Z c 2 W without fine tuning of the various vevs.

Tree-level unitarity bounds
In this section, we present an algorithm for the determination of the tree-level unitarity bounds. We check that it reproduces results available in the literature, and present its application to two new cases with three Higgs doublets. Finally, the techniques presented here can be used for faster numerical implementation of unitarity bounds in more complicated cases, not amenable to closed-form solutions. The problem of finding the constraints imposed by tree-level unitarity has been addressed in the case of the 2HDM, both with Z 2 symmetry [45] and also for most general case [46,47]. For the 2HDM with a Z 2 symmetric scalar potential, the results are simple, but for the most general case one has to compute the eigenvalues of 4 × 4 and 3 × 3 matrices [46,47]. For the case of the 3HDM a solution is known for the case of S 3 symmetry [48] where because of the symmetry the solutions are again simple, although the method to obtain them is already quite complicated.
Clearly, for the NHDM one needs an algorithm that can be easily implemented numerically. As explained in the references above, since one is interested in the high energy limit, one just needs to evaluate the scattering S-matrix for the two body scalar bosons, and these arise exclusively from the quartic part of the potential, V 4 . Then, the first part of the problem consists in finding the set of two body states that can contribute. For the cases of low N and high symmetry we can choose conveniently the sets to take advantage of this [45,46], but if we are going to solve the problem numerically it is better to have a simple algorithm of general applicability. Since the electric charge and the hypercharge are conserved in this high energy scattering, we can separate the states according to these quantum numbers. For this it is better not to separate the real and imaginary parts of the neutral components, using the following notation for Higgs doublets The relevant two body states are given in the entries of table 3, and their complex conjugates.
Q 2Y State Number of states Table 3. List of two body scalar states separated by (Q, Y).
As an example, for N = 2 we have,

4)
S 0 α ={n 1 n 1 , n 1 n 2 , n 2 n 2 }, (4.5) 16π a + 0,Y=1 αβ where we have indicated on the right-hand side the dimensionality of the resulting matrices. The compound nature of the index α should be taken in account. For instance where the set {i, j} corresponds to α and the normalization N ij is 1/ √ 2 for the 2 body states with equal particles and 1 in all the other cases. This factor can be understood in the following way. When we take the derivative with respect to the 2 body state with equal particles we should divide by the normalization 1/ √ 2. But on the right-hand side we are taking derivatives with respect to the individual fields. To avoid double counting we should divide by two in the case of identical particles, so 1 2 ( 1 For illustration, consider the case of the 2HDM. The quartic part of the potential contains Using the procedure described above one can easily get for M ++ 2 , which coincides with the result of ref. [47], up to an interchange of rows and columns. We have applied this procedure for the known cases in the literature and, to illustrate the power of the method, we also present two new cases.

An (almost) trivial case: the Standard Model
Consider the Standard Model with n H = 1. With the conventions of ref. [36], we have In this case the M matrices reduce to There are therefore two independent eigenvalues, Applying partial wave unitarity to the eigenvalues of the s-wave amplitude matrix, which is a consequence of the optical theorem, Since the right hand side of eq. (4.19) is bounded by 1 4 , it follows that [49], which translates into This in turn implies a well known result [21,22].

Other known cases
We applied the method to the complex 2HDM with Z 2 symmetry and we recover the results of ref. [45,46]. Next we studied the general complex 2HDM and we also agree with the results of ref. [46,47]. Both these cases are for N = 2. A more complicated case is the 3HDM with S 3 symmetry. The matrices are larger but we were able to recover all the results of ref. [48]. This makes us confident that the procedure can be generalized to cases where the unitarity constraints have not been studied. We consider two new cases.

Unitarity bounds for specific processes in the NHDM
In the previous sections, we have presented an algorithm which computes tree-level unitarity bounds on a given chosen model of scalar doublets. The method of Lee, Quigg and Thacker [21,22], which we have generalized and optimized, yields necessary and sufficient conditions for tree unitarity in the scalar and gauge sectors of any NHDM. A possible shortcoming of this method is the computational time needed to find all the eigenvalues for each point in parameter space in the case of a general scalar doublet theory, i.e. without symmetries, or with a large N .
A complementary approach of finding necessary, although not sufficient conditions for every NHDM can be attained by using unitarity bounds on specific processes. Therefore, we can compute the partial-wave coefficient a 0 for gauge-gauge and gauge-scalar scatterings, using that from the optical theorem | Re(a 0 )| < 1/2, and that we can use the Equivalence Theorem to simplify calculations. The Equivalence Theorem allows us to perform all calculations from the scalar potential, and therefore, we can always define where M stands for the amplitude of the process. Another important property of this method, also used in the previous section, is that we will only consider the quartic couplings contribution at high energy. We use this simple approach to work out some examples in the NHDM for an arbitrary N .

The
As a first example, we consider the process with amplitude M (W + W − → W + W − ), which we approximate at high energies to M (G + G − → G + G − ). The leading order contribution for the process will be the quartic coupling with amplitude where we have used eq. (2.96) and the notation from section 3. It is straightforward to see that this process enforces the bound It is interesting to note that in the alignment limit of the coupling [hV V ] NHDM = [hV V ] SM we cannot further constrain any masses of new neutral scalars. This observation is due to the bound in the A matrix itself, which is orthogonal, implying that β [A 1β ] 2 = 1. This bound is given as a special case of eq. (2.10c) of ref. [18], and it is valid in any NHDM.

The ZS +
a → ZS + a process We now turn to the unitarity bounds arising from the ZS + a → ZS + a scattering. Using the same reasoning as before, we use the quartic coupling of eq. (2.103). We then have where Y is the quadratic parameter of the lagrangian in the charged Higgs basis. Using the optical theorem we find that Although one cannot predict mass bounds for the charged scalars in this process, it is possible to study numerically the dependence on a given choice of Y . This parameter in the charged Higgs basis is a physical one and can, in principle, be measured.

4.5.3
The W + S − a → W + S − a process As a final example we compute the unitarity bounds for the process W + S − a → W + S − a . We follow the same arguments as before and write the quartic coupling in eq. (2.97) as where B aβ is the coupling to [S 0 β W + S − a ] and (B † ) βa is the coupling to [S 0 β W − S + a ]. We can therefore use the optical theorem to obtain, It is convenient to rewrite eq. (4.50) as, so that the dependence on cubic couplings is explicit.

Conclusions
We have studied the most general SU (2) L × U (1) Y theory with N Higgs doublets. We have stressed the importance of the charged Higgs basis, where the magnitudes and basisinvariant combinations of phases of its scalar couplings Y and Z are observables. The kinetic Lagrangian depends exclusively on a N ×2N matrix B, or equivalently the 2N ×2N real orthogonal matrixB as defined in eq. (2.64), which governs the relation between the neutral scalar components of the scalar doublets in the charged Higgs basis and the neutral Higgs mass eigenstates. The matrix B (orB) depends on N − 1 unphysical phases (corresponding to the non-uniqueness of the charged Higgs basis) and on 2(N −1) 2 physical parameters. Of these, N (N − 1) appear in the special combination A = Im(B † B).
Although new parameters beyond B appear in the scalar potential, many couplings involving the Goldstone bosons G 0 and G ± can be related to couplings involving Z 0 and W ± , as expected from consistency with gauge fixing. We use the crucial eqs. (2.46)- (2.47) to show that such relations indeed hold. This is also consistent with bounds from unitarity, which are discussed in great detail. In particular, we develop an efficient algorithm for the inclusion of such bounds in NHDM and employ it in two new 3HDM models with a Z 3 and with a order-4 CP symmetry, respectively. Some model independent necessary constraints are shown, by applying the optical theorem to selected processes.
In models where the scalar potential exhibits additional symmetries, some new parameters may appear to arise, relating the original basis to the charged Higgs basis (thereby acquiring physical significance). For example, such a case arises in the the 2HDM with a Z 2 -symmetric scalar potential, where β is the angle that rotates the basis in which the Z 2 symmetry is manifest into the charged Higgs basis. But, such parameters can always be re-expressed in terms of those discussed here [e.g., as shown for β in eqs. (51) and (52) of ref. [39]]. New parameters do arise when fermions are included, which will be addressed elsewhere [15].

Acknowledgments
We are very grateful to A. Barroso

Appendices A The Higgs basis and the charged Higgs basis
In section 2.1, we have discussed how to determine the physical charged and neutral scalar mass eigenstates starting from a generic basis of scalar fields {Φ k }. The charged and neutral components of the kth scalar doublet can be expressed as linear combinations of the physical charged and neutral scalar mass eigenstates [cf. eqs. (2.19) and (2.20), respectively], where the corresponding coefficients define the matrices U and V . To find the physical couplings, we first determine the scalar potential minimum conditions in the original (generic) basis (these stationarity conditions given in eq. (2.11) are complicated coupled cubic equations in the vevs, v k ). Substituting these conditions back into eqs. (2.7)-(2.9) and re-expressing the generic basis scalar fields in terms of the physical charged and neutral scalar mass eigenstates yields the desired expressions for the physical couplings.
Note that in determining the charged and neutral scalar mass eigenstates via (2.19) and (2.20), one must decompose the scalar doublets into their charged and neutral components, treating each component separately. This hides an important characteristic of the Higgs potential-namely, physical observables must be invariant under a unitary transformation among the N scalar doublets. This is known as basis invariance, which is discussed at length in refs. [25,26,28].
It is often more convenient to perform the analysis of the NHDM by first transforming from the generic basis to the so-called Higgs basis [23][24][25][26], in which the neutral vev resides entirely in the first scalar doublet. This is achieved through a unitary transformation X, Since X is unitary, we can invert eq. (A.1) to obtain We employ the normalization of the vevs such that Φ 0 j ≡ v j / √ 2. Hence, taking the vacuum expectation value of both sides of eq. (A.3) and making use of eq. (A.2) yields  (2.20) and eq. (A.5) is that, in the latter, we have also transformed another combination of the real components of Φ k , which yields a real field H 0 that is not a mass eigenstate. Thus, we may parameterize One particularly nice feature of the Higgs basis is the simplification obtained in the stationarity conditions of eq. (2.11), and in the masses of eq. (2.10) and eqs. (2.13)-(2.15) [11]: Any matrix X obeying eq. (A.2) will yield a Higgs basis where, by definition, only the first scalar doublet has a non-zero vev. Of the infinitely many choices for X that satisfy eq. (A.2), it is particularly interesting to consider the transformation X = U , where U is the unitary matrix defined in eq. (2.19) that yields the physical charged Higgs mass eigenstates [12]. Notice that we have left the scalar doublet structure intact; in this procedure, the neutral components transform as the charged components. Thus, in general, the corresponding neutral partner of the physical charged scalar will not be a mass eigenstate. We define such a basis as the charged Higgs basis, which is a subclass of the class of Higgs bases defined above. Now, eqs. (2.29) and (A.5) remain valid and we may parameterize The requirement that no physical observable can depend on the choice of χ provided the original motivation for the basis invariance considerations in the Higgs sector [24].
where H + 2 , . . . , H + N are the physical charged Higgs mass eigenstate fields, with corresponding masses m 2 ±,i . In this first step, the neutral components are transformed by  After this first step, eqs. (A.8)-(A.11) become Thus, In addition, 23) which implies that the N + 1 st row and column of the matrix,  vanish, corresponding to the massless Im ϕ C0 1 = G 0 . As a result, The first equality above is a consequence of eq. To determine the neutral scalar mass eigenstates, one must diagonalize the real orthogonal (2N − 1) × (2N − 1) squared-mass matrix that mixes the (2N − 1) neutral scalar fields H 0 , Re ϕ C0 2 , . . . , Re ϕ C0 N , Im ϕ C0 2 , . . . , Im ϕ C0 N (which are defined in the charged Higgs basis). This is achieved thorough the One can define the real orthogonal 2N × 2N matrixṼ C in analogy with eq. (2.21), which satisfy equations analogous to eqs. (2.22)-(2.24). Deleting the first row and column ofṼ C yields the matrix that diagonalizes the squared-mass matrix of the neutral scalars fields in the charged Higgs basis.
Performing the neutral scalar squared-mass diagonalization in the charged Higgs basis can be especially useful in some circumstances. Eqs. (2.46) and (2.47) assume very simple forms in the charged Higgs basis, As a further example, we notice that the cubic terms of the scalar potential in eq. (2.8) may be written in the charged Higgs basis as All couplings with two indices equal to 1 may be related with eqs. (A.21)-(A. 22), and hence can be related to the scalar masses. We find where there are implicit sums over repeated indices, i = 2, . . . , N and j = 2, . . . , N . As in the case of the Higgs basis, the charged Higgs basis is also not uniquely determined. Indeed, the charged Higgs basis is a class of bases, since starting in any given charged Higgs basis, one can separately rephase the N − 1 scalar doublet fields, Φ C 2 , Φ C 3 , . . . , Φ C N , while preserving the corresponding charged components as mass-eigenstate fields. It is convenient to keep track of this rephasing degree of freedom. Thus, in eq. (2.19), we will choose the unitary matrix U such that U k1 =v k and make some conventional choice for the overall phases of the column vectors U kj , for j = 2, 3, . . . , N . Such a choice picks out one of the possible charged Higgs bases. In this basis, the corresponding physical charged Higgs fields are denoted by S + a (for a = 2, 3, . . . , N ). One can of course transform to any other charged Higgs basis by an appropriate rephasing of Φ C 2 , Φ C 3 , . . . , Φ C N , in which case the corresponding physical charged Higgs fields are also rephased, S + a → e iχa S + a .

B Connections with the 2HDM
B.1 2HDM with a softly-broken Z 2 discrete symmetry The Lagrangian shown in sections 2.2 and 2.5 depends on the matrices U and V introduced in eqs. (2.19) and (2.20). This is not the notation commonly used in two Higgs doublet models (2HDM) [9,52]. Here we make the connection to the notation used in the complex two Higgs doublet model (C2HDM), where the Z 2 symmetry ( is softly broken by a complex squared-mass parameter [53][54][55][56]. In this case, one transforms to the Higgs basis [24,28] through In eq. (B.2), c β = cos β, s β = sin β, and v = v 2 1 + v 2 2 = ( √ 2G F ) −1/2 . Without loss of generality, we have taken the vevs v 1 and v 2 real. 16 The doublets in the Higgs basis may be written , where G ± and G 0 are the Goldstone bosons and H ± are the physical charged scalars. Thus, the matrix U of eqs. (2.19) is simply Let us parameterize the scalars Φ 1 and Φ 2 in the original generic basis as . (B.5)

Eqs. (B.1) and (B.3) yield for the massless would-be Goldstone boson
We define the orthogonal state The fields η 1 , η 2 , and η 3 combine into the mass eigenstates h 1 , h 2 , and h 3 as  where the orthogonal matrix may be parameterized as Here, s i = sin α i , c i = cos α i (i = 1, 2, 3), and, without loss of generality, the angles may be restricted to [57] − π/2 < α 1 ≤ π/2, −π/2 < α 2 ≤ π/2, 0 ≤ α 3 ≤ π/2. (B.9) By definition, we take the masses of the neutral scalars in increasing order: m 1 < m 2 < m 3 . We would like to recombine these expressions into the form of eqs. (2.1) and (2.21): which satisfies eqs. (2.67), as expected. As we have seen in sections 2.2 and 2.5, several important couplings, including the couplings of each neutral scalar with two vector bosons, involve the special antisymmetric combination (B.13) wheres 1 = sin (α 1 − β) andc 1 = cos (α 1 − β). Notice that, although three angles appear, there are in fact only two independent parameters in the most general orthogonal and antisymmetric 4 × 4 matrix. Indeed, such a matrix can always be parameterized as  wherec k = cos θ k ands k = sin θ k , for k = 1, 2. Thus, of the three angles in the matrix R, only two combinations can be determined by measurements involving solely the neutral scalars.
The complete set of Feynman rules for the C2HDM is presented on a webpage [58]. We have checked explicitly that the couplings in sections 2.2 and 2.5 reproduce the C2HDM Feynman rules in [58]. This is a highly non-trivial cross-check since the expressions are very complicated when written in terms of the angles α i and β. Moreover, for those couplings involving masses the equality is only obtained when using the relation which holds in the C2HDM [56]. The following relations are also useful: (B. 16) In the limit of the real 2HDM (s 2 → 0 and s 3 → 0), one finds It is also interesting to consider what happens in the SM, where Thus, V = (i 1), and confirming that Re(V † V ) = 1 2×2 and Im(V † V ) is antisymmetric.

B.2 Basis-independent treatment of the most general 2HDM
By using the basis-independent techniques introduced in refs. [26,27], one can analyze the most general CP-violating 2HDM (with no additional symmetries imposed on the scalar potential) in terms of quantities that are independent of the choice of basis for the two scalar doublet fields. All physical observables of the theory can be expressed in terms of such basis-independent quantities. It is instructive to see how this formalism is related to the treatment of the NHDM given in section 2.
Since the notation of ref. [27] differs somewhat from the notation used in this paper, we provide here a brief introduction to the basis-independent treatment of the 2HDM. We begin with the scalar potential in a generic basis given in eq. (2.2). The vevs of the scalar doublet fields are given by where v = 2m W /g ≃ 246 GeV andv = (v 1 ,v 2 ) is a complex vector of unit norm. A second unit vectorŵ can be defined that is orthogonal tov, where ǫ 12 = −ǫ 21 = 1 and ǫ 11 = ǫ 22 = 0. Indeed, the complex dot product,v * jŵ j = 0, where the sum over the repeated index j is implicit.
It is convenient to define two hermitian projection operators, Note thatv andŵ are eigenvectors of the matrix V . The matrices V and W can be used to define the following manifestly basis-invariant real quantities that depend on the scalar potential parameters defined in the generic basis of scalar fields [cf. eq. (2. 2)], In addition, we shall define the following pseudo-invariant (potentially complex) quantities, The significance of the quantities defined by eqs. (B.23)-(B.32) become clearer after rewriting the scalar potential in the Higgs basis. Using the notation of eqs. (B.20) and (B.21), the Higgs basis fields can be defined as, 17 In particular, note that the vevs of the Higgs basis fields are as required by eq. (A.5). That is, starting from the scalar potential defined in the generic basis [cf. eq. (2. 2)], we simply setv = (1, 0) andŵ = (0, 1). Applying these results to eqs. (B.23)-(B.32), we see that Y 1,2,3 are the coefficients of the squared mass terms and Z 1,2,...,7 are the coefficients of the quartic terms of the scalar potential when expressed in terms of the Higgs basis fields. In particular, In the Higgs basis, the minimization of the scalar potential yields One key observation is that the Higgs basis as defined by eq. (B.34) is unique only up to an overall phase redefinition of the Higgs basis field H 2 → e iχ H 2 . Indeed, in light of eq. (B.33), the phase freedom in defining the Higgs basis simply corresponds toŵ → e −iχŵ . It then follows that Y 3 , Z 5 , Z 6 and Z 7 also acquire a phase under H 2 → e iχ H 2 , which is why these quantities were called pseudo-invariants above. In contrast, Y 1 , Y 2 and Z 1,2,3,4 are invariant under H 2 → e iχ H 2 .
Our goal now is to evaluate the matrices B and A defined in eqs. (2.55) and (2.56), respectively. This requires the diagonalization of the charged Higgs and neutral Higgs squared-mass matrices. First, we consider the charged Higgs squared-mass matrix, which is given by eq. (2.10). It is convenient to rewrite this matrix as follows. Following ref. [27], we note that we can expand an hermitian second-ranked tensor in terms of the eigenvectors of V , Applying eq. (B.38) to the hermitian charged Higgs squared-mass matrix, it follows that The diagonalization of M 2 ± is straightforward. Defining the diagonalization matrix U as in eq. (2.19), it follows that which satisfies, That is, of the two eigenvalue of M 2 ± , one is zero, corresponding to the charged Goldstone boson, and the other is m 2 H ± = Y 2 + 1 2 Z 3 v 2 . Next, we obtain the matrix that diagonalizes the neutral Higgs squared-mass matrix. In this analysis, it will prove useful to first perform the diagonalization in the Higgs basis, since this allows us to easily identify the relevant basis-independent quantities. This has been carried out in ref. [27]. Here we summarize the main results that we need for our present analysis. The three physical neutral Higgs boson mass-eigenstates are determined by diagonalizing the 3 × 3 real symmetric squared-mass matrix that is defined in the Higgs basis [25,27], where Z 345 ≡ Z 3 + Z 4 + Re Z 5 . To identify the neutral Higgs mass-eigenstates, we diagonalize the squared-mass matrix M 2 . The diagonalization matrix is a 3 × 3 real orthogonal j q j1 q j2 0 i 0 1 c 12 c 13 −s 12 − ic 12 s 13 2 s 12 c 13 c 12 − is 12 s 13 3 s 13 ic 13 Table 4. Invariant combinations of the neutral Higgs boson mixing angles θ 12 and θ 13 , where c ij ≡ cos θ ij and s ij ≡ sin θ ij .
matrix that depends on three angles: θ 12 , θ 13 and θ 23 ,    As shown in ref. [27], the invariant angles θ 12 and θ 13 are in fact basis-independent quantities-that is, they can be expressed explicitly in terms of basis-independent combinations of quantities defined in any scalar field basis. 18 The neutral Goldstone boson and the physical neutral Higgs states (h 0 ≡ G 0 and h 1,2,3 , respectively) are then given by: where the q j1 and q j2 are invariant combinations of θ 12 and θ 13 , which are exhibited in Table 4. In particular, note that the quantities q j1 and q j2 are basis-invariants and the neutral fields h k are also invariant with respect to a rephasing of the Higgs basis field H 2 .
To identify the diagonalizing matrix V defined in eq. (2.20), 19 we make use of eq. (B.33) to rewrite eq. (B.46) as follows, for i = 1, 2. Hence, it immediately follows that It is straightforward to check that as noted in eq. (2.65).
We immediately see that B is not an invariant matrix in light of eq. (B.45). Nevertheless, in eq. (2.51) we note that the matrix B always appears along with the charged Higgs or Goldstone fields, namely B aβ S − a (and its hermitian conjugate). Recall that under the rephasing of the Higgs basis field H 2 → e iχ H 2 , we haveŵ → e −iχŵ , whereasv is invariant. Indeed, one easily checks thatB TB = 1 4×4 . In the 2HDM, the charged Higgs basis and the Higgs basis coincide. Thus, the matrix B rotates the Higgs basis fields into the neutral Higgs mass eigenstate fields (which includes the massless Goldstone field, G 0 = √ 2 Im H 0 1 ). More precisely, comparison with eq. (B.44) yields  Finally, we evaluate the matrix A, after using eq. (B.49). Once again, we can use the results of Table 4 .54) reproduce the corresponding 2HDM results given in ref. [27]. The power of the notation introduced in ref. [27] is clear in eq. (B.56), which depends on the only two invariant angles. In contrast, the notation of eqs. (B.7) and (B.8) conventionally used in the C2HDM community, leads to the matrix A given in eq. (B.13) which seems to depend on three angles. The true physical content of eq. (B.13) becomes apparent only after rewriting it as in eq. (B.14), which as in the case of eq. (B.56) depends only on two independent parameters.
Note that A is a real orthogonal antisymmetric matrix, as required by eq. (2.62). Indeed, the most general 4 × 4 real orthogonal antisymmetric matrix depends on two parameters, which we have identified with the two invariant angles of the neutral Higgs squared-mass matrix diagonalization.
In the 2HDM, the special form of the q j1 and q j2 allow us to rewrite eq. (B.55) as where 0 labels the first row and column of A, and the indices i, j, k = 1, 2, 3 (with an implicit sum over k) correspond to the second, third and fourth rows and columns of A. This allows one to simplify the expression for the Zh i h j vertex. This is special to the case of N = 2 and does not generalize to the NHDM with N > 2. Given the explicit forms for the matrices A and B given by eqs. (B.56) and (B.49) respectively, one can immediately check that eq. (2.67) is satisfied.

C Counting parameters that govern the A and B matrices
The matrices A and B enter the expression for the interaction Lagrangian that couples the Higgs mass eigenstates to the gauge bosons and Goldstone bosons. The matrix A is manifestly invariant under a change of the scalar basis used in expressing the NHDM Lagrangian in terms of interaction-eigenstate scalar fields. The matrix B is a pseudoinvariant quantity that depends on N − 1 unphysical phases. However, these phases can be completely absorbed into the definition of the N − 1 physical charged Higgs fields. In this appendix, we discuss the number of parameters that govern the A and B matrix in the NHDM.

C.1 Independent parameters of the matrix A
The key properties of the matrix A are given in eq. (2.62). Namely, A is an arbitrary real orthogonal antisymmetric 2N ×2N matrix. First, we note that A is a 2N ×2N nonsingular matrix such that det A = 1. Since A T A = 1 2N ×2N , it follows that det A = ±1, which implies that A is nonsingular. Moreover, for any even-dimensional 2N × 2N antisymmetric matrix A, the pfaffian of A, denoted by pf A, is defined by where ǫ is the rank-2N Levi-Civita tensor, and the sum over repeated indices is implied. A well-known result states that for any antisymmetric matrix A, 20 In particular, if A is also orthogonal then det A = 1, in which case pf A = ±1.
Next, we note that the eigenvalues of any real antisymmetric matrix A are purely imaginary. Moreover if λ is an eigenvalue of A then λ * is also an eigenvalue (see, e.g., ref. [60]). Thus, the eigenvalues of a 2N × 2N antisymmetric matrix A can be denoted by ±ia i , (i = 1, 2, . . . , n) where the a i are real and positive. We now exploit the real normal form of a nonsingular 2N × 2N real antisymmetric matrix A (see, e.g., appendix D.4 of ref. [61]). In particular, there exists a real orthogonal matrix Q such that is written in block diagonal form with 2 × 2 matrices appearing along the diagonal and the a i are real and positive. Note that the a i are the positive square roots of the eigenvalues of A T A. If in addition, A is a real orthogonal matrix, then we may use the fact that the eigenvalues of a real orthogonal matrix are complex numbers of unit modulus. In light of the above results, it follows that a i = 1 for all i = 1, 2, . . . , n. Thus, Hence, we conclude that any real orthogonal antisymmetric 2N × 2N matrix A can be parameterized by where J is defined in eq. (C.4) and Q is a real orthogonal matrix. We now employ the well-known property of the pfaffian that pf(QJQ T ) = pf J det Q. In light of pf J = 1, it follows that det Q = pf A , (C. 6) which determines the sign of det Q.
As discussed in appendix D.4 of ref. [61], the orthogonal matrix Q in eq. (C.5) is unique up to multiplication on the right by a 2N × 2N real orthogonal matrix S that satisfies SJS T = J. Such a matrix S is an element of Sp(N, where a proof of this isomorphism is given in ref. [62]. 21 Since O(2N ) is parameterized by N (2N − 1) continuous parameters and U(N ) is parameterized by N 2 parameters, we can use the freedom to multiply Q on the right by S to remove N 2 parameters from Q. This leaves N (2N − 1) − N 2 = N (N − 1) parameters in Q that cannot be removed.
That is, a real orthogonal antisymmetric 2N × 2N matrix A can be parameterized by N (N − 1) continuous parameters.

C.2 Independent parameters of the matrix B
To determine the number of independent parameters that govern the N × 2N matrix B, it is more convenient to consider the real orthogonal 2N × 2N matrixB. The transpose of this matrix rotates the charged Higgs basis fields into the neutral Higgs mass eigenstate fields S 0 β , where S 0 1 ≡ G 0 is the neutral Goldstone boson field, and It is convenient to remove the first column and the N +1 st row ofB to eliminate the neutral Goldstone bosons state. The resulting (2N − 1) × (2N − 1) matrix will be called R. We label R kβ with row and column indices that take on values from 1, 2, . . . , 2N , but excluding k = N + 1 and β = 1. That is, We can parameterize the matrix R as follows. First we define the 2 × 2 orthogonal matrix, We then define the (2N − 1) × (2N − 1) matrix R as a matrix whose matrix elements are given by Then, R can be written recursively as [63] (C.14) It is always possible to express R 2N −2 by making use of the decomposition given in eq. (C.28), where R c is a (2N − 2) × (2N − 2) matrix given in block form by 16) where C and D are (N − 1) × (N − 1) real antisymmetric matrices, and the matrix U R is a (2N − 2) × (2N − 2) real representation of the U(N − 1) subgroup of SO(2N − 2). Thus, our final expression for the matrix R is We have already noted in section 2.3 that under the rephasing of the physical charged Higgs fields, S + a → e iχa S + a , we must also transform B aβ → e iχa B aβ (in both cases for a = 2, 3, . . . , N ), so that the combinations B aβ S − a and its charged conjugate appearing in the interaction Lagrangian are invariant. Since U R in eq. (C.17) is a real representation of U(N − 1), which depends on (N − 1) 2 parameters, it is convenient to decompose this matrix into the product, is a real representation of an element of the diagonal U(1) × U(1) × · · · × U(1) subgroup of U(N − 1), which depends on N − 1 parameters, and U Rc incorporates the remaining (N − 1)(N − 2) degrees of freedom of U R . Then, we see that U Rd represents the freedom to perform separate rephasings of the N − 1 scalar doublets with zero vevs. That is, the degrees of freedom in the matrix U Rd represent the freedom to redefine the charged Higgs basis. Hence, the N − 1 parameters that govern the matrix U Rd are unphysical. The remaining parameters that describe the matrix R are physical. We can count these parameters as follows. First, the product R 12 R 13 · · · R 1,2N −1 consists of 2N − 2 angles θ 12 , θ 13 , . . . , θ 1,2N −1 . Second, the number of parameters that govern the (2N − 2) × (2N − 2) matrix R c is equal to the number of parameters needed to express the two (N −1)×(N −1) real antisymmetric matrices C and D. This provides (N − 1)(N − 2) additional parameters. Finally, we must include the (N − 1)(N − 2) parameters that govern the matrix U Rc . Thus, we have 2N − 2 + 2(N − 1)(N − 2) = 2(N − 1) 2 parameters. These are basis-invariant parameters, since they do not depend on the choice of the charged Higgs basis.
The end result is that the matrix B can be expressed in terms of 2(N − 1) 2 physical parameters. Indeed, this number can be obtained starting with the (N − 1)(2N − 1) parameters that describe the (2N − 1) × (2N − 1) real orthogonal matrix R [cf. eq. (C.10)], and then subtracting the N − 1 unphysical phases by absorbing them into the definition of the physical charged Higgs fields as described above.

C.3 The embedding of the U(N ) subgroup inside SO(2N )
We begin with the following theorem, which is useful in the analysis of spontaneous symmetry breaking of an SO(2N ) symmetric potential in a theory with a second-rank antisymmetric tensor multiplet of scalars [64][65][66].
Theorem: Suppose that Σ 0 is a 2N × 2N real antisymmetric matrix that satisfies Σ T 0 Σ 0 = Σ 0 Σ T 0 = c 2 1 2N ×2N for some real number c. Then, if the generators of the Lie algebra of SO(2N ), henceforth denoted by so(2N ), in the defining (2N -dimensional) representation are given by {T a , X b }, where the iT a and iX b are real antisymmetric 2N × 2N matrices that satisfy: then the T a span a u(N ) Lie subalgebra of so(2N ), while the remaining generators, X b , span elements of so(2N ) whose exponentials comprise the SO(2N )/U(N ) homogeneous space. Moreover, Tr(T a X b ) = 0.
Proof: First, we show that if Σ T 0 Σ 0 = Σ 0 Σ T 0 = c 2 1 2N ×2N and T a Σ 0 + Σ 0 T T a = 0, then the T a span an U(N ) Lie subalgebra. Note that these two conditions imply: In light of eq. (C.3), there exists a real orthogonal matrix W such that W M W T = diag(J 1 , J 2 , . . . , J n ) is block diagonal, where each block is a 2 × 2 matrix of the form J n ≡ 0 zn −zn 0 , where z n ∈ R and the z 2 n are the eigenvalues of M M T (or M T M ). Applying this result to Σ 0 , note that the eigenvalues of Σ 0 Σ T 0 are all degenerate and equal to c 2 . Moreover, since the matrixJ 1 2N ×2N , it follows that one can find real orthogonal matrices W 1 and W 2 such that W 1 Σ 0 W T 1 = cW 2J W T 2 = diag(cJ , cJ , . . . , cJ ), where J is the 2 × 2 matrix, That is, the factorization of Σ 0 and cJ both yield the same block diagonal matrix consisting of N identical 2 × 2 blocks consisting of cJ . Thus, there exists a real orthogonal matrix Likewise, one can use the same matrix V to define X b ≡ V X b V T . By an analogous computation, c 2 X T = Σ T 0 XΣ 0 , which implies that X T b = −J X bJ . Recall that T a and X b are both antisymmetric 2N × 2N matrices. Then, T a ≡ V T a V T and X a ≡ V X a V T are also antisymmetric. Hence, it follows that Using the explicit form forJ , eq. (C.26) implies that T a and X b take the following block form: where A, B, C and D are N ×N real matrices such that A, C and D are antisymmetric and B is symmetric. Thus, we have exhibited a similarity transformation (note that V T = V −1 ) that transforms the basis of the Lie algebra spanned by the T a to one that is spanned by the T a . Moreover, consider the isomorphism that maps i T a given in eq. (C.27) to the , we see that the A + iB are anti-hermitian generators (which are not generally traceless) that span a u(N ) subalgebra of so(2N ). We can check the number of u(N ) generators by counting the number of degrees of freedom of one real antisymmetric and one real symmetric matrix: . Taking the trace yields Tr T a X b = − Tr T a X b , and we conclude that Tr T a X b = 0.
To show that the {T a , X b } span the full so(2N ) Lie algebra, we have already noted above that there are N 2 generators, {T a }. In addition, there are N (N −1) generators, {X a }, corresponding to the number of parameters describing two real antisymmetric matrices [see eq. (C. 27)]. Thus, the total number of generators is N (2N − 1) which matches the total number of so(2N ) generators.
Any element of the SO(2N ) group is an exponential of an element of the corresponding so(2N ) Lie algebra. This provides many possible choices for parameterizing an arbitrary element of the SO(2N ) group. We shall choose T a and X b as generators of the so(2N ) Lie algebra. Exponentiating the appropriate linear combinations of generators [cf. eq. (C.27)] allows us to express any element R 2N ∈ SO(2N ) in the following form, where A, B, C and D are N × N real matrices such that A, C and D are antisymmetric and B is symmetric. Based on the discussion below eq. (C.27), we recognizeŨ R as the 2N -dimensional real representation of the group U(N ). That is, given an N × N unitary matrix U , one can identity,Ũ as a 2N × 2N real orthogonal matrix that provides the explicit form for the embedding of U(N ) inside SO(2N ).
Since the exponential of any element of the Lie algebra so(2N ) yields an element of SO(2N ), one can also choose a different order in the product of exponentials to parameterize an element of SO(2N ). For example, one can also express any element R 2N ∈ SO(2N ) in the following form, where C ′ and D ′ are real N × N antisymmetric matrices and W is an N × N unitary matrix. In general, C ′ = C, D ′ = D and W = U .

D Cubic couplings of the Goldsone boson and physical Higgs scalars
It is quite remarkable that the cubic scalar couplings that involve a single or two Goldstone fields can be simplified to expressions that involve either squared mass differences or squared masses of the corresponding physical Higgs scalars, as exhibited in eqs. (2.90)-(2.93). This was first noted in the context of the CP-conserving 2HDM in Ref. [20,40]. Achieving such simplified forms for these couplings is rather laborious. It is instructive to provide some of the details of the derivation. As an example, we demonstrate below how to derive the coupling between a neutral scalar (here denoted by S 0 p ) and two Goldstone bosons G 0 = S 0 1 obtained in eq. (2.93). Starting from eq. (2.85), there seem to be three relevant terms: those with (β, γ, δ) = (p, 1, 1), (β, γ, δ) = (1, p, 1), and (β, γ, δ) = (1, 1, p). However, the latter two vanish due to eqs. (2.30)-(2.31). Indeed, when β = 1, the result involves Applying eqs. (2.30)-(2.31) to the remaining term, we find for 4v thus reproducing eq. (2.93). The most crucial step is the first, where eq. (2.46) was used to relate these couplings with the mass matrices. The third line above is obtained by using eqs. (2.32) and (2.44), where the (11) entry vanishes, which shows that the charged boson masses do not contribute, as expected. The fourth line is obtained by breaking (V † ) θjvj in its real and imaginary parts, and then using eqs. (2.33) and (2.45) to show that the imaginary part involves the vanishing entries of D 2 0 . Eq. (2.34) yields the fifth line. To proceed, we break the remaining V † V matrix into its real and imaginary parts. According to eq. (2.66) the real part is the unit matrix, while the imaginary part is antisymmetric. The two terms involving the symmetric matrix Im(V † V )D 2 0 Im(V † V ) cancel each other. Given eqs. (2.45) and (2.66), we reach the last line.
It is noteworthy that it took such a long calculation to obtain such a simple result. In fact, it turns out that such proofs are simpler when performed in the charged Higgs basis.

E Generalized sum rules
In this appendix, we rederive in detail the sum rules obtained by Gunion, Haber and Wudka in ref. [19]. Following the conventions of ref. [19], we indicate vector bosons with indices a, b, c, . . . and scalars with indices i, j, k, . . .. The Feynman rules for the cubic vertices are A α a A β b A γ c : i g abc (p a − p b ) γ + (p b − p c ) α + (p c − p a ) β ≡ i g abc Γ αβγ (p a , p b , p c ), (E.1) A α a A β c φ i : i g abi g αβ , (E.2) A α a φ i φ j : i g aij (p i − p j ) α , (E. 3) with all momenta incoming, and the Feynman rules for the quarticic vertices are A α a A β b φ i φ j : i g abij g αβ , (E.4) For the sum rule involving four gauge bosons we also need the relation between the quartic and cubic term. We adopt here the conventions of Cornwall, Levin and Tiktopoulos [17]. We just need the relevant terms, This gives the following Feynman rules, A µ a A ν b A ρ c A σ d : −8i (D abcd g µν g ρσ + D acbd g µρ g νσ + D adbc g µσ g νρ ) , (E.6) where we note, for future reference, that comparing eq. (E.1) and eq. (E.7) we get, Cornwall, Levin and Tiktopoulos [17], and independently Llewellyn Smith [16]) show that, in order for unitarity to hold, the couplings C abc and D abcd must be those of a gauge theory. In particular, D abcd = 1 8 (C ace C bde − C ade C cbe ), (E.9) 0 = C abe C cde − C ace C bde − C ade C cbe , (E. 10) where the last relation is the Jacobi identity. This means that C abc are the structure constants of the gauge group. As we want to write everything in terms of the structure constants g abc of ref. [19], we use eq. (E.8) to obtain D abcd = − 1 8 (g ace g bde − g ade g cbe ), (E.11) 0 = g abe g cde − g ace g bde − g ade g cbe . The diagrams contributing to the scattering A a (p 1 )+A b (p 2 ) → A c (p 3 )+A d (p 4 ) are given in figure 1. In an obvious notation we will name the amplitudes according to the Mandelstam variables channel (s, t or u) and by the particle being exchanged. We get M 4Point = (−i) 2 8 (D abcd g αβ g γδ + D acbd g αγ g βδ + D adbc g αδ g βγ ) f αβǫδ , (E.13) To determine the coefficients of the high energy behavior (see eq. (E.31) below) we cannot use the approximate expression on the right-hand side of eq. (E.21) because we would then lose contributions from the E 4 terms that modify the E 2 terms. So we should use consistently the definitions of the left-hand side and expand the result in powers of s, t or u. For instance, for particle a we have, ǫ L a = (γ a β a , γ a β a /β a ), β a = to express the result in terms of just two independent variables. Which ones should be used will depend on the problem. Next we want to isolate the terms that grow with E 4 and E 2 . To achieve this we make the scaling s → s/x, t → t/x, u → u/x, (E. 30) and perform an expansion for small x. The terms in E 4 are the coefficients of x −2 and the terms that grow like E 2 are the coefficients of x −1 . Therefore, we can write for each amplitude where we assumed that the independent Mandelstam variables are s and t. We did this consistent expansion using FeynCalc and Mathematica for the Lorentz algebra and series expansion, respectively. To have an idea of what is involved we just write the exact amplitude for the s-channel exchange of a gauge boson: where the antisymmetry of the constants g abc was used. So the more divergent terms cancel only with the gauge part. The constraints that emerge simply imply that we must have a spontaneously broken gauge theory [16,17] as given in eqs. (E.11)-(E.12).

The E 2 terms
Having shown that a spontaneously broken gauge theory assures that the the most divergent high energy behavior cancels, we consider next the terms that diverge like E 2 . Here the gauge theory part is not enough to achieve cancellation and we get constraints on the gauge boson couplings to scalars. For convenience we definê We also note that once we use eq. (E.33) there is no contribution at this order for B st i . The results are summarized in table 6.
E.1.3 The sum rule of ref. [19] To obtain the sum rule in eq. (2.4) of ref. [19] we take as independent the Mandelstam variables s and t. The coefficients of the terms growing with s and t must vanish. If we  Table 6. where the prime in ′ indicates that the sum only runs over massive gauge bosons.

E.1.4 Another sum rule
If we take the coefficient of t, that is the sum of theB t i we obtain another sum rule: To obtain the sum rule in eq. (2.5) of ref. [19], we take as independent Mandelstam variables s and u. The coefficients of the terms growing with s and u must vanish. If we take the coefficient of s we obtain the desired sum rule. (g cik g abk − g bik g ack ) . (E.53)

E.2.4 Another sum rule
If we take the coefficient of u we get another sum rule, The diagrams contributing to the scattering A a + A b → φ i + φ j are given in figure 3.  Figure 3. Amplitudes for the scattering A a A b → φ i φ j The diagram with the triple Higgs vertex does not exhibit any bad high energy behavior, and therefore one does not have to consider it here. The other amplitudes are M 4Point = (−i)(i g abij g αβ ǫ α (p 1 )ǫ β (p 2 ), (E.55)