aainstitutetext: Institut de Physique des Deux Infinis de Lyon (IP2I),
CNRS/IN2P3, 4 rue Enrico Fermi, 69622 Villeurbanne Cedex,
Francebbinstitutetext: Université Claude-Bernard Lyon 1, Univ Lyon, Lyon,
Franceccinstitutetext: Center for AI and Natural Sciences, KIAS, Seoul 02455,
Koreaddinstitutetext: Institut für Theoretische Physik und Astrophysik, Uni
Würzburg, Emil-Hilb-Weg 22, D-97074 Würzburg, Germany

# Exploring extended Higgs sectors via pair production at the LHC

G. Cacciapaglia c T. Flacke d M. Kunkel d W. Porod d L. Schwarze
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>wuerzburg.de<EMAIL_ADDRESS><EMAIL_ADDRESS>wuerzburg.de

###### Abstract

Higgs sectors extended by electroweakly charged scalars can be explored by
scalar pair production at the LHC. We consider a fermiophobic scenario, with
decays into a pair of gauge bosons, and a fermiophilic one, with decays into
top and bottom quarks. After establishing the current bounds on simplified
models, we focus on an SU(5)/SO(5) composite Higgs model. This first
exploration demonstrates the need for dedicated searches at current and future
colliders.

††preprint: KIAS-A22008

##  1 Introduction

The Standard Model (SM) of particle physics contains a single scalar field, a
doublet of weak isospin $\text{SU}(2)_{L}$ that is responsible for the
breaking of the electroweak (EW) symmetry Englert:1964et ; Higgs:1964pj . Upon
acquiring a vacuum expectation value (VEV), a massive physical scalar particle
arises, the famous Higgs boson Higgs:1964pj  discovered in 2012 at the Large
Hadron Collider (LHC) experiments ATLAS:2012yve ; CMS:2012qbp . However, most
models of new physics feature extended Higgs sectors: for instance, minimal
supersymmetric models Martin:1997ns  and two Higgs doublet models
Branco:2011iw  feature a second doublet, type-II seesaw models
Schechter:1980gr  feature a zero hypercharge triplet, triplets appear also in
the Georgi-Machacek model Georgi:1985nv , while larger representations appear
in the custodial-preserving septet model Hisano:2013sn . In all these
scenarios, the scalar fields acquire sizeable couplings to the SM gauge bosons
and fermions via VEVs and/or via mixing with the SM Higgs boson. Hence, they
are dominantly singly-produced at colliders, and most current searches are
focused on these channels.

Single production, however, is always model dependent and it can be suppressed
by tuning small single-scalar couplings. In contrast, pair production only
depends on the gauge quantum numbers of the scalars, and cannot be tuned to be
small. The couplings of two scalars from
$\text{SU}(2)_{L}\times\text{U}(1)_{Y}$ multiplets to the EW gauge bosons
arise from the covariant derivatives in the scalar kinetic terms and are,
therefore, always present. They generate the dominant pair production via
Drell-Yan, where two initial state quarks merge via an s-channel gauge boson.
The kinetic term also yields a coupling of two scalars to two gauge bosons
that, via vector boson fusion, contributes to the scalar pair production.
However, this process is subdominant as compared to Drell-Yan
Agugliaro:2018vsu . If the EW symmetry were preserved, the Drell-Yan pair
production cross section for a given scalar multiplet would be a function of
the scalar mass, only. Via the EW breaking, the scalar gauge eigenstates can
mix through the scalar potential, hence the couplings of two (physical) mass
eigenstates to the EW gauge bosons (and thus the Drell-Yan pair production
cross sections) acquire model dependence. Nevertheless, the mass mixing cannot
reduce Drell-Yan pair production cross sections of all scalars at the same
time and some channels are guaranteed to remain sizeable. Other channels may
be present via loops or scalar self couplings: famously, the Higgs boson
itself is pair-produced via the Higgs triple coupling as well as one-loop top
box contributions. In this work, we focus on models where pair production is
the dominant mode for scalars charged under the EW gauge symmetry. This
scenario appears naturally in composite Higgs models, where the Higgs boson is
accompanied by additional light states, protected by parities internal to the
strong sector Ferretti:2016upr .

In composite Higgs models, the Higgs boson emerges as a pseudo-Nambu-Goldstone
boson (pNGB) Kaplan:1983fs  following the dynamical breaking of the EW
symmetry triggered by misalignment in a condensing strong dynamics at the TeV
scale Weinberg:1975gm ; Dimopoulos:1979es . It may well be accompanied by
additional light meson-like scalars. In fact, based only on the global
symmetries, a minimal model $\text{SO}(5)/\text{SO}(4)$ with 4 pNGBs matching
the Higgs doublet components can be constructed Agashe:2004rs  based on
holography Contino:2003ve . However, it is not easy to obtain this symmetry
pattern in an underlying gauge/fermion theory à la QCD. A fermion condensate
$\langle\psi\psi\rangle$ can only generate the following patterns:
$\text{SU}(2N)/\text{Sp}(2N)$, $\text{SU}(N)/\text{SO}(N)$ or
$\text{SU}(N)^{2}/\text{SU}(N)$ depending on whether the representation of
$\psi$ under the confining gauge symmetry is pseudo-real, real or complex
Witten:1983tw ; Kosower:1984aw , respectively. Hence, from the point of view
of the underlying gauge theory, the minimal model with custodial symmetry
Georgi:1984af ; Agashe:2006at  features $\text{SU}(4)/\text{Sp}(4)$
Ryttov:2008xe ; Galloway:2010bp ; Cacciapaglia:2014uja , and has one
additional pNGB besides the Higgs doublet. The next-to-minimal cases contain
many more pNGBs: $15$ for $\text{SU}(6)/\text{Sp}(6)$ Low:2002ws ; Cai:2019cow
, $14$ for $\text{SU}(5)/\text{SO}(5)$ Dugan:1984hq ; Arkani-Hamed:2002ikv
and $15$ for $\text{SU}(4)^{2}/\text{SU}(4)$ Ma:2015gra ; Vecchi:2015fma . For
other non-minimal patterns see Refs Mrazek:2011iu ; Bellazzini:2014yua . Note
that departure from minimality is not in contradiction with the null results
of direct Beyond-the-Standard-Model (BSM) searches at colliders: the pNGBs are
typically heavier than the Higgs and only have EW interactions, hence being
very hard to discover at hadron colliders and too heavy for past $e^{+}e^{-}$
colliders such as LEP. Other resonances, like baryon-like top partners needed
for top partial compositeness Kaplan:1991dc , can be much heavier.

Electroweak pNGBs have recently been studied in the context of exotic decays
of top partners Bizot:2018tds ; Xie:2019gya ; Benbrik:2019zdp ;
Cacciapaglia:2019zmj ; Banerjee:2022xmu , as the latter have sizeable
production cross sections at hadronic colliders like the LHC at CERN. In this
work, we instead focus on the pNGB direct production via their EW couplings.
The dominant channel is pair production via Drell-Yan: The vector boson fusion
(VBF) pair production via gauge couplings is found to be subleading to Drell-
Yan Agugliaro:2018vsu . Single production can also be phenomenologically
relevant, however it is strongly model dependent. VBF single production is
generated via topological anomalies, hence it is suppressed by a small anomaly
coupling. Drell-Yan single production could also be present if the pNGBs
couple to quarks: However, for pNGBs in models with partial compositeness,
couplings to light quark flavours are expected to be very small (roughly
proportional to the quark mass). The dominant couplings, therefore, involve
third generation quarks. In this case, neutral pNGBs can be singly-produced
via gluon fusion analogously to the Higgs boson. Finally, both neutral and
charged pNGBs can be singly-produced in association with $tt$ or $tb$,
respectively, hence providing a relevant contribution if the couplings are
large enough. In this work, we will provide the first complete analysis of how
current LHC searches probe the parameter space of the EW pNGBs via their pair
production. After providing bounds for simplified models, we will focus on a
specific model to investigate the interplay between various channels.

Among the next-to-minimal models, the $\text{SU}(5)/\text{SO}(5)$ model has
been studied since the early days of composite Higgs models Dugan:1984hq . In
the context of four-dimensional models with a microscopic description
Ferretti:2013kya ; Ferretti:2014qta ; Belyaev:2016ftv , it emerges as the
minimal symmetry pattern from the condensate $\langle\psi\psi\rangle$ of two
EW-charged fermions if $\psi$ live in a real irreducible representation
(irrep) of the confining gauge group. Initial investigations of its LHC
phenomenology were performed in Refs Agugliaro:2018vsu ; Banerjee:2022xmu ,
and a detailed description of the model can be found in Ref. Agugliaro:2018vsu
. For other models, we expect similar limits, with the caveat that they may
contain a Dark Matter candidate Wu:2017iji ; Cai:2020njb , while this is not
possible for $\text{SU}(5)/\text{SO}(5)$ due to the topological anomaly. We
conclude the introduction by recalling that the singlet pNGB in the minimal
case $\text{SU}(4)/\text{Sp}(4)$ is very hard to detect due to the small gauge
couplings Galloway:2010bp ; Arbey:2015exa , unless it is lighter than the $Z$
boson Cacciapaglia:2021agf . Finally, we recall that models with top partial
compositeness also contain QCD-coloured pNGBs Cacciapaglia:2015eqa ;
Belyaev:2016ftv ; Cacciapaglia:2020vyf ; Cacciapaglia:2021uqh  and a
ubiquitous, and potentially light, singlet associated with a global U$(1)$
symmetry Ferretti:2016upr ; Belyaev:2016ftv ; Cacciapaglia:2019bqz .

The manuscript is organised as follows: In Sec. 2 we present current bounds on
various production and decay channels of a pair of scalars, which can apply to
any model. In Sec. 3 we focus on the $\text{SU}(5)/\text{SO}(5)$ model and
investigate both the fermiophobic case in Sec. 3.2 and fermiophilic one in
Sec. 3.3. Finally, we offer our conclusions in Sec. 4.

##  2 Simplified model bounds on Drell-Yan pair-produced scalars

Many BSM models contain an extended scalar sector with
$\text{SU}(2)_{L}\times\text{U}(1)_{Y}$ multiplets beyond the Higgs doublet.
The bounds on (and signals of) these models are highly model dependent.
Yukawa-type couplings of the additional scalars are subject to constraints
from flavour physics, while the scalar potential influences the EW symmetry
breaking and is, therefore, strongly constrained. The latter mainly occurs via
VEVs of the new multiplets, while mixing with the Higgs through the scalar
potential can also influence flavour physics. At the same time, Yukawa-type
couplings and scalar VEVs and mixing patterns determine the single production
cross sections of the BSM scalars at lepton and hadron colliders. In the
following we will only focus on pair production, via the dominant Drell-Yan
channels.

###  2.1 Simplified model Lagrangian

For our phenomenological studies, we use parts of a simplified model which has
already been introduced in Banerjee:2022xmu . We extend the SM by colourless
scalar states $S^{0},S^{0\prime},S^{\pm},S^{\pm\pm}$ that are physical mass
eigenstates labelled by their electric charge. We include the minimal set of
states up to charge-2 that have all the possible couplings to the EW gauge
bosons, hence including two neutral states with opposite parity.

We consider the simplified model Lagrangian with kinetic and mass terms for
the scalars as well as interaction terms

$\displaystyle\mathcal{L}_{\rm int}$
$\displaystyle=\mathcal{L}_{SSV}+\mathcal{L}_{SVV}+\mathcal{L}_{ffS}\,,$ (1)

where the first term contains the couplings of two scalars to an EW gauge
boson, which determine the Drell-Yan pair production. The remaining terms
contain the couplings of a scalar to two EW gauge bosons or to two SM
fermions, which dictate the two-body decays into SM particles.

The first term arises from the $\text{SU}(2)_{L}\times\text{U}(1)_{Y}$
covariant derivative terms in full models and reads:

$\displaystyle\mathcal{L}_{SSV}$
$\displaystyle=\frac{ie}{s_{W}}W^{-\mu}\left(K^{S^{0}S^{+}}_{W}S^{0}\overleftrightarrow{\partial_{\mu}}S^{+}+K^{S^{0\prime}S^{+}}_{W}S^{0\prime}\overleftrightarrow{\partial_{\mu}}S^{+}+K^{S^{-}S^{++}}_{W}S^{-}\overleftrightarrow{\partial_{\mu}}S^{++}\right)+{\mathrm{h.c.}}$
$\displaystyle+\frac{ie}{s_{W}c_{W}}Z^{\mu}\left(K^{S^{0}S^{0\prime}}_{Z}S^{0}\overleftrightarrow{\partial_{\mu}}S^{0\prime}+K^{S^{+}S^{-}}_{Z}S^{+}\overleftrightarrow{\partial_{\mu}}S^{-}+K^{S^{++}S^{--}}_{Z}S^{++}\overleftrightarrow{\partial_{\mu}}S^{--}\right)$
$\displaystyle-
ieA^{\mu}\left(S^{+}\overleftrightarrow{\partial_{\mu}}S^{-}+2S^{++}\overleftrightarrow{\partial_{\mu}}S^{--}\right)\,,$
(2)

where
$\phi_{1}\overleftrightarrow{\partial_{\mu}}\phi_{2}\equiv\phi_{1}(\partial_{\mu}\phi_{2})-(\partial_{\mu}\phi_{1})\phi_{2}$.
The $K_{V}^{SS}$ parameters are determined by the
$\text{SU}(2)_{L}\times\text{U}(1)_{Y}$ representations of the scalar
multiplets as well as the mass mixing. The $K_{V}^{SS}$ coefficients for
sample models, including the model discussed in Sec. 3, are given in Appendix
B. The production cross section of each scalar pair is proportional to its
respective $|K_{V}^{SS}|^{2}$.

The second term parameterises dimension-5 operators which yield the decay of
the scalars into two EW gauge bosons, and reads

$\displaystyle\mathcal{L}_{SVV}=\frac{e^{2}}{16\pi^{2}v}\bigg{[}$
$\displaystyle
S^{0}\left(\tilde{K}^{S^{0}}_{\gamma\gamma}F_{\mu\nu}\tilde{F}^{\mu\nu}+\frac{2}{s_{W}c_{W}}\tilde{K}^{S^{0}}_{\gamma
Z}F_{\mu\nu}\tilde{Z}^{\mu\nu}+\frac{1}{s_{W}^{2}c_{W}^{2}}\tilde{K}^{S^{0}}_{ZZ}Z_{\mu\nu}\tilde{Z}^{\mu\nu}\right.$
$\displaystyle\left.\hskip
17.07164pt+\frac{2}{s_{W}^{2}}\tilde{K}^{S^{0}}_{WW}W^{+}_{\mu\nu}\tilde{W}^{-\mu\nu}\right)$
$\displaystyle+$ $\displaystyle
S^{0\prime}\left(K^{S^{0\prime}}_{\gamma\gamma}F_{\mu\nu}F^{\mu\nu}+\frac{2}{s_{W}c_{W}}K^{S^{0\prime}}_{\gamma
Z}F_{\mu\nu}Z^{\mu\nu}+\frac{1}{s_{W}^{2}c_{W}^{2}}K^{S^{0\prime}}_{ZZ}Z_{\mu\nu}Z^{\mu\nu}\right.$
$\displaystyle\left.\hskip
17.07164pt+\frac{2}{s_{W}^{2}}K^{S^{0\prime}}_{WW}W^{+}_{\mu\nu}W^{-\mu\nu}\right)$
$\displaystyle+$
$\displaystyle\left(S^{+}\left(\frac{2}{s_{W}}\tilde{K}^{S^{+}}_{\gamma
W}F_{\mu\nu}\tilde{W}^{-\mu\nu}+\frac{2}{s_{W}^{2}c_{W}}\tilde{K}^{S^{+}}_{ZW}Z_{\mu\nu}\tilde{W}^{-\mu\nu}\right)+\mathrm{h.c.}\right)$
$\displaystyle+$ $\displaystyle
S^{++}\frac{1}{s_{W}^{2}}\tilde{K}^{S^{++}}_{W^{-}W^{-}}W^{-}_{\mu\nu}\tilde{W}^{-\mu\nu}+\mathrm{h.c.}\bigg{]}\,.$
(3)

The couplings above are written assuming that all scalars are odd under
parity, except for the even state $S^{0\prime}$ in order to allow the $Z$
couplings in Eq. (2). This choice is motivated by matching to the composite
models we consider in Sec. 3, however the parity assignment can be flipped in
a straightforward manner. Note that the parity assignment does not
significantly affect the bounds we consider here, as the kinematics of the
decay is untouched. Hence we only study the case in Eq. (3).

The last term contains Yukawa-type couplings to the third generation quarks:

$\displaystyle\mathcal{L}_{ffS}$
$\displaystyle=S^{0}\bigg{[}\bar{t}\left(\kappa^{S^{0}}_{t}+i\tilde{\kappa}^{S^{0}}_{t}\gamma_{5}\right)t+\bar{b}\left(\kappa^{S^{0}}_{b}+i\tilde{\kappa}^{S^{0}}_{b}\gamma_{5}\right)b\bigg{]}+\left(S^{0}\rightarrow
S^{0\prime}\right)$
$\displaystyle+S^{+}\,\bar{t}\left(\kappa^{S^{+}}_{tb,L}P_{L}+\kappa^{S^{+}}_{tb,R}P_{R}\right)b+{\rm
h.c.}\,,$ (4)

where, motivated by the SM structure, the couplings are allowed to violate
parity. Couplings to other SM fermions could be included analogously: our
choice here is motivated by the models of top partial compositeness from Sec
3.

###  2.2 Di-scalar channels

Figure 1: Examples of di-scalar channels from pair production via Drell-Yan
processes with subsequent decays into SM particles.

We investigate all scalar pairs produced at the LHC through the Drell-Yan
processes:

$pp\to
S^{\pm\pm}S^{\mp}\,,\,S^{\pm}S^{0(\prime)}\,,\,S^{++}S^{--}\,,\,S^{+}S^{-}\,,\,S^{0}S^{0\prime}\,.$
(5)

Together with the first tier decays of the scalar pairs into SM particles,
these production processes yield many di-scalar channels, see Fig. 1 for two
examples. Charge-conjugated states belong to the same channel. For the decays
of the scalars, we consider two complementary scenarios: The fermiophobic
case, where the dominant decays are into EW gauge bosons, and the fermiophilic
case, where the scalars decay dominantly into a pair of third generation
quarks. In both cases, we only consider narrow width resonances. The two
choices are motivated by the different origins of the two sets of couplings in
Eq. (3) and Eq. (4): The former deriving from higher dimension operators or
loops, the second from Yukawa-like couplings or (small) mixing to the Higgs
boson.

fermiophobic | ${\color[rgb]{0,0,1}S^{++}}{\color[rgb]{1,0,0}S^{--}}$ | ${\color[rgb]{0,0,1}S^{\pm\pm}}{\color[rgb]{1,0,0}S^{\mp}}$ | ${\color[rgb]{0,0,1}S^{+}}{\color[rgb]{1,0,0}S^{-}}$ | ${\color[rgb]{0,0,1}S^{\pm}}{\color[rgb]{1,0,0}S^{0(\prime)}}$ | ${\color[rgb]{0,0,1}S^{0}}{\color[rgb]{1,0,0}S^{0\prime}}/{\color[rgb]{0,0,1}S^{0\prime}}{\color[rgb]{1,0,0}S^{0}}$
---|---|---|---|---|---
$WWWW$ | ${\color[rgb]{0,0,1}W^{+}W^{+}}{\color[rgb]{1,0,0}W^{-}W^{-}}$ | - | - | - | ${\color[rgb]{0,0,1}W^{+}W^{-}}{\color[rgb]{1,0,0}W^{+}W^{-}}$
$WWW\gamma$ | - | ${\color[rgb]{0,0,1}W^{\pm}W^{\pm}}{\color[rgb]{1,0,0}W^{\mp}\gamma}$ | - | ${\color[rgb]{0,0,1}W^{\pm}\gamma}{\color[rgb]{1,0,0}W^{+}W^{-}}$ | -
$WWWZ$ | - | ${\color[rgb]{0,0,1}W^{\pm}W^{\pm}}{\color[rgb]{1,0,0}W^{\mp}Z}$ | - | ${\color[rgb]{0,0,1}W^{\pm}Z}{\color[rgb]{1,0,0}W^{+}W^{-}}$ | -
$WW\gamma\gamma$ | - | - | ${\color[rgb]{0,0,1}W^{+}\gamma}{\color[rgb]{1,0,0}W^{-}\gamma}$ | - | ${\color[rgb]{0,0,1}W^{+}W^{-}}{\color[rgb]{1,0,0}\gamma\gamma}$
$WWZ\gamma$ | - | - | ${\color[rgb]{0,0,1}W^{\pm}\gamma}{\color[rgb]{1,0,0}W^{\mp}Z}$ | - | ${\color[rgb]{0,0,1}W^{+}W^{-}}{\color[rgb]{1,0,0}\gamma Z}$
$WWZZ$ | - | - | ${\color[rgb]{0,0,1}W^{+}Z}{\color[rgb]{1,0,0}W^{-}Z}$ | - | ${\color[rgb]{0,0,1}W^{+}W^{-}}{\color[rgb]{1,0,0}ZZ}$
$W\gamma\gamma\gamma$ | - | - | - | ${\color[rgb]{0,0,1}W^{\pm}\gamma}{\color[rgb]{1,0,0}\gamma\gamma}$ | -
$WZ\gamma\gamma$ | - | - | - | ${\color[rgb]{0,0,1}W^{\pm}}{\color[rgb]{1,0,0}}\\{Z\gamma\\}{\color[rgb]{0,0,1}}{\color[rgb]{1,0,0}\gamma}$ | -
$WZZ\gamma$ | - | - | - | ${\color[rgb]{0,0,1}W^{\pm}}{\color[rgb]{1,0,0}}\\{Z\gamma\\}{\color[rgb]{0,0,1}}{\color[rgb]{1,0,0}Z}$ | -
$WZZZ$ | - | - | - | ${\color[rgb]{0,0,1}W^{\pm}Z}{\color[rgb]{1,0,0}ZZ}$ | -
$\gamma\gamma\gamma\gamma$ | - | - | - | - | ${\color[rgb]{0,0,1}\gamma\gamma}{\color[rgb]{1,0,0}\gamma\gamma}$
$Z\gamma\gamma\gamma$ | - | - | - | - | ${\color[rgb]{0,0,1}Z\gamma}{\color[rgb]{1,0,0}\gamma\gamma}$
$ZZ\gamma\gamma$ | - | - | - | - | ${\color[rgb]{0,0,1}Z}{\color[rgb]{1,0,0}}\\{Z\gamma\\}{\color[rgb]{0,0,1}}{\color[rgb]{1,0,0}\gamma}$
$ZZZ\gamma$ | - | - | - | - | ${\color[rgb]{0,0,1}ZZ}{\color[rgb]{1,0,0}Z\gamma}$
$ZZZZ$ | - | - | - | - | ${\color[rgb]{0,0,1}ZZ}{\color[rgb]{1,0,0}ZZ}$
Table 1:  Classification of the 24 di-scalar channels in terms of the 5 pair
production cases (columns) and the 15 combinations of gauge bosons (rows) from
decays. In the channels, the first two and second two bosons are resonantly
produced. The notation
$\\{Z\gamma\\}={\color[rgb]{0,0,1}Z}{\color[rgb]{1,0,0}\gamma}+{\color[rgb]{0,0,1}\gamma}{\color[rgb]{1,0,0}Z}$
indicates the two permutations. Charge-conjugated states belong to the same
di-scalar channel.

In the fermiophobic case, we assume dominant decays of the scalars into EW
gauge bosons via the couplings in Eq. (3), leading to111We do not consider the
possible coupling of the neutral scalars to two gluons, as it can only be
generated if they couple to states carrying QCD charges. We remark that Drell-
Yan pair production of scalars with subsequent decays to a pair of dijets is
targeted by experimental searches CMS:2022usq , but public recasts of these
are as of now not available.

$\displaystyle S^{++}$ $\displaystyle\to W^{+}W^{+}\,,$ (6a) $\displaystyle
S^{+}$ $\displaystyle\to W^{+}\gamma,\,W^{+}Z\,,$ (6b) $\displaystyle
S^{0(\prime)}$ $\displaystyle\to W^{+}W^{-},\,\gamma\gamma,\,\gamma Z,\,ZZ.$
(6c)

Combining the different Drell-Yan scalar pairs with the above decay channels
leads to 24 di-scalar channels – each containing four gauge bosons – for which
we present bounds in Sec. 2.4. One sample process is shown in the left diagram
of Fig. 1, while a complete list of all channels is shown in Table 1.

fermiophilic | ${\color[rgb]{0,0,1}S^{++}}{\color[rgb]{1,0,0}S^{--}}$ | ${\color[rgb]{0,0,1}S^{++}}{\color[rgb]{1,0,0}S^{-}}$ | ${\color[rgb]{0,0,1}S^{+}}{\color[rgb]{1,0,0}S^{-}}$ | ${\color[rgb]{0,0,1}S^{+}}{\color[rgb]{1,0,0}S^{0(\prime)}}$ | ${\color[rgb]{0,0,1}S^{0}}{\color[rgb]{1,0,0}S^{0\prime}}/{\color[rgb]{0,0,1}S^{0\prime}}{\color[rgb]{1,0,0}S^{0}}$
---|---|---|---|---|---
$tttt$ | - | - | - | - | ${\color[rgb]{0,0,1}t\bar{t}}{\color[rgb]{1,0,0}t\bar{t}}$
$tttb$ | - | - | - | ${\color[rgb]{0,0,1}t\bar{b}}{\color[rgb]{1,0,0}t\bar{t}}$ | -
$ttbb$ | - | - | ${\color[rgb]{0,0,1}t\bar{b}}{\color[rgb]{1,0,0}b\bar{t}}$ | - | ${\color[rgb]{0,0,1}t\bar{t}}{\color[rgb]{1,0,0}b\bar{b}}$
$tbbb$ | - | - | - | ${\color[rgb]{0,0,1}t\bar{b}}{\color[rgb]{1,0,0}b\bar{b}}$ | -
$bbbb$ | - | - | - | - | ${\color[rgb]{0,0,1}b\bar{b}}{\color[rgb]{1,0,0}b\bar{b}}$
$Wttbb$ | - | ${\color[rgb]{0,0,1}W^{+}t\bar{b}}{\color[rgb]{1,0,0}b\bar{t}}$ | - | - | -
$WWttbb$ | ${\color[rgb]{0,0,1}W^{+}t\bar{b}}{\color[rgb]{1,0,0}W^{-}b\bar{t}}$ | - | - | - | -
Table 2:  Classification of the 8 di-scalar channels in terms of the 5 pair
production cases (columns) and the 5 combinations of top and bottom from
decays (rows). In cases with one or two doubly charged scalars, one always
obtains $ttbb$ with one or two additional $W$’s, respectively. The charge-
conjugated states are not shown.

In the fermiophilic scenario we assume dominant couplings of the scalars to
third family quarks. Note that doubly charged scalars cannot decay to two
quarks due to their charge, but if they are part of an $\text{SU}(2)_{L}$
multiplet, the three-body decay $S^{++}\rightarrow W^{+}S^{+*}\rightarrow
W^{+}t\bar{b}$ is allowed. The dominant decay channels we consider for the
fermiophilic scenario are thus222Note that top and bottom loops generate
effective couplings to gluons and EW gauge bosons, however they lead to
subleading decay channels.

$\displaystyle S^{++}$ $\displaystyle\to W^{+}t\bar{b},$ (7a) $\displaystyle
S^{+}$ $\displaystyle\to t\bar{b},$ (7b) $\displaystyle S^{0(\prime)}$
$\displaystyle\to t\bar{t}\,\,\,\mbox{ or }\,\,\,b\bar{b}.$ (7c)

For pair-produced scalars, this yields 8 possible di-scalar channels in the
fermiophilic scenario. One sample process is shown in the right diagram of
Fig. 1, while a complete list is showcased in Table 2.

###  2.3 Simulation setup and determination of LHC bounds

For the simulation of signal events, we use the publicly available eVLQ model
first presented in Ref. Banerjee:2022xmu , which implements the simplified
models in Eq. 1 as a FeynRules Alloul:2013bka  model at next-to-leading order
in QCD. The implementation contains one doubly charged, one singly charged and
one neutral scalar, and we expanded it by another neutral scalar to allow for
$S^{0}S^{0\prime}$ production.

All events are generated at a centre-of-mass energy of $13$ TeV in proton-
proton collisions. For each di-scalar channel, we perform a parameter scan
over the scalar mass $m_{S}$, starting at the decay mass threshold and up to
$\left.m_{S}\right|_{\rm max}=1$ TeV. For channels involving two different
scalars, we assume them to be mass degenerate. For each scan point, we
generate $10^{5}$ events of Drell-Yan scalar pairs with decay into the target
channel, using MadGraph5_aMC@NLO Alwall:2014hca  version 3.3.2 at NLO
(including patches that were incorporated in version 3.4.0 after the
completion of this work), in association with the parton densities in the
NNPDF 2.3 set Ball:2012cx ; Buckley:2014ana . We then interface the events
with Pythia8 Sjostrand:2014zea  for SM particle decays, showering and
hadronisation. The resulting showered signal events are analysed with
MadAnalysis5 Conte:2012fm ; Conte:2014zja ; Dumont:2014tja ; Conte:2018vmg
version 1.9.60 and CheckMATE Drees:2013wra ; Dercks:2016npn  version 2.0.34
(commit number 8952e7). Both tools reconstruct the events using Delphes 3
deFavereau:2013fsa  and the anti-$k_{T}$ algorithm Cacciari:2008gp
implemented in FastJet Cacciari:2011ma . The exclusion associated with the
events is calculated with the CLs prescription Read:2002hq . We also run the
events against the SM measurements implemented in Rivet Bierlich:2019rhm
version 3.1.5 and extract exclusions from the respective YODA files using
Contur Butterworth:2019wnt ; Buckley:2021neu  version 2.2.1.

To present simplified model bounds, we determine the signal cross section
$\sigma_{95}$ which is excluded at 95% CL. The procedure differs between the
tools:

  * •

MadAnalysis5 explicitly calculates the upper limit on the signal cross section
(both expected sig95exp and observed sig95obs) from the simulated signal
events. We use the observed bound from the signal region to which MadAnalysis5
ascribes the highest sensitivity (best).

  * •

CheckMATE quotes upper limits on the signal, $S_{95}^{\mathrm{exp}}$
(expected) and $S_{95}^{\mathrm{obs}}$ (observed). From these, the input cross
section $\sigma_{\mathrm{in}}$ and the signal $S$ that passed the cuts, we
calculate the upper limit on the cross section as

$\sigma_{95}=\frac{S_{95}^{\mathrm{obs}}}{S}\sigma_{\mathrm{in}}.$ (8)

We follow the default procedure recommended by the CheckMATE collaboration for
determining the best signal region, i.e. we use the observed bound of the
signal region with the strongest expected bound.

  * •

Contur does not calculate bounds on the cross section, so we determine them
manually by running Rivet and Contur multiple times on the same events and
dynamically adjust the input cross section until we obtain
$\mathrm{CL}_{s}=0.05\pm 0.01$.

For each channel and each parameter point, we take the minimal value for
$\sigma_{95}$ obtained from MadAnalysis5, CheckMATE, and Contur as the final
bound, i.e. we do not attempt to combine them.

###  2.4 Simplified model results and discussion

(a) Scalar pair with decays to quarks

(b) $S^{++}S^{--}$ and $S^{\pm\pm}S^{\mp}$ with di-boson decays

(c) $S^{+}S^{-}$ with di-boson decays

(d) $S^{\pm}S^{0}$ with di-boson decays

(e) $S^{0}S^{0\prime}$ with di-boson decays with $\geq 2$ photons

(f) $S^{0}S^{0\prime}$ with di-boson decays with $\leq 1$ photons

Figure 2: Upper limits on the cross section of the di-scalar channels from
Drell-Yan pair production. The scalars decay to: (a) third generation quarks
or (b)-(f) two vector bosons. Both scalars are assumed to have the same mass.
The analyses contributing to the bounds are Refs. ATLAS:2021jol ; CMS:2019xud
; ATLAS:2018nud ; ATLAS:2021mbt ; ATLAS:2021fbt ; ATLAS:2021twp ;
ATLAS:2018zdn ; ATLAS:2019gdh ; CMS:2019xjf ; CMS:2017abv ; ATLAS:2020qlk ;
ATLAS:2021kog ; ATLAS:2019gey ; ATLAS:2018nci ; ATLAS:2016uwq ; CMS:2017moi
(see Tab. 3 and Tab. 4 in Appendix A for details). The numerical values of the
limits are available on https://github.com/manuelkunkel/scalarbounds.

To the best of our knowledge, only very few of these di-scalar channels have
been explicitly targeted by LHC searches. The doubly charged scalar
$S^{\pm\pm}$ is the only one that cannot be singly-produced at a sizable rate,
hence a direct search for its pair production or single production in
association with a singly charged scalar is available ATLAS:2021jol . The
ATLAS search ATLAS:2021jol  can be directly applied to the channels
$S^{++}S^{--}\rightarrow WWWW$ and $S^{\pm\pm}S^{\mp}\rightarrow WWWZ$.
Searches for di-Higgs production also contain final states of interest, like
$WWWW$ and $bbbb$, however those searches strongly focus on scalar masses of
$125$ GeV with production via SM processes CMS:2022cpr ; CMS:2022nmn ;
CMS:2022kdx ; ATLAS:2022ycx  or via a massive resonance ATLAS:2016paq ;
ATLAS:2018hqk . Searches for pair production of two neutral scalars decaying
to $WWWW$ and $bbbb$ final states are also available, but only via resonant
production ATLAS:2018ili ; CMS:2022qww . These searches have limited
applicability to the di-scalar channels discussed in this article because of
the different kinematics.

The di-scalar channels we obtained in Sec. 2.2, however, do populate the
signal region of many BSM searches at the LHC as well as those used for SM
cross section measurements. To determine reliable bounds, it is required to
recast the searches as the topology and/or kinematics can be very different.
As a first step towards determining appropriate bounds we therefore simulate
the processes and determine bounds from all LHC search recasts and
measurements which are publicly available in MadAnalysis5, CheckMATE, and
Contur. The results are showcased in Fig. 2, where we present the simplified
model bounds on the cross section for each of the 32 di-scalar channels. For
each channel, we simulate Drell-Yan produced pairs of bosons with subsequent
decay into the target state (four EW gauge bosons, or four fermions with 0, 1,
or 2 additional $W$ bosons) which are then further decayed, hadronised and
analysed using the procedure described in the previous subsection. Further
details on the dominant analyses are given in Sec. A.2.

Figure 2(a) shows the bounds on the 8 di-scalar channels in the fermiophilic
scenario, consisting of third generation quarks plus one additional $W$ boson
per doubly charged scalar due to the 3-body decay of $S^{\pm\pm}$. In channels
with multiple top quarks, dominant bounds arise from a search for $R$-parity
violating supersymmetry ATLAS:2021fbt , while various supersymmetric searches
CMS:2019xjf ; ATLAS:2018zdn ; ATLAS:2021twp ; ATLAS:2019gdh  and the generic
search in Ref. CMS:2017abv  are relevant for the multi-bottom channels.

Figures 2(b) to 2(f) show the bounds for channels of the fermiophobic scenario
that, for readability, are split into 5 figures and regrouped according to the
charges of the di-scalar states. In case of $S^{0}S^{0\prime}$, the channels
are further sub-grouped according to the number of photons in the final state.
Fig. 2(b) is dedicated to di-scalar channels with at least one doubly-charged
scalar, leading to at least 3 $W$ bosons plus a $W$, $Z$, or photon. The
photon channel $WWW\gamma$ can be constrained using measurements of the
$Z\gamma$ production cross section ATLAS:2019gey ; ATLAS:2018nci . The main
searches for the $WWWW$ and $WWWZ$ channels look for multi-lepton final states
CMS:2019xud ; CMS:2017moi . For these two channels, the results of the ATLAS
search for doubly and singly charged Higgs bosons decaying into vector bosons
in multi-lepton final states ATLAS:2021jol  apply, and they are shown as blue
and orange dashed lines. As is to be expected, the bounds from the ATLAS
search dedicated to these final states are stronger than the bound we obtain
from recasts of a large number of BSM searches targeting different signatures
and scenarios. This also suggests that dedicated searches for the other di-
scalar channels discussed in this article can lead to substantial improvement
in covering their signatures. In Fig. 2(c) we show the di-scalar channels from
$S^{+}S^{-}$ production. The bounds on $W\gamma W\gamma$ are by far the
strongest, coming from a search for gauge-mediated supersymmetry in final
states containing photons and jets ATLAS:2018nud . The main bounds for the
channels $WZWZ$ and $W\gamma WZ$ stem from a multi-lepton search CMS:2017abv
and the $Z\gamma$ cross section measurements ATLAS:2019gey ; ATLAS:2018nci ,
respectively. Fig. 2(d) is dedicated to the di-scalar channels from
$S^{\pm}S^{0}$ production. As for the previous panel, the searches can be
split by the number of photons, leading to the multi-lepton search CMS:2017abv
for channels containing 0 photons, measurements of the $Z\gamma$-cross section
ATLAS:2019gey ; ATLAS:2018nci  for single photon channels and Ref.
ATLAS:2018nud  for multi-photon channels. In Fig. 2(e) we present the
$S^{0}S^{0\prime}$ channels that contain at least 2 photons. The
$\gamma\gamma\gamma\gamma$ channel is constrained by the generic search
ATLAS:2018zdn  and the measurement of the $\gamma\gamma$-production cross
section ATLAS:2021mbt . For the remaining channels, the dominant analysis is
the (multi-)photon search ATLAS:2018nud . Finally, Fig. 2(f) contains the
remaining $S^{0}S^{0\prime}$ channels with at most one photon, which are less
strongly constrained than the multi-photon channels. The main searches
contributing to the bounds are the multi-lepton search CMS:2017abv  and the
$Z\gamma$-cross section measurements ATLAS:2018nci ; ATLAS:2019gey  for the
channels with 0 and 1 photon, respectively.

###  2.5 Applicability and limitations of simplified model bounds

The bounds presented in Fig. 2 are based on recasts of other searches and SM
measurements, apart from the ATLAS direct searches for $S^{\pm\pm}$ production
ATLAS:2021jol . As the most important and generic limitation, we wish to re-
emphasise that our study is based only on searches and measurements by ATLAS,
CMS, and LHCb, for which recasts in MadAnalysis5, CheckMATE, or Contur are
available. This represents only a fraction of (in particular, the newest) LHC
searches and measurements, implying that including recasts of additional
searches will improve the bounds. Performing all these recasts is beyond the
scope of this article.

Another limitation of the simplified model approach stands in the fact that
limits are extracted for each specific channel, in our case applying to 32 (24
fermiophobic and 8 fermiophilic) di-scalar channels. However, realistic models
with an extended Higgs sector contain several scalar mass eigenstates, which
can decay into more than one final state. How can the limits in Fig. 2 be used
to extract reliable limits on a more complex extended Higgs sector? To answer
this question, we consider below three template scenarios, which cover
exhaustively all possibilities.

  1. 1.

_Single scalar, Drell-Yan, several decay channels:_
If only a single particle is produced with a single decay mode, the bounds on
the mass of this particle can be immediately read off from Fig. 2. If the
Drell-Yan produced scalar has several decay channels, it is required to
compute the cross section times branching ratio for each matching di-scalar
channel, and compare them to the corresponding limit in Fig. 2. The most
conservative limit comes from the channel that has the strongest bound. As
different channels may contribute to the same signal region of the leading
search, the actual bound on $m_{S}$ can be further improved by simulating the
complete signal from the scalar pair production.

  2. 2.

_Several scalars, Drell-Yan, several decay channels:_
If the model contains several scalars of similar masses, even more di-scalar
channels can be matched. Besides the most conservative bound described above,
one can extract a more realistic bound by combining various channels. This is
feasible if the scalars are relatively close in mass, so that the acceptances
remain similar. Hence, the procedure would consist of summing the cross
sections times branching ratios of all processes that contribute to the same
di-scalar channel. An even more aggressive approach is to sum all the channels
that contribute to the same search, as we will illustrate with an explicit
example in the next section.

  3. 3.

_Non Drell-Yan and/or new decay channels:_
Finally, there are models where the dominant production is not Drell-Yan, in
which case the limits in Fig. 2 cannot be directly applied. However, the
impact of the different kinematics on the bound is typically limited because
the searches are not dedicated to the specific final state and production
mechanism. Hence, we expect the limits in Fig. 2 to provide a good estimate,
while a full simulation is needed to extract a more reliable bound. The same
consideration applies if additional decay channels are available, like for
instance cascade or three-body decays.

More generally, dedicated searches for the di-scalar channels could give
stronger bounds than the ones we obtained in Fig. 2, and detailed studies are
needed to determine the most promising final states. We leave this
investigation for future work.

##  3 Bounds on the SU(5)/SO(5) pNGBs

The simplified model approach is very useful as the limits can be applied to a
broad class of models, at least to a certain extent. In this section, we
investigate a specific full model with an extended EW scalar sector, study the
bounds on the full model and compare the results to estimates one can very
quickly obtain by using the simplified model approach of Sec. 2.

We focus here on composite Higgs models based on gauge/fermionic underlying
dynamics Ferretti:2013kya ; Ferretti:2014qta ; Ferretti:2016upr . Minimal
models feature one of the following cosets in the EW sector:
$\text{SU}(4)/\text{Sp}(4)$, $\text{SU}(5)/\text{SO}(5)$ or
$\text{SU}(4)\times\text{SU}(4)/\text{SU}(4)_{D}$. A first rough sketch of the
LHC phenomenology of the pNGBs can be found in Ferretti:2016upr . We focus on
the $\text{SU}(5)/\text{SO}(5)$ coset Agugliaro:2018vsu  in the following as
it features a doubly charged scalar.

###  3.1 The electroweak pNGBs and their LHC phenomenology

The pNGBs from the $\text{SU}(5)/\text{SO}(5)$ coset have been investigated in
detail in Ref. Agugliaro:2018vsu  (see also Refs. Dugan:1984hq ;
Ferretti:2014qta ). Here, we summarise the key elements and discuss in some
detail the underlying LHC phenomenology. A complete summary of the pNGB
couplings to vector bosons, which are relevant for this study, can be found in
Ref. Banerjee:2022xmu . The pNGBs of the EW sector form a $\mathbf{14}$ of
$\text{SO}(5)$, which decomposes with respect to the custodial
$\text{SU}(2)_{L}\times\text{SU}(2)_{R}$ as

$\mathbf{14}\to(\mathbf{3},\mathbf{3})+(\mathbf{2},\mathbf{2})+(\mathbf{1},\mathbf{1})\,.$
(9)

We identify the $(\mathbf{2},\mathbf{2})$ with the Higgs doublet (bi-doublet
of the custodial symmetry). Following the notation of Ref. Agugliaro:2018vsu ,
the bi-triplet can be decomposed under the custodial
$\text{SU}(2)_{D}\subset\text{SU}(2)_{L}\times\text{SU}(2)_{R}$ as

$(\mathbf{3},\mathbf{3})\to\mathbf{1}+\mathbf{3}+\mathbf{5}\equiv\eta_{1}+\eta_{3}+\eta_{5}\,,$
(10)

where

$\eta_{1}=\eta_{1}^{0},\quad\eta_{3}=(\eta_{3}^{+},\eta_{3}^{0},\eta_{3}^{-}),\quad\,\eta_{5}=(\eta_{5}^{++},\eta_{5}^{+},\eta_{5}^{0},\eta_{5}^{-},\eta_{5}^{--}).$
(11)

This basis is suggested by the fact that the vacuum of the strong sector
preserves the custodial $\text{SU}(2)_{D}$.333Note that
$(\mathbf{2},\mathbf{2})\rightarrow\mathbf{3}+\mathbf{1}$, where $\mathbf{3}$
are the longitudinal $W$ and $Z$ components and $\mathbf{1}$ is the Higgs
boson. Nevertheless, a mixing among the states is induced by the terms in the
scalar potential that violate it. To simplify the analysis, in the following
we neglect the mixing and assume that the three multiplets have common masses
$m_{1}$, $m_{3}$ and $m_{5}$, respectively. Mass differences are due to the EW
symmetry breaking, hence one naively expects a relative mass split of the
order $v/m_{i}$ ($i=1,3,5$) where $v$ is the VEV of the Higgs boson. The
precise values depend on the details of the scalar potential: here, we
consider the mass differences as free parameters, and allow them to vary up to
$200$ GeV. Besides the masses, there is an additional parameter that is
important for the phenomenology: $\sin\theta=v/f_{\psi}$, with $f_{\psi}$
being the decay constant of this pNGB sector. Electroweak precision data give
a lower bound of about 1 TeV on $f_{\psi}$ Agugliaro:2018vsu . Last but not
least, we assume that the vacuum is only misaligned along the Higgs direction
in order to avoid large breaking of the custodial symmetry. We remark that,
while the EW quantum numbers of the scalars are similar to those of the
Georgi-Machacek model Georgi:1985nv , all states in Eq. 11 are parity-odd,
except for $\eta_{3}$ which is parity even. Hence, in the composite model only
the custodial triplet can develop a VEV without breaking CP.

In composite Higgs models with an extended pNGB sector, there are three types
of couplings that determine the phenomenology of the scalars:

  * (i)

Gauge interactions due to the EW quantum numbers of the pNGBs. In absence of
VEVs, they lead to couplings of two scalars with one (and two) gauge boson(s),
along the lines of Eq. 2. For the $\text{SU}(5)/\text{SO}(5)$ coset, a
complete list of these couplings is reported in Ref. Banerjee:2022xmu .

  * (ii)

Couplings of one pNGB to two EW gauge bosons generated by the topological
anomaly of the coset, in the form of Eq. 3. They correspond to dimension-5
operators and are suppressed by one loop. For $\text{SU}(5)/\text{SO}(5)$, the
coefficients are listed in Refs. Dugan:1984hq ; Ferretti:2014qta . Note that
the parity-even state $\eta^{0}_{3}$ lacks these couplings, gluons do not
appear as the underlying fermions are only charged under the EW symmetry, and
the model dependence is contained in a pre-factor that depends on the gauge
group of the underlying confining dynamics.

  * (iii)

Couplings of one pNGB to SM fermions, in the form of Eq. 4, where only top and
bottom appear following top partial compositeness. These couplings depend on
the properties of the top partners, and they are classified in Ref.
Agugliaro:2018vsu .

The couplings (i) are responsible for Drell-Yan pair production, which
dominate as (ii) and (iii) lead to very small cross sections. The cross
sections of all pNGB pairs as a function of a common mass are shown in Fig. 3,
which include a K-factor of 1.15 arising from QCD corrections Fuks:2013vua .
Finally, all types of couplings determine the decay patterns of the scalar
pair. We illustrate an example in Fig. 4. Besides the cascade decays, which
are relevant for large enough mass splits between multiplets, the final states
match the di-scalar channels discussed in Sec. 2. In particular, when
couplings to fermions are present, they tend to dominate over the decays to
gauge bosons.

Figure 3: Cross sections for the Drell-Yan production of
$\text{SU}(5)/\text{SO}(5)$ pNGBs at the LHC with $\sqrt{s}=13\,$TeV, assuming
the same mass for all states of the custodial singlet, triplet, and
quintuplet. Note that the $\eta_{1}^{0}\eta_{5}^{0}$ combination is not
allowed as they are both parity-odd.

The LHC signatures of pNGB pair production depend strongly on whether the
pNGBs are fermiophilic or fermiophobic. We start the discussion with the
fermiophobic case, in which case interactions to the EW gauge bosons are
relevant. The corresponding branching ratios are shown in Figs. 5 and 6. For
the lightest multiplet and near-degenerate masses, the anomaly couplings
determine decays into a pair of EW gauge bosons, with the exception of
$\eta_{3}^{0}$. At the leading order in $v/f_{\psi}$, only decays involving
neutral gauge bosons appear.444This is due to the fact that the only gauge-
invariant operator appears for the neutral triplets,
$\phi^{a}W^{a}_{\mu\nu}\tilde{B}^{\mu\nu}$, where $B$ contains the hypercharge
gauge boson. Couplings with only $W^{\pm}$ need two insertion of the Higgs
VEV, hence they are suppressed by $v^{2}/f_{\psi}^{2}$. Hence, the singly
charged states decay as

$\eta_{3,5}^{+}\to W^{+}\gamma,\,W^{+}Z\,,$ (12)

with dominant photon channel as $\mathrm{Br}(\eta_{3,5}^{+}\to
W^{+}\gamma)\approx\cos^{2}\theta_{W}\approx 78\%$ Agugliaro:2018vsu  for both
multiplets, as shown in Figs. 5(c), 6(a) and 6(b) for small mass split. The
neutral singlet and quintuplet can decay as

$\eta_{1,5}^{0}\to\gamma\gamma,\,\gamma Z,\,ZZ\,,$ (13)

with comparable branching ratios, see for example the $\eta^{0}_{1}$ decays in
Fig. 5(a) and the $\eta^{0}_{5}$ decays in Fig. 5(d) for small mass split.
Couplings to charged $W^{\pm}$ are suppressed by the EW scale, hence they lead
to $\eta_{1,5}^{0}\to W^{+}W^{-}$ branching ratios suppressed by
$(v/f_{\psi})^{4}\raise
1.29167pt\hbox{$\;<$\kern-7.5pt\raise-4.73611pt\hbox{$\sim\;$}}10^{-3}$, which
we neglect. Instead, while still suppressed, this provides the only available
decay channel for the doubly charged pNGB in the quintuplet:

$\eta_{5}^{++}\to W^{+}W^{+}.$ (14)

Finally, the $\eta_{3}^{0}$ is CP-even and thus has no couplings to the
anomaly. It therefore undergoes three-body decays via off-shell pNGBs:

$\displaystyle\eta_{3}^{0}\to W^{+}W^{-}\gamma,\,W^{+}W^{-}Z\quad$
$\displaystyle\text{ via }\eta_{3,5}^{\pm(*)}\,,\mbox{ and }$ (15a)
$\displaystyle\eta_{3}^{0}\to Z\gamma\gamma,\,ZZ\gamma,\,ZZZ\quad$
$\displaystyle\text{ via }\eta_{1,5}^{0(*)}\,.$ (15b)

These processes contribute to the upper tier in Fig. 4. There is an
interesting cancellation taking place in the three-body decays: In the limit
$\theta\to 0$, the contributions to Eq. 15a cancel exactly if $m_{3}=m_{5}$.
The same holds for Eq. 15b if $m_{1}=m_{3}=m_{5}$. Thus, if the pNGBs are
mass-degenerate, the $\eta_{3}^{0}$ becomes rather long-lived and leaves the
detector before it decays. In practice, however, we expect at least a small
split, so $\eta_{3}^{0}$ decays promptly to three vector bosons. The main
effect on the phenomenology is that the decays through the charged channel Eq.
15a are suppressed if $m_{1}\gg m_{5}\gtrsim m_{3}$, which we explore further
in Sec. 3.2 and which is illustrated in Fig. 6(d).

Figure 4: Example of Drell-Yan production of two pNGBs with cascade and
anomaly decays. If the triplet is the lightest multiplet, the $\eta_{3}^{0}$
undergoes three-body decays via off-shell pNGBs.

(a) Decays of $\eta_{1}^{0}$ for $m_{1}=600~{}\mathrm{GeV}>m_{3}$

(b) Decays of $\eta_{5}^{++}$ for $m_{5}=600~{}\mathrm{GeV}>m_{3}$

(c) Decays of $\eta_{5}^{+}$ for $m_{5}=600~{}\mathrm{GeV}>m_{3}$

(d) Decays of $\eta_{5}^{0}$ for $m_{5}=600~{}\mathrm{GeV}>m_{3}$

Figure 5: Overview of the pNGB decays in the fermiophobic case. The mass of
the decaying particles is set to 600 GeV. The heavier state decays either via
the anomaly into di-boson final states or via an (off-shell) gauge boson into
a lighter pNGB.

(a) Decays of $\eta_{3}^{+}$ for $m_{5}\gg m_{3}=600~{}\mathrm{GeV}>m_{1}$

(b) Decays of $\eta_{3}^{+}$ for $m_{1}\gg m_{3}=600~{}\mathrm{GeV}>m_{5}$

(c) Decays of $\eta_{3}^{0}$ for $m_{5}\gg m_{1}>m_{3}=600$ GeV

(d) Decays of $\eta_{3}^{0}$ for $m_{1}\gg m_{5}>m_{3}=600$ GeV

Figure 6: Overview of the pNGB decays in the fermiophobic case (continued from
Fig. 5). The neutral triplet component decays into three gauge bosons, as it
does not couple to the anomaly.

The discussion so far applies to the lightest multiplet and also covers the
case where the multiplets are very close in mass. However, there can be a
sizeable mass split, in which case cascade decays from one multiplet into a
lighter one and a (potentially off-shell) vector boson become important.
Assuming for example $m_{5}>m_{3}>m_{1}$, we have

$\displaystyle\eta_{5}^{++}\to W^{+(*)}\eta_{3}^{+},\qquad\eta_{5}^{+}\to
Z^{(*)}\eta_{3}^{+},\,W^{+(*)}\eta_{3}^{0},\qquad\eta_{5}^{0}\to
W^{\pm(*)}\eta_{3}^{\mp},\,Z^{(*)}\eta_{3}^{0},$ (16a)
$\displaystyle\eta_{3}^{+}\to W^{+(*)}\eta_{1}^{0},\qquad\eta_{3}^{0}\to
Z^{(*)}\eta_{1}^{0}.$ (16b)

We find that both classes of decays are of similar importance once the mass
split is between $30$ and $50$ GeV, see Figs. 5 and 6, while cascade decays
dominate for larger mass splits. The two exceptions to this rule of thumb are
$\eta_{5}^{++}$ as shown in Fig. 5(b), whose anomaly coupling is suppressed by
$v^{2}/f_{\psi}^{2}$, and $\eta_{3}^{0}$, for which the anomaly-induced three-
body decays are irrelevant as soon as any cascade decay is accessible. We
note, for completeness, that the quintuplet does not couple to the singlet in
the model considered.

We turn now to the fermiophilic case. We assume here that only couplings to
quarks are present. One expects that the couplings in Eq. 4 scale like the
quark masses, e.g.

$\displaystyle\kappa^{\eta^{0}_{i}}_{t}=c_{t}^{i}\,\frac{m_{t}}{f}\quad,\quad\kappa^{\eta^{0}_{i}}_{b}=c_{b}^{i}\,\frac{m_{b}}{f}\quad\text{and}\quad\kappa^{\eta^{+}_{j}}_{tb}=c_{tb}^{j}\,\frac{m_{t}}{f}\,,$
(17)

where the $c$ coefficients are of order one. In this case the decays to third
generation quarks dominate over the loop-level anomaly-induced decays into two
vector bosons or the three-body decays discussed above. Hence, we consider for
this scenario the decays

$\displaystyle\eta_{3,5}^{+}\to t\bar{b},\qquad\eta_{1,3,5}^{0}\to
t\bar{t},\,b\bar{b}\,.$ (18)

From Eq. 17, the $t\bar{t}$ channel dominates over $b\bar{b}$ above threshold.
In the case of $\eta^{++}_{5}$, it turns out that the three-body decay

$\eta_{5}^{++}\to W^{+}t\bar{b}$ (19)

via an off-shell $\eta_{3,5}^{+}$ is dominant over the decay to $W^{+}W^{+}$.
In case of $m_{5}>m_{3}$ also the decay $\eta^{++}_{5}\to
W^{+(*)}\eta^{+}_{3}$ becomes important. We have checked that for mass
differences below $25$ GeV the decay into quarks clearly dominates and for a
mass difference of $50$ GeV the modes $W^{+}t\bar{b}$ and
$W^{+(*)}\eta^{+}_{3}$ are of equal importance. For larger mass differences
the latter mode is the most important one. Here we have assumed that the
coefficients $c$ are equal to one.

###  3.2 LHC bounds in the fermiophobic case

(a) Bounds from individual channels

(b) Bounds from sum of multiphoton channels

Figure 7: Application of the model-independent bounds to a specific model, the
custodial quintuplet $\eta_{5}$ from the $\text{SU}(5)/\text{SO}(5)$ coset. In
(a) we determine the bounds from the dominant individual channels by comparing
the cross section time branching ratio from the model (solid) with the upper
limits from Fig. 2 (dashed). In green we show the results of a full
simulation. The blue line in (b) is the sum of the individual multi-photon
cross sections shown in (a). Further details are given in the text.

As a first step, we consider only the quintuplet $\eta_{5}$ and apply the
simplified model bounds from Sec. 2, where we found that final states with
multiple photons and at least one $W/Z$ yield the strongest constraints. In
Fig. 7(a) we compare the cross section times branching ratio of all multi-
photon final states (solid lines) with the corresponding bounds from Fig. 2
(dashed lines). From the individual channels we find that masses below $340$
GeV are excluded, with the strongest bound coming from
$\eta_{5}^{\pm}\eta_{5}^{0}\to W\gamma\gamma\gamma$. In addition, we perform a
full simulation in which all states contained in the quintuplet are pair-
produced and decayed according to the specific model under study. The solid
green line denotes the sum over all pair production cross sections of the
quintuplet. The dashed green line shows the corresponding bound, i.e. the sum
of scalar pair production cross sections that would be needed in order to
exclude the convolution of all decay channels from quintuplet states. As can
be seen, the bound on the mass $m_{S}$ is $485$ GeV and thus significantly
stronger than the bounds obtained from individual channels. The apparent
discrepancy between simplified models and the full simulation stems from the
fact that all multi-photon channels populate the same signal region of the
search ATLAS:2018nud  that yields the dominant bound. Also, all multi-photon
channels have a similar upper limit, indicating that the signal acceptances
are comparable. Adding up the various signal cross sections with two or more
photons results in the blue line shown in Fig. 7(b). Comparing this summed
cross section with the bounds from different multi-photon channels (see the
shaded area in Fig. 7(b)) yields an estimated bound on $m_{S}$ of $460-500$
GeV, in agreement with the result of the full simulation. This example shows
the usefulness (and limitations) of the simplified model bounds and how they
can be combined in the context of a particular model.

In a second step, we take all multiplets into account and consider scenarios
with fixed mass differences. We study the following benchmark scenarios,
characterised by varying a single mass scale $m_{S}$:

S-eq: $\displaystyle m_{3}=m_{S}-2~{}\mathrm{GeV},\quad$ $\displaystyle
m_{5}=m_{S},\quad$ $\displaystyle m_{1}=m_{S}+2~{}\mathrm{GeV}\,;$ (20a)
S-135: $\displaystyle m_{1}=m_{S}-50~{}\mathrm{GeV},\quad$ $\displaystyle
m_{3}=m_{S},\quad$ $\displaystyle m_{5}=m_{S}+50~{}\mathrm{GeV}\,;$ (20b)
S-531: $\displaystyle m_{5}=m_{S}-50~{}\mathrm{GeV},\quad$ $\displaystyle
m_{3}=m_{S},\quad$ $\displaystyle m_{1}=m_{S}+50~{}\mathrm{GeV}\,;$ (20c)
S-351: $\displaystyle m_{3}=m_{S}-50~{}\mathrm{GeV},\quad$ $\displaystyle
m_{5}=m_{S},\quad$ $\displaystyle m_{1}=m_{S}+50~{}\mathrm{GeV}\,.$ (20d)

The choice of $50$ GeV is motivated by the fact that the mass splits are
expected to be a fraction of the Higgs VEV. The phenomenology differs in each
case: In S-eq, all particles decay via the anomaly and $\eta_{3}^{0}$ exhibits
three-body decays. We introduce a small mass split of 2 GeV to avoid the
cancellation for some $\eta_{3}^{0}$ decays discussed below Eq. 15. In S-135
and S-531, the heavier states decay into the next lighter states or di-bosons,
while the lightest states only have anomaly decays. Finally, in S-351 both
$\eta_{1}$ and $\eta_{5}$ decay into the triplet, and $\eta_{3}^{0}$ decays
into three vector bosons.

(a) Scenario S-eq

(b) Scenario S-135

(c) Scenario S-531

(d) Scenario S-351

Figure 8: Bounds on the pNGB masses for the Drell-Yan production of the full
bi-triplet for multiple benchmark mass spectra defined in Eq. 20. In (a), all
masses are approximately equal. In the remaining panels, there is a 50 GeV
mass split between the multiplets.

We present the bounds on the mass parameter $m_{S}$ for the four benchmark
scenarios in Fig. 8. In orange, we show the sums over all scalar pair
production cross sections $\sigma_{95}$ that would be needed to exclude the
model at 95% CL at each parameter point. As discussed above, the strongest
bounds come from multi-photon channels, with Ref. ATLAS:2018nud  being the
dominant analysis, cf. Tab. 4 in Sec. A.2. The kink in $\sigma_{95}$ is due to
a change in dominant signal region within the same analysis. The actual sum
over all pair production cross sections is drawn in blue. The bounds range
from $640$ GeV for S-135 to $720$ GeV for S-153. The case S-eq can be
understood by adding the additional channels due to the triplet and using the
same procedure as in case of the pure quintuplet. The fact that the
$\eta^{0}_{3}$ decays only via three-body modes is of lesser importance for
final states containing photons. The different bounds for the other scenarios
considered are due to the relative size of the cross section for the triplet
and quintuplet. In the case where the quintuplet is heavier than the triplet,
the decay $\eta^{++}_{5}\to W^{+*}\eta^{+}_{3}$ leads to additional photons
stemming from the $\eta^{+}_{3}$ decays that increase the bound compared to
the scenarios in which $\eta^{++}_{5}$ decays only into $W^{+}W^{+}$.

Figure 9: Bounds on the pNGB masses for the Drell-Yan production of the
custodial triplet $\eta_{3}$ and singlet $\eta_{1}$ with the quintuplet
$\eta_{5}$ decoupled (scenario S-31). Depending on the mass hierarchy, the
pNGBs decay either into di-bosons or into one vector boson and a lighter pNGB.
The heatmap and the dotted contours show the total cross section. The bounds
are obtained from Ref. ATLAS:2018nud , with the dominant signal region
indicated by the marker symbol. The 95% and 68% CL exclusion contours are
drawn in solid black and gray, respectively.

Finally, we consider a third case where one of the multiplets is effectively
decoupled, and define two benchmarks:

$\text{S-31}:\quad m_{5}\gg m_{3,1}\,;\qquad\qquad\text{S-35}:\quad m_{1}\gg
m_{3,5}\,.$ (21)

The case $m_{3}\gg m_{1,5}$ is already covered by our first example of this
section since the singlet and quintuplet do not couple and only the quintuplet
members are produced via Drell-Yan processes. For both scenarios, we scan over
the two light masses with a mass split of up to $200$ GeV and simulate the
Drell-Yan production of two pNGBs. In Fig. 9, we show the results for S-31 in
the $m_{3}$-$\Delta m_{13}$ plane, where $\Delta m_{13}=m_{1}-m_{3}$. In
addition to the exclusion contours at 95% CL (solid black) and 68% CL (solid
gray), we also show the sum over pair production cross sections as a heatmap
with dotted contours. This highlights interesting features in the form of
regions where the bounds deviate from the cross section contours. Following
the 95% CL bound, we identify three such regions: In the lower half, the
triplets decay to the singlet and the final state is determined by the anomaly
decays of $\eta_{1}^{0}$, see Fig. 5(a). From $\Delta m_{13}=-200$ GeV to
$-100$ GeV, the bounds grow weaker as the $W$ and $Z$ bosons from
$\eta_{3}$-cascade decays get softer, followed by an increase towards $\Delta
m_{13}=0$ as the $\eta_{3}^{+}\to W^{+}\gamma$ decay sets in. Finally, when
the singlet is heavier than the triplet, the bounds are weaker again due to
the decreasing $\mathrm{Br}(\eta_{3}^{0}\to\gamma\gamma Z)$.

Figure 10: Bounds on the pNGB masses for the Drell-Yan production of the
custodial triplet $\eta_{3}$ and quintuplet $\eta_{5}$ with the singlet
$\eta_{1}$ decoupled (scenario S-35). Depending on the mass hierarchy, the
pNGBs decay either into di-bosons or into one vector boson and a lighter pNGB.
The heatmap and the dotted contours show the sum over all scalar pair
production cross sections. The bounds are obtained from Ref. ATLAS:2018nud ,
with the dominant signal region indicated by the marker symbol. The 95% and
68% CL exclusion contours are drawn in solid black and grey, respectively.

In Fig. 10, we show the bounds on $m_{3}$ and $m_{5}$ with the singlet
decoupled, S-35. To understand the features, we again follow the 95% CL
exclusion contour. For negative $\Delta m_{53}=m_{5}-m_{3}$, the quintuplet
states decay via the anomaly. The $\eta_{3}^{+}$ dominantly decays into the
$\eta_{5}^{++}$, which cannot contribute photons to the final state. Thus, the
bounds increase relative to the cross section from $\Delta m_{53}=-50$ GeV as
the $\eta_{3}^{+}$ anomaly decays become relevant. For a positive mass split,
the bounds rapidly increase. The reason for this is that the $\eta_{5}^{++}$,
which is produced with a large cross section, has a very small anomaly
coupling so that already at $25$ GeV mass split the branching ratio of the
cascade decay to $\eta_{3}^{+}$ is almost 100%, resulting in a large photon
production. With increasing mass split, the bounds become weaker again. This
is mainly due to the dependence of the $\eta^{0}_{3}$ decays on the mass
difference, see Fig. 6(c). In the nearly mass degenerate case the decays into
$W$ bosons are strongly suppressed, leading to an enhancement of photons from
the $\eta^{0}_{3}$ decays.

###  3.3 LHC bounds in the fermiophilic case

(a) $S^{+}tb\to tbtb$

(b) $S^{0}tt\to tttt$

(c) $S^{0}\to tt$

(d) $S^{0}\to\gamma\gamma$

Figure 11: Bounds on the single production of pNGBs in association with two
third- generation quarks. The blue shaded area indicates typical cross
sections assuming that only one scalar is present. The coupling to the quarks
is given by $c\,\frac{m_{t}}{f_{\psi}}$ and the side band maps the blue shade
to the corresponding value of $c$ (for $f_{\psi}=1$ TeV) and $f_{\psi}$ (for
$c=1$). From a model-building point of view, it is always possible to choose a
larger scale, $f_{\psi}>1$ TeV. Hence, as $c\sim\mathcal{O}(1)$, the darker
shaded regions correspond to theoretically more favourable parameter points
than the light and white regions.

We turn now to the scenarios in which the pNGBs couple to quarks, where decays
via the anomaly are strongly suppressed and can be neglected, as already
discussed in Sec. 2. In these scenarios, one has single scalar production via
the processes

$\displaystyle pp\to S^{0}t\bar{t}\quad\text{ and }\quad pp\to S^{\pm}tb\,.$
(22)

Moreover, the couplings of the neutral scalars to bottom and top quarks induce
couplings to gluons and photons at the one-loop level. This leads to processes
like

$\displaystyle pp\to S^{0}\to t\bar{t}\quad\text{ and }\quad pp\to
S^{0}\to\gamma\gamma\,.$ (23)

We show in Fig. 11 bounds on various processes for $f_{\psi}=1$ TeV and three
different values of the factors $c$ defined in Eq. (17): 1/5, 1 and 5. Note,
that different values of $c$ and $f_{\psi}$ can be obtained by a simple
rescaling of the $c=1$ line by a factor $(c/f_{\psi})^{2}$, with $f_{\psi}$ in
TeV. We compare available searches for $pp\to S^{\pm}tb\to\bar{t}bt\bar{b}$
CMS:2019rlz ; ATLAS:2021upq , $pp\to H/At\bar{t}\to t\bar{t}t\bar{t}$
ATLAS:2022ohr , $pp\to Z^{\prime}\to t\bar{t}$ ATLAS:2018rvc ; CMS:2018rkg
and $pp\to S^{0}\to\gamma\gamma$ ATLAS:2021uiz , from which we extract a limit
on the respective signal cross section. For Fig. 11(a) we use the
renormalisation and factorisation scales
$\mu_{R}=\mu_{F}=(m_{t}+m_{b}+m_{S^{+}})/3$ as this gives a K-factor very
close to 1 Dittmaier:2009np . For the other plots we have taken the cross
sections from the Higgs Xsection working group Higgs:Xsec ;
LHCHiggsCrossSectionWorkingGroup:2016ypw  and have rescaled the Yukawa
couplings accordingly. We see that currently we do not get any bounds except
for $c=5$ and $f_{\psi}=1$ TeV in the 4 $t$ channel, Fig. 11(b), which gives a
bound of about 640 GeV on $m_{S^{0}}$. This corresponds to a rather small
fraction of the available parameter space and if one reduces
$(c=5)/(f_{\psi}=1$ TeV) by a factor $\simeq 1/\sqrt{3}$ one does not get any
bound.

Figure 12: Bounds on the pNGB masses for the Drell-Yan production of the full
bitriplet with decays to third-generation quarks.

We now turn to Drell-Yan pair production, for which we give our results in
Fig. 12. Here we have assumed that all pNGBs have the same mass and all
factors $c=1$ (neither branching ratios nor production cross sections depend
on $f_{\psi}$). The blue line gives the total cross section summing over all
pNGBs irrespective of their decay modes. The orange lines give the exclusion
when considering all possible channels. They are dominated by Ref.
ATLAS:2021fbt  implemented in CheckMATE. Note that CheckMATE uses the signal
region with the strongest expected bound and reports the corresponding
observed bound as the final result. Using this standard procedure, one obtains
the bound given by the solid orange line. However, this can lead to
difficulties if observed and expected bounds differ significantly leading to
the kinks at $m_{S}=350$ GeV and $450$ GeV. Modifying the procedure such that
always the strongest observed bound is taken, one obtains a smoother curve for
the limit, shown by the dashed orange line. This yields a somewhat stronger
bound of about $500$ GeV. We detail the differences of these procedures in
Sec. A.1.

##  4 Conclusions and outlook

In this work, we investigate the bounds on the Drell-Yan pair production of
scalar bosons that carry electroweak charges at the LHC. We first consider all
possible channels in a simplified model approach, leading to 32 distinct
channels: 24 containing four vector bosons, and 8 with top and bottom quarks.
The two scenarios arise from fermiophobic and fermiophilic models,
respectively. The only channels that have dedicated searches contain doubly
charged scalars decaying into a pair of same-charge $W$ bosons. For other
channels, we use all the available recast searches for new physics and
measurements of SM cross sections. These limits, showcased in Fig. 2, can be
applied to any model with an extended Higgs sector dominated by pair
production.

As a concrete example, we focus on a composite Higgs model based on the coset
$\text{SU}(5)/\text{SO}(5)$, which features a custodial bi-triplet. We show
that the limits on individual channels lead to relatively weak bounds on the
scalar masses. Instead, stronger bounds can be obtained by combining various
pair production channels. Considering several benchmark scenarios, we
establish limits on the scalar mass scale around $500-700$ GeV in the
fermiophobic case. For decays into top and bottom quarks, the bounds are
around $500$ GeV.

The main limitation of the simplified model approach is the restriction to
searches and measurements that have been recast. By determining limits from
MadAnalysis5, CheckMATE and Rivet/Contur, we cover a considerable amount of
analyses. Still, there are many searches, not yet implemented, that have the
potential to significantly improve these bounds. Another limitation is in the
combination of different searches, which is not possible without detailed
knowledge of the experimental correlations between the various signal regions.
Designing simple combination procedures, like the one proposed in Ref.
Araz:2022vtr , could mitigate this issue. Furthermore, the combination removes
the ambiguity in choosing the most sensitive signal region.

Within its limitations, our analysis proves that current non-dedicated
searches and standard model measurements impose significant bounds on extended
Higgs sectors, which contain many scalar bosons with electroweak charges.
Nevertheless, the variety of production channels and available final states
leaves open the possibility to improve the coverage of this large class of
models by means of dedicated searches. In the fermiophobic case, final states
with photons are remarkable, as they already lead to the best bounds via
generic multi-photon searches. The reach could be improved by searches
targeting photon resonances or kinematic features related to the mass of the
decaying scalar bosons. These channels are particularly relevant for composite
models, where final states with photons are naturally abundant. In the
fermiophilic case, final states with multi tops could be targeted by resonance
searches. A dedicated experimental search programme by ATLAS, CMS, and LHCb
could immediately improve the experimental coverage of extended Higgs sectors.

The reach can clearly be extended further by prospective future colliders. We
show in Fig. 13 the cross section of a typical process, the pair production of
the doubly charged $\eta^{++}_{5}$, for various colliders: Besides the High-
Luminosity LHC at $14$ TeV, we consider a $100$ TeV proton-proton collider and
a muon-collider with a centre-of-mass energy of $3$, $10$, and $14$ TeV. For a
$100$ TeV $pp$-collider with a typical integrated luminosity of $30$ ab-1
Hinchliffe:2015qma ; Arkani-Hamed:2015vfh , by naive re-scaling we estimate
its reach to cover pNGB masses up to $4$ TeV. For the various muon-collider
options, the reach should be close to $m_{S}\sim\sqrt{s}/2$ assuming that the
integrated luminosity scales as $(\sqrt{s}/(10~{}\mathrm{TeV}))^{2}\ 10^{4}$
ab-1 Delahaye:2019omf ; Han:2022edd . Dedicated studies will be necessary to
obtain more realistic values for the reach of the different collider options.

Figure 13: Cross sections for the $\eta^{++}_{5}\eta^{--}_{5}$ pair production
at the $14$ TeV HL-LHC and some of the proposed future collider options.

## Acknowledgement

We thank A. Banerjee, J. Butterworth, G. Ferretti, K. Rolbiecki and R.
Ströhmer for useful discussions. This work has been supported by the “DAAD,
Frankreich” and “Partenariat Hubert Curien (PHC)” PROCOPE 2021-2023, project
number 57561441 as well as the international cooperation program “GenKo”
managed by the National Research Foundation of Korea (No. 2022K2A9A2A15000153,
FY2022) and DAAD, P33 - projekt-id 57608518. G.C. and T.F. acknowledge support
from the Campus-France STAR project n. 43566UG, “Higgs and Dark Matter
connections”. T.F. is supported by a KIAS Individual Grant (AP083701) via the
Center for AI and Natural Sciences at Korea Institute for Advanced Study.

##  Appendix A Technical Notes

###  A.1 Choosing the best signal region

When choosing the most sensitive signal region for a given analysis, CheckMATE
uses the signal region with the strongest expected bound but reports the
corresponding observed bound as the final result. This can lead to some
unintuitive results when there is a large difference between the expected and
observed bound, such as the sudden increase in $\sigma_{95}$ in Fig. 12 at
$m_{S}=450$ GeV. To illustrate what causes this behaviour, we show the bounds
from all relevant signal regions in Fig. 14: The dominant analysis for the
decays to quarks is Ref. ATLAS:2021fbt , and there are three important signal
regions, SR11 (blue), SR13 (orange) and SR15 (green). For each signal region,
we show the observed bounds as a solid line and the expected bounds as a
dashed line. The black dotted line indicates the “strongest” bound
$\sigma_{95}$ using the default method described above, while the grey dotted
line is the $\sigma_{95}$ that is obtained by choosing the minimum of the
observed bounds for each parameter point.

When there is one signal region that clearly dominates, such as SR11 for small
masses, the default and minimum procedures coincide. However, for $m_{S}\geq
450$ GeV, the expected bounds from SR13 and SR15 are very similar with SR15
being marginally more sensitive. The default procedure then dictates using the
observed bounds from SR15 for $\sigma_{95}$, although they are significantly
weaker than the ones from SR13. Given that the difference in the expected
significance is small, we find it justified to use SR13 instead.

Figure 14: Bounds on the pNGB masses for the Drell-Yan production of the full
bitriplet with decays to third-generation quarks.

###  A.2 List of dominant analyses

In Fig. 2 in the main text, we present upper limits on the Drell-Yan
production cross section of electroweak scalars for a variety of decay
channels. Due to the different topologies of the resulting final states, the
analyses that yield the strongest constraints differ among the various
channels. In this appendix we break down which analyses contribute to which
decay channel. Tab. 3 gives a brief description of the relevant analyses,
including the recasting tool they are implemented in and their respective
tool-internal name. In Tab. 4, we then list for each channel the analyses that
give the dominant bound for at least one mass point. The full information is
available on https://github.com/manuelkunkel/scalarbounds.

Analysis | Description | Recast
---|---|---
|  ATLAS JHEP ATLAS:2021jol
---
$139~{}\mathrm{fb}^{-1}$
|  $S^{++}S^{--}\to 4W$, $S^{++}S^{-}\to WWWZ$;
---
2, 3 or 4 leptons, MET and jets
–
|  CMS PAS EXO-19-002 CMS:2019xud
---
$137~{}\mathrm{fb}^{-1}$
| Type-III seesaw and light scalars;
---
at least 3 charged leptons
| MadAnalysis5
---
cms_exo_19_002
|  ATLAS PRD 97 ATLAS:2018nud
---
$36.1~{}\mathrm{fb}^{-1}$
| Gauge mediated SUSY breaking;
---
(multi)photon and jets
| CheckMATE
---
atlas_1802_03158
|  ATLAS JHEP ATLAS:2021mbt
---
$139~{}\mathrm{fb}^{-1}$
| Measurement of prompt photon-pair
---
production
| Rivet/Contur
---
ATLAS_2021_I1887997
|  ATLAS EPJ C 81 ATLAS:2021fbt
---
$139~{}\mathrm{fb}^{-1}$
| RPV SUSY; many jets,
---
$\geq 1$ leptons and 0 or $\geq 3$ $b$-jets
| CheckMATE
---
atlas_2106_09609
|  ATLAS EPJ C 81 ATLAS:2021twp
---
$139~{}\mathrm{fb}^{-1}$
| Squarks and gluinos;
---
1 lepton, jets and MET
| CheckMATE
---
atlas_2101_01629
|  ATLAS EPJ C 79 ATLAS:2018zdn
---
$3.2~{}\mathrm{fb}^{-1}$
General search for new phenomena |  | CheckMATE
---
atlas_1807_07447
|  ATLAS JHEP ATLAS:2019gdh
---
$139~{}\mathrm{fb}^{-1}$
| Bottom-squark pair production;
---
no leptons, $\geq 3$ $b$-jets and MET
| CheckMATE
---
atlas_1908_03122
|  CMS PAS SUS-19-006 CMS:2019xjf
---
$137~{}\mathrm{fb}^{-1}$
| Gluinos and squarks;
---
no leptons, multiple jets and MET
| MadAnalysis5
---
cms_sus_19_006
|  CMS-SUS-16-033 CMS:2017abv
---
$35.9~{}\mathrm{fb}^{-1}$
| Gluinos and stops;
---
no leptons, multiple jets and MET
| MadAnalysis5
---
cms_sus_16_033
|  ATLAS JHEP ATLAS:2020qlk
---
$139~{}\mathrm{fb}^{-1}$
| Chargino-neutralino production;
---
MET and $h\to\gamma\gamma$
| CheckMATE
---
atlas_2004_10894
|  ATLAS JHEP ATLAS:2021kog
---
$139~{}\mathrm{fb}^{-1}$
| Measurements of four-lepton
---
differential cross sections
| Rivet/Contur
---
ATLAS_2021_I1849535
|  ATLAS JHEP ATLAS:2019gey
---
$139~{}\mathrm{fb}^{-1}$
|  Measurement of the $Z(\to\ell^{+}\ell^{-})\gamma$
---
production cross section
| Rivet/Contur
---
ATLAS_2019_I1764342
|  ATLAS JHEP ATLAS:2018nci
---
$36.1~{}\mathrm{fb}^{-1}$
|  Measurement of the $Z(\to\nu\bar{\nu})\gamma$
---
production cross section
| Rivet/Contur
---
ATLAS_2018_I1698006
|  ATLAS-CONF-2016-096 ATLAS:2016uwq
---
$13.3~{}\mathrm{fb}^{-1}$
| Electroweakino production;
---
2 to 3 leptons, MET and no jets
| CheckMATE
---
atlas_conf_2016_096
|  CMS PAS SUS-16-039 CMS:2017moi
---
$35.9~{}\mathrm{fb}^{-1}$
| Electroweakino production;
---
$\geq 2$ leptons and MET
| CheckMATE
---
cms_sus_16_039
Table 3: Summary of the analyses that contribute to the simplified model bounds in Fig. 2. Production | Channel | MadAnalysis5 | CheckMATE | Rivet/Contur
---|---|---|---|---
$S^{++}S^{--}$ | $WWWW$ | CMS:2019xud  | CMS:2017moi  |
$S^{\pm\pm}S^{\mp}$ | $WWWZ$ | CMS:2019xud  | CMS:2017moi  |
$WWW\gamma$ |  |  | ATLAS:2019gey ; ATLAS:2018nci
$S^{+}S^{-}$ | $WZWZ$ |  | CMS:2017moi  |
$WZW\gamma$ |  |  | ATLAS:2019gey ; ATLAS:2018nci
$W\gamma W\gamma$ |  | ATLAS:2018nud ; ATLAS:2020qlk  |
$S^{\pm}S^{0}$ | $WZWW$ | CMS:2019xud  | CMS:2017moi  |
$W\gamma WW$ |  |  | ATLAS:2019gey ; ATLAS:2018nci
$WZZZ$ |  | CMS:2017moi  |
$(W\gamma)(ZZ)$ |  |  | ATLAS:2019gey ; ATLAS:2021kog ; ATLAS:2018nci
$(WZ)(Z\gamma)$ |  |  | ATLAS:2019gey ; ATLAS:2018nci
$(WZ)(\gamma\gamma)$ |  | ATLAS:2018nud  |
$(W\gamma)(Z\gamma)$ |  | ATLAS:2020qlk ; ATLAS:2018nud  |
$W\gamma\gamma\gamma$ |  | ATLAS:2020qlk ; ATLAS:2018zdn ; ATLAS:2018nud  |
$S^{0}S^{\prime 0}$ | $WWWW$ | CMS:2019xud  | CMS:2017moi  |
$WWZZ$ |  | CMS:2017moi  |
$WW\gamma Z$ |  |  | ATLAS:2019gey ; ATLAS:2018nci
$WW\gamma\gamma$ |  | ATLAS:2018nud  |
$ZZZZ$ |  | CMS:2017moi  | ATLAS:2021kog
$\gamma ZZZ$ |  |  | ATLAS:2018nci ; ATLAS:2019gey ; ATLAS:2021kog
$(\gamma Z)(\gamma Z)$ |  | ATLAS:2018nud ; ATLAS:2020qlk  | ATLAS:2021kog ; ATLAS:2019gey
$(\gamma\gamma)(ZZ)$ |  | ATLAS:2018nud  |
$\gamma\gamma\gamma Z$ |  | ATLAS:2020qlk ; ATLAS:2018nud  |
$\gamma\gamma\gamma\gamma$ |  | ATLAS:2018zdn  | ATLAS:2021mbt
$S^{++}S^{--}$ | $WtbWtb$ |  | ATLAS:2021fbt  |
$S^{\pm\pm}S^{\mp}$ | $Wtbtb$ |  | ATLAS:2021fbt  |
$S^{+}S^{-}$ | $tbtb$ |  | ATLAS:2021fbt ; ATLAS:2021twp  |
$S^{\pm}S^{0}$ | $tbtt$ |  | ATLAS:2021fbt  |
$tbbb$ | CMS:2019xjf  | ATLAS:2018zdn ; ATLAS:2021fbt ; ATLAS:2021twp ; ATLAS:2019gdh  |
$S^{0}S^{\prime 0}$ | $tttt$ |  | ATLAS:2021fbt  |
$ttbb$ | CMS:2017abv  | ATLAS:2021fbt ; ATLAS:2021twp  |
$bbbb$ | CMS:2017abv  | ATLAS:2018zdn  |
Table 4: Experimental analyses contributing to the simplified model bounds in
Fig. 2. More details are available on
https://github.com/manuelkunkel/scalarbounds.

##  Appendix B $K^{SS}_{V}$ coefficients from the scalar kinetic term

As outlined in the main text, the Drell-Yan pair production process of two
gauge eigenstate scalars of an $\text{SU}(2)_{L}$ multiplet arises from a
coupling in the kinetic term of the scalar, and as such, it depends only on
the $\text{SU}(2)_{L}\times\text{U}(1)_{Y}$ quantum numbers of the scalar.

As a first example, we review the calculation for a complex scalar
$\text{SU}(2)_{L}$ triplet with hypercharge $Y$ which we denote by
$\phi_{3,Y}$. We write $\phi_{3,Y}$ as

$\displaystyle\phi_{3,Y}=\frac{1}{2}\begin{pmatrix}\phi_{3,Y}^{3}&\phi_{3,Y}^{1}-i\phi_{3,Y}^{2}\\\
\phi_{3,Y}^{1}+i\phi_{3,Y}^{2}&-\phi_{3,Y}^{3}\\\
\end{pmatrix}=\frac{1}{2}\begin{pmatrix}\phi_{3,Y}^{0}&\sqrt{2}\phi_{3,Y}^{+}\\\
\sqrt{2}\phi_{3,Y}^{-}&-\phi_{3,Y}^{0}\\\ \end{pmatrix}.$

The covariant derivative is

$\displaystyle D_{\mu}\phi_{3,Y}$
$\displaystyle=\partial_{\mu}\phi_{3,Y}-ig[W_{\mu},\phi_{3,Y}]-ig^{\prime}YB_{\mu}\phi_{3,Y}.$

The kinetic term reads

$\displaystyle\mathcal{L}_{\phi_{3,Y,\rm{kin}}}=$ $\displaystyle
2\text{Tr}(D_{\mu}\phi_{3,Y})^{\dagger}(D^{\mu}\phi_{3,Y})$ $\displaystyle=$
$\displaystyle\partial_{\mu}\phi_{3,Y}^{a*}\partial^{\mu}\phi^{a}+ieq^{a}A^{\mu}\phi_{3,Y}^{a*}\overleftrightarrow{\partial_{\mu}}\phi_{3,Y}^{a}+i\frac{e}{s_{w}c_{w}}(t_{3}^{a}-q^{a}s_{w}^{2})Z^{\mu}\phi_{3,Y}^{a*}\overleftrightarrow{\partial_{\mu}}\phi_{3,Y}^{a}$
$\displaystyle+igW^{\mu,-}\left(\phi_{3,Y}^{-*}\overleftrightarrow{\partial_{\mu}}\phi_{3,Y}^{0}-\phi_{3,Y}^{0*}\overleftrightarrow{\partial_{\mu}}\phi_{3,Y}^{+}\right)+\text{h.c.}+\mathcal{O}\left({g^{2}}\right)$

where $q^{a}=t_{3}^{a}+Y$ is the electric charge of the field $\phi_{3,Y}^{a}$
and $t_{3}^{a}$ is the corresponding eigenvalue of $T^{3}$. Here,
$a\in\\{+,0,-\\}$ indicates the $T^{3}$ quantum numbers. Note that for a
complex representation, $\phi^{+*}\neq\phi_{3,Y}^{-}$ and
$\phi_{3,Y}^{*0}\neq\phi_{3,Y}^{0}$. Comparing this with Eq. (2) and
expressing the fields via the electric charge yields the coefficients

$\displaystyle K_{W}^{\phi_{3,Y}^{1-Y}\phi_{3,Y}^{0+Y}}$ $\displaystyle=1\,,$
(24) $\displaystyle K_{W}^{\phi_{3,Y}^{0-Y}\phi_{3,Y}^{1+Y}}$
$\displaystyle=-1\,,$ (25) $\displaystyle
K_{Z}^{\phi_{3,Y}^{(t_{3}^{a}+Y)}\phi_{3,Y}^{-(t_{3}^{a}+Y)}}$
$\displaystyle=-(t_{3}^{a}-s^{2}_{w}q^{a}).$ (26)

Note that for a real representation with $Y=0$ we have $\phi^{+*}=\phi^{-}$
and $\phi^{0*}=\phi^{0}$. The kinetic term is then given by
$\mathcal{L}\phi_{3,Y}=\text{Tr}(D_{\mu}\phi_{3,Y})^{\dagger}D^{\mu}\phi_{3,Y}$,
which yields the same coefficients as Eq. 24 with $Y=0$.

As a second example, we present the $K^{SS}_{V}$ coefficients of the pNGBs
from $\text{SU}(5)/\text{SO}(5)$ breaking. They are determined analogously,
have been presented in Banerjee:2022xmu  and are listed in the following for
completeness:

$\begin{array}[]{cccccccccc}\hline\cr\hline\cr&\lx@intercol\hfil
K^{S^{0}_{i}S^{+}_{j}}_{W}\hfil\lx@intercol&K^{S^{-}_{i}S^{++}_{j}}_{W}\\\
\hline\cr&\eta^{+}_{3}&\eta^{+}_{5}&\eta_{5}^{++}\\\ \hline\cr
h&0&0&\hbox{\multirowsetup$-$}\\\
\eta^{0}_{3}&-\frac{i}{2}&\frac{c_{\theta}}{2}&\\\
\eta^{0}_{5}&-\frac{c_{\theta}}{2\sqrt{3}}&\frac{i\sqrt{3}}{2}&\\\
\eta^{0}_{1}&\sqrt{\frac{2}{3}}c_{\theta}&0&\\\ \eta&0&0&\\\
\eta^{-}_{3}&\lx@intercol\hfil\hbox{\multirowsetup$-$}\hfil\lx@intercol&\frac{c_{\theta}}{\sqrt{2}}\\\
\eta^{-}_{5}&&&-\frac{i}{\sqrt{2}}\\\ \vphantom{--}\\\
\hline\cr\hline\cr\end{array}\quad\begin{array}[]{ccccccccccccc}\hline\cr\hline\cr&\lx@intercol\hfil
K^{S^{0}_{i}S^{0}_{j}}_{Z}\hfil\lx@intercol&\lx@intercol\hfil
K^{S^{+}_{i}S^{-}_{j}}_{Z}\hfil\lx@intercol&K^{S^{++}_{i}S^{--}_{j}}_{Z}\\\
\hline\cr&h&\eta^{0}_{3}&\eta^{0}_{5}&\eta^{0}_{1}&\eta&\eta^{-}_{3}&\eta^{-}_{5}&\eta_{5}^{--}\\\
\hline\cr
h&0&0&0&0&0&\lx@intercol\hfil\hbox{\multirowsetup$-$}\hfil\lx@intercol&\hbox{\multirowsetup$-$}\\\
\eta^{0}_{3}&&0&\frac{ic_{\theta}}{\sqrt{3}}&i\sqrt{\frac{2}{3}}c_{\theta}&0\\\
\eta^{0}_{5}&&&0&0&0\\\ \eta^{0}_{1}&&&&0&0\\\ \eta&&&&&0\\\
\eta^{+}_{3}&\lx@intercol\hfil\hbox{\multirowsetup$-$}\hfil\lx@intercol&-\frac{c_{2w}}{2}&-\frac{ic_{\theta}}{2}\\\
\eta^{+}_{5}&&&&&&&-\frac{c_{2w}}{2}\\\
\eta^{++}_{5}&\lx@intercol\hfil-\hfil\lx@intercol&-c_{2w}\\\
\hline\cr\hline\cr\end{array}$

where $c_{\theta}={\rm cos}(\theta)$ with $\sin(\theta)=v/f_{\psi}$ and
$c_{2w}=\rm{cos}(2\theta_{W})$ is the cosine of twice the Weinberg angle.

## References

  * (1) F. Englert and R. Brout, “Broken Symmetry and the Mass of Gauge Vector Mesons,” Phys. Rev. Lett. 13 (1964) 321–323.
  * (2) P. W. Higgs, “Broken Symmetries and the Masses of Gauge Bosons,” Phys. Rev. Lett. 13 (1964) 508–509.
  * (3) ATLAS Collaboration, G. Aad et al., “Observation of a new particle in the search for the Standard Model Higgs boson with the ATLAS detector at the LHC,” Phys. Lett. B 716 (2012) 1–29, arXiv:1207.7214 [hep-ex].
  * (4) CMS Collaboration, S. Chatrchyan et al., “Observation of a New Boson at a Mass of 125 GeV with the CMS Experiment at the LHC,” Phys. Lett. B 716 (2012) 30–61, arXiv:1207.7235 [hep-ex].
  * (5) S. P. Martin, “A Supersymmetry primer,” Adv. Ser. Direct. High Energy Phys. 18 (1998) 1–98, arXiv:hep-ph/9709356.
  * (6) G. C. Branco, P. M. Ferreira, L. Lavoura, M. N. Rebelo, M. Sher, and J. P. Silva, “Theory and phenomenology of two-Higgs-doublet models,” Phys. Rept. 516 (2012) 1–102, arXiv:1106.0034 [hep-ph].
  * (7) J. Schechter and J. W. F. Valle, “Neutrino Masses in SU(2) x U(1) Theories,” Phys. Rev. D 22 (1980) 2227.
  * (8) H. Georgi and M. Machacek, “Doubly charged Higgs bosons,” Nucl. Phys. B 262 (1985) 463–477.
  * (9) J. Hisano and K. Tsumura, “Higgs boson mixes with an SU(2) septet representation,” Phys. Rev. D 87 (2013) 053004, arXiv:1301.6455 [hep-ph].
  * (10) A. Agugliaro, G. Cacciapaglia, A. Deandrea, and S. De Curtis, “Vacuum misalignment and pattern of scalar masses in the SU(5)/SO(5) composite Higgs model,” JHEP 02 (2019) 089, arXiv:1808.10175 [hep-ph].
  * (11) G. Ferretti, “Gauge theories of Partial Compositeness: Scenarios for Run-II of the LHC,” JHEP 06 (2016) 107, arXiv:1604.06467 [hep-ph].
  * (12) D. B. Kaplan and H. Georgi, “SU(2) x U(1) Breaking by Vacuum Misalignment,” Phys. Lett. B 136 (1984) 183–186.
  * (13) S. Weinberg, “Implications of Dynamical Symmetry Breaking,” Phys. Rev. D 13 (1976) 974–996. [Addendum: Phys.Rev.D 19, 1277–1280 (1979)].
  * (14) S. Dimopoulos and L. Susskind, “Mass Without Scalars,” Nucl. Phys. B 155 (1979) 237–252.
  * (15) K. Agashe, R. Contino, and A. Pomarol, “The Minimal composite Higgs model,” Nucl. Phys. B 719 (2005) 165–187, arXiv:hep-ph/0412089.
  * (16) R. Contino, Y. Nomura, and A. Pomarol, “Higgs as a holographic pseudoGoldstone boson,” Nucl. Phys. B 671 (2003) 148–174, arXiv:hep-ph/0306259.
  * (17) E. Witten, “Global Aspects of Current Algebra,” Nucl. Phys. B 223 (1983) 422–432.
  * (18) D. A. Kosower, “Symmetry Breaking Patterns in Pseudoreal and Real Gauge Theories,” Phys. Lett. B 144 (1984) 215–216.
  * (19) H. Georgi and D. B. Kaplan, “Composite Higgs and Custodial SU(2),” Phys. Lett. B 145 (1984) 216–220.
  * (20) K. Agashe, R. Contino, L. Da Rold, and A. Pomarol, “A Custodial symmetry for $Zb\bar{b}$,” Phys. Lett. B 641 (2006) 62–66, arXiv:hep-ph/0605341.
  * (21) T. A. Ryttov and F. Sannino, “Ultra Minimal Technicolor and its Dark Matter TIMP,” Phys. Rev. D 78 (2008) 115010, arXiv:0809.0713 [hep-ph].
  * (22) J. Galloway, J. A. Evans, M. A. Luty, and R. A. Tacchi, “Minimal Conformal Technicolor and Precision Electroweak Tests,” JHEP 10 (2010) 086, arXiv:1001.1361 [hep-ph].
  * (23) G. Cacciapaglia and F. Sannino, “Fundamental Composite (Goldstone) Higgs Dynamics,” JHEP 04 (2014) 111, arXiv:1402.0233 [hep-ph].
  * (24) I. Low, W. Skiba, and D. Tucker-Smith, “Little Higgses from an antisymmetric condensate,” Phys. Rev. D 66 (2002) 072001, arXiv:hep-ph/0207243.
  * (25) C. Cai, H.-H. Zhang, G. Cacciapaglia, M. Rosenlyst, and M. T. Frandsen, “Higgs Boson Emerging from the Dark,” Phys. Rev. Lett. 125 no. 2, (2020) 021801, arXiv:1911.12130 [hep-ph].
  * (26) M. J. Dugan, H. Georgi, and D. B. Kaplan, “Anatomy of a Composite Higgs Model,” Nucl. Phys. B 254 (1985) 299–326.
  * (27) N. Arkani-Hamed, A. G. Cohen, E. Katz, and A. E. Nelson, “The Littlest Higgs,” JHEP 07 (2002) 034, arXiv:hep-ph/0206021.
  * (28) T. Ma and G. Cacciapaglia, “Fundamental Composite 2HDM: SU(N) with 4 flavours,” JHEP 03 (2016) 211, arXiv:1508.07014 [hep-ph].
  * (29) L. Vecchi, “A dangerous irrelevant UV-completion of the composite Higgs,” JHEP 02 (2017) 094, arXiv:1506.00623 [hep-ph].
  * (30) J. Mrazek, A. Pomarol, R. Rattazzi, M. Redi, J. Serra, and A. Wulzer, “The Other Natural Two Higgs Doublet Model,” Nucl. Phys. B 853 (2011) 1–48, arXiv:1105.5403 [hep-ph].
  * (31) B. Bellazzini, C. Csáki, and J. Serra, “Composite Higgses,” Eur. Phys. J. C 74 no. 5, (2014) 2766, arXiv:1401.2457 [hep-ph].
  * (32) D. B. Kaplan, “Flavor at SSC energies: A New mechanism for dynamically generated fermion masses,” Nucl. Phys. B 365 (1991) 259–278.
  * (33) N. Bizot, G. Cacciapaglia, and T. Flacke, “Common exotic decays of top partners,” JHEP 06 (2018) 065, arXiv:1803.00021 [hep-ph].
  * (34) K.-P. Xie, G. Cacciapaglia, and T. Flacke, “Exotic decays of top partners with charge 5/3: bounds and opportunities,” JHEP 10 (2019) 134, arXiv:1907.05894 [hep-ph].
  * (35) R. Benbrik et al., “Signatures of vector-like top partners decaying into new neutral scalar or pseudoscalar bosons,” JHEP 05 (2020) 028, arXiv:1907.05929 [hep-ph].
  * (36) G. Cacciapaglia, T. Flacke, M. Park, and M. Zhang, “Exotic decays of top partners: mind the search gap,” Phys. Lett. B 798 (2019) 135015, arXiv:1908.07524 [hep-ph].
  * (37) A. Banerjee et al., “Phenomenological aspects of composite Higgs scenarios: exotic scalars and vector-like quarks,” arXiv:2203.07270 [hep-ph].
  * (38) G. Ferretti and D. Karateev, “Fermionic UV completions of Composite Higgs models,” JHEP 03 (2014) 077, arXiv:1312.5330 [hep-ph].
  * (39) G. Ferretti, “UV Completions of Partial Compositeness: The Case for a SU(4) Gauge Group,” JHEP 06 (2014) 142, arXiv:1404.7137 [hep-ph].
  * (40) A. Belyaev, G. Cacciapaglia, H. Cai, G. Ferretti, T. Flacke, A. Parolini, and H. Serodio, “Di-boson signatures as Standard Candles for Partial Compositeness,” JHEP 01 (2017) 094, arXiv:1610.06591 [hep-ph]. [Erratum: JHEP 12, 088 (2017)].
  * (41) Y. Wu, T. Ma, B. Zhang, and G. Cacciapaglia, “Composite Dark Matter and Higgs,” JHEP 11 (2017) 058, arXiv:1703.06903 [hep-ph].
  * (42) H. Cai and G. Cacciapaglia, “Singlet dark matter in the SU(6)/SO(6) composite Higgs model,” Phys. Rev. D 103 no. 5, (2021) 055002, arXiv:2007.04338 [hep-ph].
  * (43) A. Arbey, G. Cacciapaglia, H. Cai, A. Deandrea, S. Le Corre, and F. Sannino, “Fundamental Composite Electroweak Dynamics: Status at the LHC,” Phys. Rev. D 95 no. 1, (2017) 015028, arXiv:1502.04718 [hep-ph].
  * (44) G. Cacciapaglia, A. Deandrea, A. M. Iyer, and K. Sridhar, “Tera-$Z$ stage at future colliders and light composite axionlike particles,” Phys. Rev. D 105 no. 1, (2022) 015020, arXiv:2104.11064 [hep-ph].
  * (45) G. Cacciapaglia, H. Cai, A. Deandrea, T. Flacke, S. J. Lee, and A. Parolini, “Composite scalars at the LHC: the Higgs, the Sextet and the Octet,” JHEP 11 (2015) 201, arXiv:1507.02283 [hep-ph].
  * (46) G. Cacciapaglia, A. Deandrea, T. Flacke, and A. M. Iyer, “Gluon-Photon Signatures for color octet at the LHC (and beyond),” JHEP 05 (2020) 027, arXiv:2002.01474 [hep-ph].
  * (47) G. Cacciapaglia, T. Flacke, M. Kunkel, and W. Porod, “Phenomenology of unusual top partners in composite Higgs models,” JHEP 02 (2022) 208, arXiv:2112.00019 [hep-ph].
  * (48) G. Cacciapaglia, G. Ferretti, T. Flacke, and H. Serôdio, “Light scalars in composite Higgs models,” Front. in Phys. 7 (2019) 22, arXiv:1902.06890 [hep-ph].
  * (49) CMS Collaboration, “Search for resonant and nonresonant production of pairs of dijet resonances in proton-proton collisions at $\sqrt{s}$ = 13 TeV,” arXiv:2206.09997 [hep-ex].
  * (50) A. Alloul, N. D. Christensen, C. Degrande, C. Duhr, and B. Fuks, “FeynRules 2.0 - A complete toolbox for tree-level phenomenology,” Comput. Phys. Commun. 185 (2014) 2250–2300, arXiv:1310.1921 [hep-ph].
  * (51) J. Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, O. Mattelaer, H. S. Shao, T. Stelzer, P. Torrielli, and M. Zaro, “The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations,” JHEP 07 (2014) 079, arXiv:1405.0301 [hep-ph].
  * (52) R. D. Ball et al., “Parton distributions with LHC data,” Nucl. Phys. B 867 (2013) 244–289, arXiv:1207.1303 [hep-ph].
  * (53) A. Buckley, J. Ferrando, S. Lloyd, K. Nordström, B. Page, M. Rüfenacht, M. Schönherr, and G. Watt, “LHAPDF6: parton density access in the LHC precision era,” Eur. Phys. J. C 75 (2015) 132, arXiv:1412.7420 [hep-ph].
  * (54) T. Sjöstrand, S. Ask, J. R. Christiansen, R. Corke, N. Desai, P. Ilten, S. Mrenna, S. Prestel, C. O. Rasmussen, and P. Z. Skands, “An introduction to PYTHIA 8.2,” Comput. Phys. Commun. 191 (2015) 159–177, arXiv:1410.3012 [hep-ph].
  * (55) E. Conte, B. Fuks, and G. Serret, “MadAnalysis 5, A User-Friendly Framework for Collider Phenomenology,” Comput.Phys.Commun. 184 (2013) 222–256, arXiv:1206.1599 [hep-ph].
  * (56) E. Conte, B. Dumont, B. Fuks, and C. Wymant, “Designing and recasting LHC analyses with MadAnalysis 5,” Eur. Phys. J. C74 no. 10, (2014) 3103, arXiv:1405.3982 [hep-ph].
  * (57) B. Dumont, B. Fuks, S. Kraml, S. Bein, G. Chalons, et al., “Toward a public analysis database for LHC new physics searches using MADANALYSIS 5,” Eur.Phys.J. C75 no. 2, (2015) 56, arXiv:1407.3278 [hep-ph].
  * (58) E. Conte and B. Fuks, “Confronting new physics theories to LHC data with MADANALYSIS 5,” Int. J. Mod. Phys. A33 no. 28, (2018) 1830027, arXiv:1808.00480 [hep-ph].
  * (59) M. Drees, H. Dreiner, D. Schmeier, J. Tattersall, and J. S. Kim, “CheckMATE: Confronting your Favourite New Physics Model with LHC Data,” Comput. Phys. Commun. 187 (2015) 227–265, arXiv:1312.2591 [hep-ph].
  * (60) D. Dercks, N. Desai, J. S. Kim, K. Rolbiecki, J. Tattersall, and T. Weber, “CheckMATE 2: From the model to the limit,” Comput. Phys. Commun. 221 (2017) 383–418, arXiv:1611.09856 [hep-ph].
  * (61) DELPHES 3 Collaboration, J. de Favereau, C. Delaere, P. Demin, A. Giammanco, V. Lemaître, A. Mertens, and M. Selvaggi, “DELPHES 3, A modular framework for fast simulation of a generic collider experiment,” JHEP 02 (2014) 057, arXiv:1307.6346 [hep-ex].
  * (62) M. Cacciari, G. P. Salam, and G. Soyez, “The anti-$k_{t}$ jet clustering algorithm,” JHEP 04 (2008) 063, arXiv:0802.1189 [hep-ph].
  * (63) M. Cacciari, G. P. Salam, and G. Soyez, “FastJet User Manual,” Eur. Phys. J. C 72 (2012) 1896, arXiv:1111.6097 [hep-ph].
  * (64) A. L. Read, “Presentation of search results: The CL(s) technique,” J. Phys. G 28 (2002) 2693–2704.
  * (65) C. Bierlich et al., “Robust Independent Validation of Experiment and Theory: Rivet version 3,” SciPost Phys. 8 (2020) 026, arXiv:1912.05451 [hep-ph].
  * (66) J. Butterworth, “BSM constraints from model-independent measurements: A Contur Update,” J. Phys. Conf. Ser. 1271 no. 1, (2019) 012013, arXiv:1902.03067 [hep-ph].
  * (67) A. Buckley et al., “Testing new physics models with global comparisons to collider measurements: the Contur toolkit,” SciPost Phys. Core 4 (2021) 013, arXiv:2102.04377 [hep-ph].
  * (68) ATLAS Collaboration, G. Aad et al., “Search for doubly and singly charged Higgs bosons decaying into vector bosons in multi-lepton final states with the ATLAS detector using proton-proton collisions at $\sqrt{s}$ = 13 TeV,” JHEP 06 (2021) 146, arXiv:2101.11961 [hep-ex].
  * (69) CMS Collaboration, “Search for new physics in multilepton final states in pp collisions at $\sqrt{s}=13~{}\mathrm{TeV}$,”. CMS PAS EXO-19-002.
  * (70) ATLAS Collaboration, M. Aaboud et al., “Search for photonic signatures of gauge-mediated supersymmetry in 13 TeV $pp$ collisions with the ATLAS detector,” Phys. Rev. D 97 no. 9, (2018) 092006, arXiv:1802.03158 [hep-ex].
  * (71) ATLAS Collaboration, G. Aad et al., “Measurement of the production cross section of pairs of isolated photons in $pp$ collisions at 13 TeV with the ATLAS detector,” JHEP 11 (2021) 169, arXiv:2107.09330 [hep-ex].
  * (72) ATLAS Collaboration, G. Aad et al., “Search for R-parity-violating supersymmetry in a final state containing leptons and many jets with the ATLAS experiment using $\sqrt{s}=13$ TeV proton–proton collision data,” Eur. Phys. J. C 81 no. 11, (2021) 1023, arXiv:2106.09609 [hep-ex].
  * (73) ATLAS Collaboration, G. Aad et al., “Search for squarks and gluinos in final states with one isolated lepton, jets, and missing transverse momentum at $\sqrt{s}=13$ TeV with the ATLAS detector,” Eur. Phys. J. C 81 no. 7, (2021) 600, arXiv:2101.01629 [hep-ex]. [Erratum: Eur.Phys.J.C 81, 956 (2021)].
  * (74) ATLAS Collaboration, M. Aaboud et al., “A strategy for a general search for new phenomena using data-derived signal regions and its application within the ATLAS experiment,” Eur. Phys. J. C 79 no. 2, (2019) 120, arXiv:1807.07447 [hep-ex].
  * (75) ATLAS Collaboration, G. Aad et al., “Search for bottom-squark pair production with the ATLAS detector in final states containing Higgs bosons, $b$-jets and missing transverse momentum,” JHEP 12 (2019) 060, arXiv:1908.03122 [hep-ex].
  * (76) CMS Collaboration, “Search for supersymmetry in proton-proton collisions at 13 TeV in final states with jets and missing transverse momentum,”. CMS PAS SUS-19-006.
  * (77) CMS Collaboration, A. M. Sirunyan et al., “Search for supersymmetry in multijet events with missing transverse momentum in proton-proton collisions at 13 TeV,” Phys. Rev. D 96 no. 3, (2017) 032003, arXiv:1704.07781 [hep-ex].
  * (78) ATLAS Collaboration, G. Aad et al., “Search for direct production of electroweakinos in final states with missing transverse momentum and a Higgs boson decaying into photons in pp collisions at $\sqrt{s}$ = 13 TeV with the ATLAS detector,” JHEP 10 (2020) 005, arXiv:2004.10894 [hep-ex].
  * (79) ATLAS Collaboration, G. Aad et al., “Measurements of differential cross-sections in four-lepton events in 13 TeV proton-proton collisions with the ATLAS detector,” JHEP 07 (2021) 005, arXiv:2103.01918 [hep-ex].
  * (80) ATLAS Collaboration, G. Aad et al., “Measurement of the $Z(\rightarrow\ell^{+}\ell^{-})\gamma$ production cross-section in $pp$ collisions at $\sqrt{s}=13$ TeV with the ATLAS detector,” JHEP 03 (2020) 054, arXiv:1911.04813 [hep-ex].
  * (81) ATLAS Collaboration, M. Aaboud et al., “Measurement of the $Z\gamma\to\nu\overline{\nu}\gamma$ production cross section in pp collisions at $\sqrt{s}=13$ TeV with the ATLAS detector and limits on anomalous triple gauge-boson couplings,” JHEP 12 (2018) 010, arXiv:1810.04995 [hep-ex].
  * (82) ATLAS Collaboration, “Search for supersymmetry with two and three leptons and missing transverse momentum in the final state at $\sqrt{s}=13$ TeV with the ATLAS detector,”. ATLAS-CONF-2016-096.
  * (83) CMS Collaboration, A. M. Sirunyan et al., “Search for electroweak production of charginos and neutralinos in multilepton final states in proton-proton collisions at $\sqrt{s}=$ 13 TeV,” JHEP 03 (2018) 166, arXiv:1709.05406 [hep-ex].
  * (84) CMS Collaboration, A. Tumasyan et al., “Search for Higgs Boson Pair Production in the Four b Quark Final State in Proton-Proton Collisions at $\sqrt{s}$=13 TeV,” Phys. Rev. Lett. 129 no. 8, (2022) 081802, arXiv:2202.09617 [hep-ex].
  * (85) CMS Collaboration, “Search for nonresonant pair production of highly energetic Higgs bosons decaying to bottom quarks,” arXiv:2205.06667 [hep-ex].
  * (86) CMS Collaboration, “Search for Higgs boson pairs decaying to WWWW, WW$\tau\tau$, and $\tau\tau\tau\tau$ in proton-proton collisions at $\sqrt{s}$ = 13 TeV,” arXiv:2206.10268 [hep-ex].
  * (87) ATLAS Collaboration, “Search for non-resonant pair production of Higgs bosons in the $b\bar{b}b\bar{b}$ final state in $pp$ collisions at $\sqrt{s}=13$ TeV with the ATLAS detector,”. ATLAS-CONF-2022-035.
  * (88) ATLAS Collaboration, M. Aaboud et al., “Search for pair production of Higgs bosons in the $b\bar{b}b\bar{b}$ final state using proton–proton collisions at $\sqrt{s}=13$ TeV with the ATLAS detector,” Phys. Rev. D 94 no. 5, (2016) 052002, arXiv:1606.04782 [hep-ex].
  * (89) ATLAS Collaboration, M. Aaboud et al., “Search for Higgs boson pair production in the $\gamma\gamma WW^{*}$ channel using $pp$ collision data recorded at $\sqrt{s}=13$ TeV with the ATLAS detector,” Eur. Phys. J. C 78 no. 12, (2018) 1007, arXiv:1807.08567 [hep-ex].
  * (90) ATLAS Collaboration, M. Aaboud et al., “Search for Higgs boson pair production in the $WW^{(*)}WW^{(*)}$ decay channel using ATLAS data recorded at $\sqrt{s}=13$ TeV,” JHEP 05 (2019) 124, arXiv:1811.11028 [hep-ex].
  * (91) CMS Collaboration, A. Tumasyan et al., “Search for new particles in an extended Higgs sector with four b quarks in the final state at $\sqrt{s}$ = 13 TeV,” arXiv:2203.00480 [hep-ex].
  * (92) B. Fuks, M. Klasen, D. R. Lamprea, and M. Rothering, “Precision predictions for electroweak superpartner production at hadron colliders with Resummino,” Eur. Phys. J. C 73 (2013) 2480, arXiv:1304.0790 [hep-ph].
  * (93) CMS Collaboration, A. M. Sirunyan et al., “Search for a charged Higgs boson decaying into top and bottom quarks in events with electrons or muons in proton-proton collisions at $\sqrt{\mathrm{s}}$ = 13 TeV,” JHEP 01 (2020) 096, arXiv:1908.09206 [hep-ex].
  * (94) ATLAS Collaboration, G. Aad et al., “Search for charged Higgs bosons decaying into a top quark and a bottom quark at $\sqrt{\mathrm{s}}$ = 13 TeV with the ATLAS detector,” JHEP 06 (2021) 145, arXiv:2102.10076 [hep-ex].
  * (95) ATLAS Collaboration, “Search for $t\bar{t}H/A\rightarrow t\bar{t}t\bar{t}$ production in the multilepton final state in proton-proton collisions at $\sqrt{s}=13$ TeV with the ATLAS detector,”. ATLAS-CONF-2022-008.
  * (96) ATLAS Collaboration, M. Aaboud et al., “Search for heavy particles decaying into top-quark pairs using lepton-plus-jets events in proton–proton collisions at $\sqrt{s}=13$ TeV with the ATLAS detector,” Eur. Phys. J. C 78 no. 7, (2018) 565, arXiv:1804.10823 [hep-ex].
  * (97) CMS Collaboration, A. M. Sirunyan et al., “Search for resonant $\mathrm{t}\overline{\mathrm{t}}$ production in proton-proton collisions at $\sqrt{s}=13$ TeV,” JHEP 04 (2019) 031, arXiv:1810.05905 [hep-ex].
  * (98) ATLAS Collaboration, G. Aad et al., “Search for resonances decaying into photon pairs in 139 fb-1 of $pp$ collisions at $\sqrt{s}$=13 TeV with the ATLAS detector,” Phys. Lett. B 822 (2021) 136651, arXiv:2102.13405 [hep-ex].
  * (99) S. Dittmaier, M. Kramer, M. Spira, and M. Walser, “Charged-Higgs-boson production at the LHC: NLO supersymmetric QCD corrections,” Phys. Rev. D 83 (2011) 055005, arXiv:0906.2648 [hep-ph].
  * (100) LHC Higgs Cross Section Working Group Collaboration. https://twiki.cern.ch/twiki/bin/view/LHCPhysics/CERNYellowReportPageAt13TeV.
  * (101) LHC Higgs Cross Section Working Group Collaboration, D. de Florian et al., “Handbook of LHC Higgs Cross Sections: 4. Deciphering the Nature of the Higgs Sector,” arXiv:1610.07922 [hep-ph].
  * (102) J. Y. Araz, A. Buckley, B. Fuks, H. Reyes-Gonzalez, W. Waltenberger, S. L. Williamson, and J. Yellen, “Strength in numbers: optimal and scalable combination of LHC new-physics searches,” arXiv:2209.00025 [hep-ph].
  * (103) I. Hinchliffe, A. Kotwal, M. L. Mangano, C. Quigg, and L.-T. Wang, “Luminosity goals for a 100-TeV pp collider,” Int. J. Mod. Phys. A 30 no. 23, (2015) 1544002, arXiv:1504.06108 [hep-ph].
  * (104) N. Arkani-Hamed, T. Han, M. Mangano, and L.-T. Wang, “Physics opportunities of a 100 TeV proton–proton collider,” Phys. Rept. 652 (2016) 1–49, arXiv:1511.06495 [hep-ph].
  * (105) J. P. Delahaye, M. Diemoz, K. Long, B. Mansoulié, N. Pastrone, L. Rivkin, D. Schulte, A. Skrinsky, and A. Wulzer, “Muon Colliders,” arXiv:1901.06150 [physics.acc-ph].
  * (106) T. Han, S. Li, S. Su, W. Su, and Y. Wu, “BSM Higgs Production at a Muon Collider,” in 2022 Snowmass Summer Study.  5, 2022.  arXiv:2205.11730 [hep-ph].