Higgs–Higgs scattering and the (non-)existence of the Higgsonium

We study the Higgs–Higgs scattering process and the possible emergence of a Higgs–Higgs bound state (Higgsonium) in any Higgs potential with the vacuum expectation value and second derivative matching the corresponding values from the Standard Model (SM). From the tree-level Higgs–Higgs scattering amplitude, we construct the unitarized amplitude using two different unitarization schemes (the well-known on-shell and N/D methods). We reproduce the known result that there is no Higgsonium state in the SM and, in addition, we determine the S-, D-, and G-wave SM scattering lengths, both at tree-level and upon unitarization. In doing so, we refine previous results by checking the convergence of the N/D approach. Next, we extend the calculation for non-SM potentials and investigate under which conditions a formation of a bound state close to the Higgs–Higgs threshold is possible. In this way, the assumption that no Higgsonium exist, imposes certain bounds on the values of the self-interaction parameters that complement those imposed by the vacuum stability condition.


Introduction
The electroweak sector is considered apt to search for the signatures of departures from the standard model (SM).Aiding this belief are the various (but inconclusive) inconsistencies observed in the experiments.Examples include, the R D ( * ) anomalies [1,2,3,4,5], interactions of the W ± and Z 0 bosons [6,7], the muon g − 2 problem [8,9,10], etc.Most recently, the CDF collaboration reported a tension in the mass of the W -boson compared to the value predicted by the SM [11].Another facet of the electroweak physics is the possibility of discovering new forces of nature.The spontaneous breaking of the EW symmetry results in the W ± and Z 0 bosons acquiring a large mass.This mechanism can also be implemented in scenarios where there are more than one Higgs-like state which acquire large masses upon breaking of the symmetry.Many such scenarios have been explored theoretically leading to various models [12].The largeness of the masses of these additional Higgs-like states makes them -at present-not observable in experiments.Indirect searches for such states, however, have lead to null results [13].There is no doubt that this sector of the SM has elicited a lion's share of the interest from the theoretical as well as the experimental community.
A long standing topic of interest is the stability of the EW vacuum, which is related to the shape of the Higgs potential [14,15,16,17,18,19].According to the SM, the Higgs potential is a "Mexican hat" potential with the La-grangian containing a mass term, cubic term (tri-Higgs coupling, in short 3HC), and a quartic term (4HC).However, any Higgs potential that reproduces the vacuum expectation value (vev, the value of the field at the minimum of the potential, denoted as v) and the Higgs mass (fixed by the second derivative at the minimum and called m H ) is in principle acceptable, since three-or four-leg (or larger) self-interactions could not be measured yet.In this respect, the SM Higgs is the "most economic" choice, since there is no term beyond the four-leg one.
An important consequence of the uncertainty in the shape of the Higgs potential is the possibility of the vacuum decay.Simply put, if the electroweak vacuum is the global minimum of the Higgs potential, then the Universe is in a stable vacuum.Otherwise, if there are minima deeper than the ones predicted by the SM, then our Universe is not in the global minimum of the Higgs potential, and in principle, can transition into the global minimum.As such, this has been a hotly debated topic [20,21,22,23,24,25,26,27,28,29,30,31,32].The precise shape of the potential will depend entirely on the Higgs self-coupling constants, and understanding them is of foremost concern towards the confirmation of the SM or for finding new processes.
Many different versions of the Higgs potential have been proposed in the literature (see, e.g, Ref. [12] for a review 1 ).These potentials include the likes of Coleman-Weinberg Higgs [35], Nambu-Goldstone Higgs [36,37], tadpole-induced Higgs [38,39], etc which behave like the SM Higgs potential only in the vicinity of the SM minima, and the extensions of the SM -like the additional singlet [40,41,42,43,44], two Higgs doublets (2HDM) [45,46,47], minimally supersymmetric SM Higgs (MSSM) [48,49], etc.Another common feature of all these potentials is that the vev of the field (in case of doublet models allowing more than one vev's, then the sum of the squares of the vev's, see e.g., Ref. [12] for exceptions to this rule) agrees with that of the SM.Furthermore, in models with more than one Higgs-like field, the physical Higgs either emerges as an admixture of the fields [50,51] or is decoupled from the rest (see e.g, Ref. [52,53]).In the former case, the mixing angle(s) can be constrained using the mass of the physical Higgs.However, the Higgs self-coupling are (largely) unconstrained.It is pertinent to note that searches for additional Higgs have so far provided null results [54,55,56].
The experimental data on the Higgs self-couplings are scarce.This is partly due to the fact that the current state-of-art is insufficient to perform Higgs-Higgs scattering experiment, and the accessible probes like di-Higgs and tri-Higgs production through gluon-gluon fusion, the top quark-Higgs (t thh) coupling, and so on have large backgrounds as well as non-perturbative effects arising from strong interactions that limit the accuracy of the measurements [57,58,59].The possibility of measurements of the self-couplings at various future experiments has been studied in various works [58,59,60,61,62,63,64,65,66,67,68] and is a topic of great interest for the community [69,70,71].Two alternate avenues -pair production of the Higgs [72,73,74,75] and the single Higgs production process [76,77,78] -were explored to study the self-coupling of the Higgs.The tri-Higgs coupling constant was extracted from the resonant and non-resonant pair production of the Higgs boson in the b b channels, as well as the pp → HH channel.The resonant production of Higgs studied by AT-LAS provided the bounds of −5 < d 3 < 12 [72] which was later revised to −0.4 < d 3 < 6.3 using additional data and by including the single Higgs production channels [79].The non-resonant production studied by CMS provided the estimate −1.24 < d 3 < 6.49 [80].Hopefully, the knowledge about this part of the SM can be improved and resolved in more detail in future colliders [81,82,83,84,85,86].
Moreover, if the extended Higgs potential has a suitable form, with the 3HC and 4HC coupling constants taking suitable values, one can expect to have a di-Higgs 1 An important related issue of the Higgs theory is the socalled triviality, which implies that the Higgs self-interaction coupling constant(s) as well as the mass go to zero under the renormalization group (RG) evolution if the RG cut-off is taken to be infinite [33].Due to this, the SM Higgs boson is expected to exhibit self-duality [34].The triviality problem requires that the Higgs field be associated with a finite scale.This scale (and the possible new degrees of freedom) is in itself a topic of intense debate [12].
bound state or a Higgsonium (also called the Higgsium).A study of the relationship between the coupling constants and the Higgsonium can shed light on the shape of the potential around the minima.One of the earliest studies of the formation of the bound state of two Higgs was performed by Cahn and Suzuki using a model similar to the linear sigma model [87].This model is essentially the same as the non-interacting Higgs doublet model.They find that a (shallow) two-Higgs bound state could exist only if the Higgs has a mass > 1 TeV (which is presently excluded).Similar results were obtained with interacting Higgs using the N/D unitarization scheme [88] and with the SM Higgs using Bethe-Salpeter equations [89].An analogous outcome was found through a relativistic calculation using the minimal SM Lagrangian [90].An extension of that work to the SM Lagrangian showed that a bound state can be observed for a light Higgs (in the vicinity of the physical mass of the Higgs boson), only if the Higgs self-coupling constants are large [91].A nonrelativistic study of interactions between two Higgs in a two-Higgs doublet model also showed the presence of Higgsonium even when the mass of Higgs is smaller than the presently known value [92,93].
A systematic study based on the linear realization of the Higgs effective field theory (HEFT) showed that a Higgsonium is possible in the strong regime of Higgs selfcoupling, reaffirming the earlier results [94].The nonrelativistic HEFT derived in Ref. [94] was extended to the two-Higgs doublet model (2HDM) in Ref. [95], that showed that the heavy partners of the Higgs -called the "obese Higgs" -can form bound states in the Type-I 2HDMs and not in the type-II 2HDMs, provided the mass of the "obese Higgs" is less than the scale of expansion of the HEFT (and greater than the mass of the SM Higgs).That work also explored the decays of such a bound state into various channels.
This paper is organized as follows.In Sec. 2 we detail the formalism used in the study.We present the results of the calculation in Sec. 4 and conclude in Sec. 5.
In the present work, we are concerned with a related possible avenue for the studying the Higgs self-coupling.Namely, we are interested in the study of the Higgs-Higgs scattering and the consequent determination of the (S-, D-, and G-wave) scattering lengths and phase shifts in both the SM Higgs potential and in extensions of the latter.An important part of this work concerns the unitarization of the scattering results of the SM Higgs Lagrangian (and extensions thereof).To this end, we use the well-known onshell and the already mentioned N/D schemes to improve the tree-level results for Higgs-Higgs scattering (scattering length and phase shift).The use of two schemes allows for a direct comparison, thus checking the independence of the result on the adopted unitarization.In particular, for the N/D case, we also go beyond the usually employed lowest order approximation and show (to our knowledge for the first time) the convergence of the results.Moreover, we also revisit the possibility of the formation of the bound state of two Higgs bosons by investigating in a systematic way the whole parameter space for the three-and four-leg interaction strengths, upon keeping the vev (v) and the mass (m H ) fixed at their physical values.We find that there is a minimum value of the 3HC coupling constant below which no bound state can exist and a lower bound for the 4HC coupling constant above which the repulsion does not allow for a bound state.In this way, assuming that no Higgsonium exists, a further constrain of the 3H-4H parameter space that improves previously known bounds is determined.
Quite interestingly, the study of scattering and the possible emergence of bound states has been recently also tackled in other areas of physics: in Ref. [96] the bound states of two dilatons/glueballs has been studied, finding that it is quite possible to have such a state in YM theory of QCD.This is particularly interesting since the dilaton potential resembles the Higgs one (indeed, it is similar to the Coleman-Weinberg Higgs mentioned above [35]).Bound states involving glueballs have been later on investigated in Ref. [97].An interesting novel development in the study of bound states has been recently applied to the scattering of gravitons [98,99], leading to the graviball hypothesis in Ref. [99].

The Higgs potential 2.1 General considerations
Consider a potential V (H) which has the basic features of the Mexican hat potential in the relevant part of the domain, but can vary arbitrarily at large values of the field or when the field is close to zero.The potential must, 1. possess minima at H = ±v, where v ∼ 246 GeV is the vev [13], and 2. have a value equal to m H for the second derivative at the minima (H = ±v).If the potential is well behaved around H = ±v, we can expand it in a Taylor series as, where (dV /dH) | H=v = 0 because of the extremum.In the last line, h = H − v has been defined.Note that the coupling constants g nH represent n−Higgs interactions, with g 2H = m 2 H , g 3H = g, and g 4H = λ.The assumptions listed above restrict the mass of the Higgs boson and the vev of the Higgs potential.However, the coupling constants g and λ remain unrestricted.Thus, the current experimental data can fix only the position of the minimum of the Higgs field and its second derivative around the minimum.The shape of the potential in general is not probed.We use this observation to express the Higgs self-coupling constants as in Ref. [19], The numerical values d 3 = d 4 = 1 correspond to the "standard model" self-coupling and any deviation from these values would hint a disagreement with it.In the discussions below, we alternate between (g, λ) and (d 3 , d 4 ) based on convenience.(This is also true for the notations h and H representing the Higgs field.)

Higgs-Higgs scattering
At tree-level, only the cubic and quartic terms contribute to the HH → HH scattering, as shown in the Fig. 1.The amplitudes for these diagrams are given by, where s, t, u are the Mandelstam variables.Upon projecting the amplitudes to the s-channel and later onto various angular momentum channels, we get the partial wave amplitudes as, where ℓ is the orbital angular momentum quantum number, Q ℓ (s) are the Legendre functions, and k = Fig. 2. Assuming the non-existence of the Higgsonium, the allowed parameter space in the (d3, d4)-plane lies above the dotteddashed (on-shell) and the dashed (N/D) curves.The requirement that the vacuum of the Higgs potential corresponds to an absolute minimum implies the region above the solid line [19].As it is visible, the non-existence of the Higgsonium further constraints the parameter space in the right-upper part of the (d3, d4)-plane.Inset: The region where the different curves cross each other (see text for details).s − 4m 2 H is the 3-momentum of any of the final states in the center of mass frame of reference.The odd-ℓ partial waves vanish due to bosonic nature of the Higgs.One crucial observation is that the ℓ = 0 partial wave contains a branch point at s = 3m 2 H in addition to the pole at s = m 2 H .This feature is used in the so-called "on-shell" unitarization scheme followed in the present work.
Finally, the s-wave tree-level scattering length (a ℓ ) reads 2 For the SM parameters, we get: which is in good agreement with the result from nonrelativistic effective theory calculations [100].Such a small value for the s-wave scattering length is a sign of the weakness of the interaction between the Higgs states.For comparison, the s-wave scattering length for W ± W ∓ , and W ± Z 0 scattering is nearly 3 orders of magnitude higher 3 2 We follow the convention where δ0 = k a0 [101].Thus, our scattering length varies by a sign compared to that of [100]. 3It should be noted that the scattering between the W ± and/or Z 0 states are inelastic and hence the scattering lengths have complex values.Thus, in the comparison above one refers to their absolute values.[102,103,104].The smallness of the SM HH-scattering length shall be confirmed later on by the very small unitarity corrections.For the higher partial waves, we can evaluate the scattering lengths (or better, hyper-volumes).The general expressions for the scattering length is, Using Eq. 9, we get the values for ℓ = 2 and ℓ = 4 partial waves as, Interestingly, the scattering amplitudes in higher partial wave channels involve only the attractive three-leg interaction.Later on, we shall reevaluate the scattering within the introduced unitarization schemes.Remarkably, for the SM case, the modification turn out to be completely negligible.
In the future, it would be also interesting to calculate how additional SM processes (as for instance the box diagrams with top quarks) would modify, within the SM, the effective parameters d 3 and d 4 .Due to the suppression of the involved diagrams, one should expect only small deviations from d 3 = d 4 = 1.Note that the effects of the quantum corrections on the Higgs self-coupling constants have been studied up to the next-to-next-to-leading order in various extensions of the SM.
The estimates vary depending on the type of extensions considered as well as the parameters of the respective models [105,106,107,108] but is of (maximally) order 10 2 for d 3 .Thus, the qualitative conclusions arrived at in the present work remain valid and Eq. 10 could be still used to determine the scattering length for the updated 'dressed' values of d 3 and d 4 .

Unitarization
In this section we discuss the main features of a generic unitarization scheme and present two concrete ones: the on-shell and the N/D schemes.Then, we shall apply them to the Higgs-Higgs scattering process.

General considerations
The amplitude calculated in the previous section was obtained at tree-level, which, as it is known, fulfills unitarity only order-by-order in perturbation theory.An effective way to respect unitarity (as well as including the effect of quantum fluctuations) is possible via an appropriate unitarization scheme that provides a prescription to replace the amplitude where A U l (s) is the unitarized amplitude.A unitarized amplitude typically consists of (some) contributions coming from all orders of perturbation involving the vertices and degrees of freedom present in the tree-level amplitude.In the literature, different unitarization schemes have been investigated (see e.g.Refs.[96,109,110,111,112,113,114,115,116,117,118,119,120] and refs therein).These schemes have been used to study a wide variety of topics including, but not limited to ππ, πK [109], ρρ scattering [110], and lately to the scattering of two scalar glueballs and the possibilities of a "glueballonium" [96],etc.In fact, the unitarization is a common tool in the theoretical search for short-lived (a.k.a, broad) resonances like the scalar mesons [120].
Before discussing the two chosen approaches (on-shell and N/D) in more in detail, we provide some relevant general properties of a scattering process.For simplicity, we concentrate on the S-wave only (ℓ = 0), whose index is omitted.The phase-shift δ(s) is linked to the unitarized amplitude A U (s) by the following relation: where the phase-space function ρ(s) reads Since for elastic scattering δ(s) is purely real, Eq. 16 implies that: This is a necessary condition for any unitarization approach.
The emergence of the bound state is also connected to the value of the unitarized scattering length, which, for the S-wave, reads When a bound state forms right at threshold, the scattering length diverges.Moreover, the scattering length is positive when no bound state forms and negative otherwise.

On-shell unitarization
In the on-shell unitarization approach [109,111,121,122], one introduces the loop function Σ(s) entering into the resummation where Im (Σ(s)) = ρ(s) for s > 4m 2 H is required to fulfill Eq. 18.The main feature of this approach is that the treelevel amplitude and the loop factorize.
The loop function is fixed by considering two subtractions [96]: whose analytical form is: Note, while a single subtraction is enough for convergence, two subtractions are requested for avoiding the emergence of unphysical ghost states [123,124].Here, the first subtraction is taken at s = m 2 H in such a way to maintain the Higgs mass at the same physical value even after performing the unitarization.The second subtraction is chosen at s = 3m 2 H , in order to preserve the divergence of the tree-level amplitude at the branch point of the left-hand cut.In fact, this branch point is also caused by the single particle pole at t = m 2 H and u = m 2 H projected onto the s-channel.In other words, for s ∼ 3m 2 H , the unitarized S-wave amplitude is given by  Finally, a bound state with mass M 2H is a solution of the equation that implies that the S-wave scattering length diverges when M 2H = 2m H .This is due to the fact that a bound state just at threshold generates a divergent pole of the amplitude.A feature of the on-shell method is that M 2H ∈ ( √ 3m H , 2m H ), but the lower limit √ 3m H should not be regarded as a physical constraint.In turn, it implies that this unitarization scheme should not be trusted when M 2H is far below the threshold.
The unitarized scattering length in the on-shell approach reads where is the total contribution of the loops to the scattering length.For the SM parameters (d 3 = d 4 = 1), this value is extremely small compared to the inverse of the tree-level amplitude at threshold, A(s = 4m 2 H ) and hence can be neglected.Thus, the value of the unitarized scattering length is nearly the same as that obtained at the tree-level.This conclusion applies with even better accuracy to the higher D-and G-waves, as explicitly shown in Ref. [96].
The irrelevance of the unitarization (and thus of quantum fluctuations) does not, however, apply in general: the closer one gets to the generation of a bound state, the more important are the loops.

The N/D unitarization
In the N/D approach [125,126,110,112,87,113,96], the unitarized S-wave amplitude takes also into account the left-hand cut behavior.Namely, at tree-level, one has This is due to the single-particle pole at t = u = m 2 H . Here, Eq. 26 is fulfilled by requiring that the unitarized amplitude is expressed as a ratio where the numerator N contains the left-hand cut and the denominator D the right-hand cut: and By construction, Eq. 18 and Eq. 26 are used.We require also that D(m 2 H ) = 1, thus the residue of the amplitude at the single-particle pole is unchanged by quantum fluctuations.
The previous constraints lead to the following coupled integral equations for N and D: where we have taken into account that the numerator reduces to the tree-level amplitude at lowest order (see below).The mass of the Higgsonium bound state is obtained by A useful method to solve the system above is by iterations.By denoting N (k) (s) and D (k) (s) as the functions N and D at the k th iteration, at lowest order we set: thus the amplitude (tree-level result).At the next order, the denominator takes the form (34) which leads to the first non-trivial N/D result [87,110] with where D (1) (s) fulfills Im D (1) (s) = −ρ(s)N (0) (s), while Im N (0) = σ(s).Note, the single subtraction in the denominator is taken at s = m 2 H giving D (1) (m 2 H ) = 1, as desired.Here, there is no need for further subtractions, since no ghost appears in this approach.Moreover, the N/D amplitude automatically keeps the branch point ) (even though the constant in front of it is modified by the unitarization).
At this order, the bound state equation is D (1) (s = M 2 2H ) = 0. Quite remarkably, it turns out that the results are quite accurate, since higher orders cause only small changes, as we show in Appendix A. Due to its relatively simple implementation, we shall use this approximation in the following plots.
Going further, one needs to evaluate the chain: leading to the function N (1) (s) and D (2) (s) given by: At an arbitrary iteration: with Im Finally, the iterative N/D amplitudes read: (41) Clearly, the bound states are searched as poles of A N/D(2k+1) (s), alias the zeroes of D (k+1) need to be determined.Within the present method, the mass of the bound state belongs to the interval (m H , 2m H ), thus the lowest limit approaches the mass of the Higgs.In this respect, this method allows to study situations in which the mass of the bound state is also well below the threshold, thus represents an improvement w.r.t. the on-shell scheme.Nevertheless, the results should be comparable for what concerns the emergence of the bound state.This turns out to be the case, as we will show in the next section.

Results and Discussion
In this section, we present the results of the present study.To this end, the tree-level amplitude obtained in Eq. 8 was unitarized using the two schemes discussed in Sec. 3. The unitarized amplitudes provide bounds for the values of the coupling constants that do not allow for a bound state, as summarized in Fig. 2.
Namely, let us first start with the on-shell scheme by assuming (in agreement with present experimental data) that no Higgsonium exists, in the on-shell scheme the corresponding parameter space in the (d 3 , d 4 )-plane lies below the dot-dashed blue line shown in Fig. 2.This curve corresponds to the equation which can be (also) obtained from the divergence of the scattering length.Furthermore, the following condition holds true.For a fixed d 3 , Next, the N/D unitarization scheme delivers qualitatively similar results as can be also seen from Fig. 2. The dashed red curve in Fig. 2 represents the upper bound for the coupling constants obtained from the 1 st iteration of the N/D scheme (i.e., D (1) (s)).Successive iterations do not alter this curve by any large amount, as discussed in detail in Appendix A.
For completeness, we also discuss the bounds to the coupling constants d 3 and d 4 given by the demand that h = 0 corresponds to the global minimum of the potential.(These results have been presented in, e.g., Ref. [19]; we briefly recall them here for the sake of clarity.)Referring to Eq. 2, the extrema of the potential are given by the zeroes of the first derivative of the potential with respect to the field.Considering only the first three terms, Thus, the extrema of the potential are, The first root (h = 0) corresponds to the EW vacuum (more precisely, the point H = v), and is a local minimum (as can be verified from the second derivative of the potential).The other two roots form a pair of minimum and maximum of the potential.In the SM (d 3 = 1 = d 4 ) these points reduce to h = −2v and h = −v respectively (equivalently, H = −v and H = 0).These two roots are real only if d 3 ≥ 8d4 9 .Thus, if we demand that the root h = 0 represents a global minimum, we get the lower bound for the parameters [19].The corresponding curve is represented by the solid black line in Fig. 2. It should be also stressed that the above derivation does not require that the potential is Z 2 -symmetric around h = −v.Note that, the curves, including the on-shell and N/D, are symmetric under d 3 → −d 3 .This is due to the fact that the tree-level amplitude is quadratic in d 3 .But, this is to be expected as d 3 characterises the attractive part of the Higgs potential and a change in the sign will flip the positions of the maxima and the minima given in Eq. 45.
The bounds on the value of d 3 obtained here lie close to the ones measured by the experiments [72,79,80].Particularly, the value obtained from the double-Higgs pair production processes is consistent with the lower limit obtained in the present work [72].However, additional data and the inclusion of the single Higgs production measurements tightens the experimental bounds [79], which has also been corroborated by the resonant production channel measurements [80].In this case, we conclude from our work that the Higgsonium cannot exist, in the SM as well as beyond it.However, it would serve useful to continue the search for Higgs-Higgs self interactions, as it can provide rich information about the SM and eventual extensions of it.Summarizing, it is visible from Fig. 2 that the nonexistence of the Higgsonium further constrains the admissible region of the (d 3 , d 4 )-plane and constitutes a more stringent constraint w.r.t. the global minimum one of Ref. [19] in the upper part of the parameter plane.
Next, we present other quantities relevant quantities for the the Higgs-Higgs scattering.Using the on-shell scheme, we plot the scattering phase shift in Fig. 3.The SM case is basically not distinguishable from the horizontal axis, implying that the phase-shift is always very small (in agreement with the fact that the tree-level result offers a very good approximation).In general, the phase-shift is positive when no Higgsonium is present, but is negative otherwise.Note, in the latter case the S-wave scattering phase shift goes to −2π as s increases, thereby confirming the presence of a bound state [127] (in fact, two states are present below threshold: a bound state and the Higgs particle, each of them generating a −π shift).
In Fig. 4 we plot, for the N/D case, the scattering length as a function of the 4HC coupling d 4 for two illustrative values of d 3 .For d 3 = d 4 = 1 (SM value, left plot), the scattering length is positive; then, for increasing d 4 it decreases and eventually becomes negative.The smooth behavior of the scattering length indicates the absence of any bound state4 .For the illustrative large value of d 3 = 15 (right plot), the scattering length is negative for any value of d 4 420 signalling the presence of the Higgsonium.The divergence at d 4 ∼ 420 implies that the bound state is realized at threshold.For larger values of d 4 there is no bound state.
Finally, we show in Fig. 5 the mass of the Higgsonium as derived from the pole of the unitarized amplitudes for a fixed value of d 4 = 272, as a function of d 3 .As the strength of attractive interactions (d 3 ) is increased keeping the repulsion (d 4 ) fixed, the depth of the bound state increases.The converse is true for d 4 , the depth of the bound state for a given d 3 decreases as repulsion is increased.In the end, it should also be noted that, by construction, the mass of the Higgsonium goes to 3m 2 H (m 2 H ) in the onshell (N/D) unitarization scheme as d 3 → ∞, implying that an agreement is only possible when the bound state is not too deep.Yet, the results are quite similar for all the range d 3 20 under investigation.

Summary and Conclusions
In this work, we have studied the scattering of two Higgs bosons.The tree-level scattering amplitudes are derived from the generalised Higgs potential and the scattering lengths for the S-, D-, and G-waves have been evaluated (the particular case of the SM parameter has been also outlined).
Then, the scattering amplitudes are unitarized using the well-known on-shell and the N/D unitarization schemes.The SM scattering lengths are basically left unchanged, thus offering a clear prediction for the SM outcomes.Away from the SM values, the unitarization renders the investigation of an eventual bound state possible.Both procedures indicate the presence of a bound state of two Higgs bosons (Higgsonium) provided the the 3H self-interactions of the Higgs boson is at least 10 times larger than the value predicted by the SM.Indeed, the requirement that the Higgsonium does not exist generates some constraints for the 3H and 4H coupling that are more stringent than the ones obtained by the global minimum requirement in the upper part of the 3H-4H coupling plane, and in line with the current experimental bounds for d 3 within the uncertainty levels.
In the future, experimental searches for Higgs-Higgs scattering (and, in a broader sense, for a Higgsonium) could turn up some interesting results and help us further constrain the SM results.It is also possible (although at present quite speculative) that a two-Higgs resonance state exists well above threshold [128] and the constraints one can derive from the study of such a state could be used in conjunction with those from the present work to further narrow down the 3HC and 4HC parameters.

A N/D at higher orders
In this appendix, we explore the effects of going to higher orders in the N/D unitarization scheme.The variation of the mass of the bound state as a function of the order (k) of the denominator (D (k) (s)) is shown in Fig. 6.A slight variation is seen in the mass when going from k = 1 to k = 2 and k = 3 in the numerical evaluation of the N/Dunitarized amplitude.However, the value stabilizes soon after and for k ≥ 4 the mass of the bound state is fixed.
Plugging Eq. 39 and 40 into each other iteratively, one finds that the denominator at order k can be written as, D (k) (s) = D (1) (s) + I(s) − C (k−1) (s); k ≥ 1 , (46) where,  .( 48) Above, D (1) (s) is as given in Eq. 34.In Eq. 46, the first two terms are independent of the order of approximation of the denominator (note, I(s) = C (0) (s)).In the integrand of I(s), 0 < ρ(s 1 ) ≤ 1/32π, |σ(s H owing to the domains of integrations.Thus, we can impose an upper bound for the integral given in Eq.48 as, Using the same arguments, we can show from Eq. 34 that, |D (1) (s)| < 1.These bounds can then be used to show that |D (k) − 1| ≪ 1 as long as the tree level amplitude is finite in the region around the expected pole.It then follows that the contribution of C (k−1) (s) is small in the vicinity of s = s pole , provided s pole ∈ (3m 2 H , 4m 2 H − δ).Thus, close to the pole of the possible bound state, as the order of the denominator increases the corrections keep reducing in magnitude.This behavior carries over to the position of the zero of the denominator D(s), as demonstrated in Fig. 6.This also means that the upper bound imposed on the value of the 4HC parameter d 4 does not vary significantly from that derived from the 1 st iteration (shown in Fig. 2).

Fig. 1 .
Fig. 1.Tree level diagrams representing the HH → HH scattering process.The red square and black circles represent the vertices coming from the H 4 and H 3 terms respectively.

Fig. 3 .
Fig. 3. Plots of the S-wave scattering phase shift for various values of the parameters.

Fig. 4 .
Fig. 4. The Higgs-Higgs scattering length as function of d4 for two different values of d3.The starting point on the left graph gives the value for the SM parameters.Numerically, it is very close to the tree-level value of Eq. 11 .

Fig. 5 .
Fig. 5. Plot of the mass of the Higgsonium as a function of d3 for a representative value of d4(= 272).The black solid line represents the value from the N/D unitarization and the blue dashed line gives the value from the on-shell unitarization.The first two vertical lines represent (from left to right) the critical values for having an Higgsonioum in the on-shell and N/D schemes, and the third vertical line marks the region excluded by the global-minimum requirement.

Fig. 6 .
Fig. 6.Plot showing the mass of the Higgsonium at different orders of the N/D unitarization.The vertical axis is in the units of mH , see text in Appendix A for details.