Relaxion and light (pseudo)scalars at the HL-LHC and lepton colliders
Claudia Frugiuele Elina Fuchs Gilad Perez Matthias Schlaffer




Deep nested level sets: Fully automated segmentation of cardiac MR images in patients with pulmonary hypertension
    Jinming Duan^1,2, Jo Schlemper^1, Wenjia Bai^1, Timothy J W Dawes^2, Ghalib Bello^2, Georgia Doumou^2, Antonio De Marvao^2, Declan P O'Regan^2, Daniel Rueckert^1
    December 30, 2023
=====================================================================================================================================================================





§ INTRODUCTION


The discovery of the Higgs boson in 2012 was an extraordinary success for both the Standard Model
(SM) and the Large Hadron Collider (LHC).  Six years later, the LHC is forcing us to reconsider our
views regarding the existence of New Physics (NP) at or near the electroweak scale.  Despite a rich
program of dedicated searches in many channels, no conclusive signs of NP between the electroweak
and the TeV scale have been found.  Consequently, the idea is questioned that NP solving the
hierarchy problem is at the few-TeV scale and traditional symmetry-based approaches are being
challenged.

Generally speaking, constraints on new particles coupled to the SM via either the gauge interactions
or other couplings of similar strength are getting increasingly strong, and push models with such
couplings to non-generic corners of parameter space.  As a consequence, a growing attention is now
given to an alternative point of view where the additional degrees of freedom are light and have
ultra-weak couplings to the SM, and hence are hard to find.  Searching for NP of this kind can be
motivated based on purely phenomenological reasonings, but it is also motivated by generic
theoretical concepts that support the existence of new weakly coupled light fields.  Such moduli or
pseudo-Nambu-Goldstone-Bosons (pNGBs) appear in many extensions of the SM, see
e.g. Ref. <cit.> for a recent review and references therein.  The possibility for the
presence of these new light fields is further supported by theoretical frameworks that address other
conceptual and observational open questions such as the hierarchy and the strong  problems or
the origin of dark matter (see e.g. Refs. <cit.>).

Recently, a new mechanism has been proposed that addresses the hierarchy problem in a way that goes
beyond the conventional paradigm of a symmetry-based solution to fine-tuning.  It belongs to the
abovementioned class of models where the solution is associated with the existence of a new and
special kind of pNGB, and it is worth to be examined in detail.  It is denoted as the relaxion
mechanism <cit.>, where the pNGB—the relaxion—stabilizes the Higgs mass dynamically. The
Higgs mass depends on the classical relaxion field value which evolves in time. The relaxion rolls
down a potential, eventually stopping its rolling at a special field value where the Higgs mass is
much smaller than the theory's cutoff, hence addressing the fine-tuning problem. Since the rolling
needs to be very slow, a frictional force is necessary.  The main sources of friction that were
proposed thus far are either the Hubble friction during inflation, thereby linking the solution of
the hierarchy problem to cosmology <cit.>, or alternatively through particle production, see
Ref. <cit.>. Relaxion models do not require TeV-scale top, gauge or Higgs
partners. However, they do generically lead to an interesting phenomenology. In addition to
cosmological signatures <cit.>, relaxion models can also leave fingerprints at the
low energy precision frontier, the intensity frontier and colliders <cit.>.

As already discussed in Refs. <cit.>, the relaxion parameter space spans many
orders of magnitude in mass and coupling, see Fig. <ref>.  This is somewhat similar to
the case of the axion, but also different in ways discussed below.  The relaxion can be as light as
10^-20 eV (e.g. for Hubble-friction based models <cit.>, and ignoring the
relaxion's quality problem <cit.>) up to tens of
GeV <cit.>. This is a huge parameter space to probe and it cannot be
scanned by a single experiment, or a single frontier of new physics searches. Even more so, most
experiments do not even touch the physical parameter space of relaxion models. It is left to
consider which of the High-Luminosity (HL)-LHC and future lepton colliders have the power to probe
relaxion dynamics.

In this paper we mostly focus on searching for GeV-scale relaxions, which are ultra-heavy with
respect to the relaxion framework, but at the same time at the low-mass end of the collider reach
and therefore pose an experimental challenge.  Admittedly, this region covers only a small part of
the parameter space of the framework, however, as we shall see, future experiments will have the
power to significantly probe physical regions of it. Furthermore, this region, which addresses the
little hierarchy problem <cit.>, is typically free of some of the model's theoretical flaws
such as an exponential number of e-folds <cit.>, the need for significant
clockworking <cit.>, and finally a quality
problem <cit.>.


As for the present, low-mass resonances are already being searched for at the LHC, but the further
development of reconstruction and analysis techniques could boost its power in this regime.  Indeed,
an increased effort has been invested to show that the LHC is powerful in probing several
interesting BSM scenarios featuring light particles <cit.>.

We focus in particular on an interesting property of the relaxion, namely that it mixes with the
Higgs due to the fact that  is spontaneously broken by its vacuum expectation value (vev), see
e.g. Refs. <cit.>.  The presence of a light new scalar mixed
with the SM-like Higgs is not only a feature common to several new approaches to solve the hierarchy
problem <cit.>, but it is ubiquitous in many other more
traditional scenarios for physics beyond the SM (BSM) such as the singlet
extension <cit.>, the Higgs portal <cit.>, 2-Higgs
doublet models (2HDM) <cit.>, supersymmetry
(SUSY) <cit.> and several others.

As a theoretically motivated framework with a light pNGB with both -even and -odd couplings,
the relaxion provides us with a useful benchmark model to examine the capabilities of the HL-LHC and
future lepton colliders to probe the corresponding parameter space in direct and indirect channels.
While our discussion mainly focuses on the relaxion framework, most of the resulting bounds also
apply to the singlet extension, Higgs portal and other models.

Our paper is structured as follows. In Sect. <ref> we summarize the relaxion
mechanism and highlight its couplings that are relevant for collider phenomenology. In addition, we
collect in Sect. <ref> the existing bounds on the model's rather wide parameter
space.  In Sect. <ref> we discuss in detail direct and indirect probes of
-even and -odd couplings, both at hadron and lepton colliders.  Finally we conclude in
Sect. <ref>.




§ RELAXION PHENOMENOLOGY




In the following we will briefly review the relaxion mechanism. The effective scalar potential of
the theory depends both on the Higgs doublet H and the relaxion ϕ,

    V(H,ϕ)   =μ^2(ϕ) H^† H+λ (H^ † H)^2 + V_ sr(ϕ)+ V_br(h,ϕ)  ,
    μ^2(ϕ)   =-Λ^2+g Λϕ + …  ,

where h is the physical component of the Higgs doublet (before mixing with ϕ, see below), and
Λ is the cutoff scale of a Higgs loop and. As discussed in Ref. <cit.>, Λ
needs to be significantly smaller than the Planck mass. Hence at energies above Λ,
additional structure is required to fully stabilize the Higgs mass. See Refs. <cit.> for attempts in UV completing relaxion models by supersymmetry and
composite Higgs.  During its evolution, the relaxion scans the Higgs mass parameter μ^2(ϕ)
from a large and positive cutoff energy Λ^2 down to negative values because of the slow-roll
potential

    V_ sr(ϕ)=r gΛ^3 ϕ ,

where g is a (small) dimensionless coupling which breaks the shift symmetry of the backreaction
potential V_br and r > 1/16 π^2 due to naturalness requirements. Once
μ^2(ϕ) becomes negative, the Higgs gets a vev v^2(ϕ)=-μ^2(ϕ)/λ. This
non-zero vev activates a backreaction potential V_br which eventually stops the rolling of
the relaxion at a value ϕ_0, where v(ϕ_0)=246.

The properties of the backreaction mechanism are model-dependent. In the minimal relaxion model
discussed in Ref. <cit.>, V_ br is generated by low-energy QCD and thus the relaxion
is identified with the QCD axion. Yet, this setup typically predicts a too large phase
θ_QCD and is therefore ruled out by the upper bound on the neutron electric dipole
moment. Possible ways to suppress the  violation associated with the relaxion mechanism are
discussed in Refs. <cit.>.  Alternatively one can introduce a new sector with
strongly <cit.> or weakly <cit.> coupled vector-like fermions, which generates a
Higgs-dependent backreaction potential of the form

    V_br( h,ϕ)= -M̃^4-j(v(ϕ)+h/√(2))^j cos( ϕ/f),
 
  with 1≤ j≤ 4 and M̃ being a mass parameter <cit.>. For later
  convenience we define the backreaction scale
  Λ_br(v(ϕ))^4 ≡M̃^4-j v(ϕ)^j/√(2)^j.  This scale is not
  predicted. It can be close to (or larger than) than the electroweak scale or as small as sub-GeV
  as in the QCD axion case <cit.>, which corresponds to j=1.  In this paper we consider
  only the case of j=2.  In order to avoid fine-tuning, the backreaction scale is limited
  by <cit.>
 
    Λ_ br(v(ϕ_0))^2 ≲𝒪(2 π v^2)
      .


In the following we will focus on the class of relaxion models where the backreaction has the form
of Eq. (<ref>) and investigate its phenomenological consequences. Other realizations of the
relaxion mechanism have been introduced for instance in
Refs. <cit.>, where our analysis does not apply.
  
 



 §.§ Relevant parameters



In Refs. <cit.> the couplings of the relaxion to the SM particles and
their phenomenological consequences were studied.  It was pointed out that generically the
backreaction and/or the relaxion vev <cit.> lead to  violation. As a result the
relaxion mixes with the Higgs and inherits its couplings to SM fields.

It is convenient to write

    Λ_ br = r_ br v ,

where r_ br≲√(2π) from Eq. (<ref>).

From the diagonalization of the 2× 2 mass matrix (see Ref. <cit.> for the full
expressions) and imposing the largest eigenvalue to be the Higgs mass of m_h=125 we obtain

    m^2_ϕ   = ^4 v^4 [f^2 (c_0 m_h^2 - 4^4 v^2 s_0^2) - c_0^2 ^4 v^4 ]/f^2 (m_h^2 f^2 -c_0 ^4 v^4)
       ≃ r_ br^4 v^4 [f^2 (c_0 - 16  r_ br^4 s_0^2 ) - 
        4 c_0^2  r_ br^4 v^2] /f^2 (f^2 - 4  c_0  r_ br^4 v^2) ,

where the approximation holds for m_h≈ v/2. Here and in the following we denote
s_0≡sin(ϕ_0/f), c_0≡cos(ϕ_0/f), where ϕ_0 is the endpoint of the rolling
of the relaxion. Moreover, to keep the notation simple, we also denote the mixed mass eigenstates by
h and ϕ. In the limit of f ≫ r^2_ br v, the expression for the relaxion mass is
simplified to

    m_ϕ≃ r_ br^2 v^2 /f√(c_0 -16 r_ br^4 s_0^2 ) .

This limit requires 16^4 s_0^2 < c_0 which demands either a small backreaction scale or a
suppressed -violating angle s_0.  This condition is easily met since for
Λ_ br > v the end-point of the relaxion rolling has a suppressed s_0 as clarified in
Ref. <cit.>.  The mixing with the Higgs is given by

    sinθ= | 8 f r_ br^4  s_0 v/√(
    64 f^2 r_ br^8 s_0^2 v^2 + (f^2 - 4  c_0 r_ br^4 v^2)^2)|
     ≃| 8   r_ br^4 s_0 v/ f  | ,

where the approximation holds in the limit of f ≫ r^2_ br v. For later use we define
s_θ≡sinθ and c_θ≡cosθ. When 16 ^4 s_0^2 ≪ c_0,
Eq. (<ref>) further simplifies to

    m_ϕ≃ r_ br^2 v^2 /f√(c_0 ) .

Moreover, using Eqs. (<ref>) and (<ref>) we find that the mixing is bounded
by

    sinθ≤ 2/v.

This bound holds within the above approximations, which are fulfilled in most of the parameter
space.




 §.§ The relaxion-Higgs sector




As discussed in Ref. <cit.>, via its mixing with the Higgs, the relaxion inherits the
Higgs couplings to SM particles g_hX, where X=f,V suppressed by a universal factor
i.e. the mixing angle sinθ. This is precisely the case as for Higgs-portal models.
The couplings of the relaxion ϕ to fermions f and vector bosons V are thus given by

    g_ϕ X = sinθ g_h X .

The upper bound on the Higgs-relaxion mixing is of the same strength as in general Higgs portal
models (see e.g. Refs. <cit.>), where the mixing of ϕ with
the Higgs as a function of the mass  is bounded by naturalness as
sinθ≤ 2/v, cf. Eq. (<ref>).


Besides the couplings to SM fermions and gauge bosons, the relaxion also couples to the Higgs boson
via c_hϕϕ  hϕϕ which is given by (see Ref. <cit.>, but here with a
general s_0, c_0)

    c_ϕϕ h = ^4 v^3/f^2 c_0 ^3
    - 2 ^4  v^2/f s_0 ^2 
    - ^4 v^4/2f^3 s_0 ^2 
    - 2^4 v^3/f^2 c_0 ^2
    + 3v λ^2
    + ^4 v^2/f s_0 ^3 .

Note that only two of the parameters f, r_br and θ are independent. In the limit
of small relaxion-Higgs mixing, i.e. sinθ→0, the coupling c_ϕϕ h does
not vanish, but is reduced to the first term with →1, originating from the backreaction
potential.  Furthermore, the comparison with the expression for  in Eq. (<ref>)
allows for c_0≫ 16^4 s_0^2 to express the coupling in terms of the mass as

    c_ϕϕ h|_θ→0≃^4 v^3/f^2 c_0 ^3
      ≃^2/v .

Hence the coupling becomes independent of θ. This observation is reflected in the limits on
(, sinθ) derived from bounds on the c_ϕϕ h coupling from
h→ϕϕ decays, see Sects. <ref> and <ref> as well as
Figs. <ref> and <ref>.

Furthermore, the mixing with the relaxion modifies the Higgs self-coupling λ with respect to
its SM value of λ_ SM=m_h^2/2 v^2 (at tree level). By demanding the heavier
mass eigenvalue to correctly reproduce the observed Higgs mass, λ can be expressed as a
function of f,  and s_0 as

    λ   = -f^2 m_h^4 + c_0  m_h^2 ^4 v^4 + 4^8  s_0^2  v^6/-2f^2 m_h^2 v^2 + 2c_0 ^4 v^6≃f^2 - 4^4 (c_0+16^4 s_0^2) v^2/8(f^2 - 4c_0 ^4 v^2) ,

where the simplification holds for the approximation of m_h≃ v/2. In addition, in order to
avoid a negative λ and other overly large contributions to Higgs-couplings, we shall take
f>1 in our numerical analysis.




 §.§ Relaxion–gauge-boson interactions



An important ingredient for bounds on the relaxion parameter space are the couplings to the SM gauge
bosons.  The relaxion couples to the gauge bosons via two different classes of couplings, one that
arises through its mixing with the Higgs and one that can arise due to its pNGB nature and is
generically dictated by symmetry principles. The former couplings are found when  violation
(CPV) is present, as discussed above, while the latter are  conserving (CPC).

The leading terms in the effective Lagrangian that describes the -even relaxion–gauge-boson
interactions are

    ℒ⊃ -ϕ sinθ(   α/16 π v  F^μνF_μν + α/8π v tanθ_ W  Z^μνF_μν
               + α_s/16  π v  G^μν,aG^a_μν
       + 2m_W^2/v  W^+ μW^-_μ 
               +m_Z^2/v Z^μZ_μ)  ,

where an order one factor that depends on m_ϕ is omitted for simplicity (see
e.g. <cit.>).  If some of the SM fields couple directly to the relaxion in a
manner that breaks the shift symmetry, then the above couplings will also contain non-universal
pieces that cannot be described only by sinθ, a fact that leads to model dependence.

As for the pseudoscalar couplings of the relaxion to the gauge fields, these arise if the
backreaction sector is anomalously charged under the corresponding gauge group, we thus write the
effective Lagrangian as

    ℒ⊃ϕ/4π f(
      c̃_γγ/4
      F_μνF^μν
      + c̃_Zγ/2Z_μνF^μν + c̃_ZZ/4
      Z_μνZ^μν + c̃_WW/4 W_μνW^μν +
      c̃_GG/4G_μνG^μν) .

In order to compare the relevant importance of the pseudoscalar couplings to the ones induced by the
mixing with the Higgs let us consider for instance the couplings to photons.  For example, in cases
where the coupling to photons is induced via a weakly coupled set of fermion fields, we expect
c̃_γγ∼α Q_ eff^2, where Q_ eff stands for the chiral sum of
the fermion charges that induces the anomaly.  We note that if the backreaction sector is
anomaly-free, e.g. consists of a SM gauge-neutral sector (see for example the models
described in Refs. <cit.>), then Q_ eff is suppressed by at
least one additional loop factor and is expected to be rather small <cit.>.

Two cases can be distinguished.


  * There are two relevant couplings that are unique to the CPV sector: the triple scalar
  coupling, c_ϕϕ h hϕϕ∼^2/v, given in Eq. (<ref>); and
  the relaxion coupling to WW (ZZ), defined in Eq. (<ref>) as
  c_ϕ WW (ZZ)∼ s_θ 2 m_W^2/v (s_θ 2 m_Z^2/v).  These couplings
  generically dominate the CPC ones at colliders, in particular when one can produce at least one
  electroweak/Higgs boson on shell. Thus play a crucial role below.

  * As for the coupling to photons, in this case both the CPV and CPC couplings are irrelevant in our effective description, and their ratio is given by 

    sinθ×f/v×1/Q_eff^2∼r_br^2/Q_eff^2 .

 We find that the CPV coupling would dominate only for a relatively large backreaction scale, or in cases where the backreaction sector is non-anomalous.
 
 


 



§ STATUS OF EXPERIMENTAL PROBES FOR THE RELAXION



As just summarized in Sect. <ref>, the relaxion mixes with the SM-like Higgs, leading to a
scalar interaction between the relaxion ϕ and the SM particles. This provides the possibility
to search for such a particle in various experimental setups <cit.>,
depending on its mass.  Its possible mass range spans several orders of magnitude from sub-eV to
several tens of GeV. Hence, it can be directly produced at high-energy colliders or be indirectly
detected through its long-range, spin-independent interaction between matter constituents.



 §.§ CP-even couplings


Refs. <cit.> studied the probes of the relaxion parameter space in
detail. We summarize their results in the overview plot in Fig. <ref>.  The diagonal
line represents the maximal mixing as allowed by Eq. (<ref>).  One can classify the
following types of probes that are relevant for the parameter space below or near this line:


  * Fifth force experiments for m_ϕ <10 eV. They comprise bounds from the Casimir
  effect <cit.> (shown here) and the Eöt-Wash-type
  experiments <cit.> for the equivalence
  principle <cit.>
  and the inverse square law <cit.> (for lower masses than
  shown here, see Ref. <cit.>).

  * Astrophysical bounds for ≲ 300. They stem from red giants,
  horizontal branch stars <cit.> and the
  Supernova 1987A <cit.>.

  * Beam dump experiments for 1<≲ 300. In particular proton fixed
  target experiments are sensitive. The presented bound is from
  CHARM <cit.>.

  * Rare meson decays for ≲ 5. Relaxions can mediate rare decays of
  K,B and Υ mesons, hence their branching ratios constrain the relaxion-Higgs mixing
  angle, see e.g. Refs. <cit.> and
  references therein.

  * LEP for 0.3<<116. Searches for Higgs-like particles produced in the
  Z decay Z→ ff̅ϕ and via Higgs-strahlung e^+ e^- →ϕ Z have been
  performed <cit.>.

  * Higgs decays for <m_h/2. The LHC Run-1 set an upper limit on the untagged
  branching ratio of the Higgs boson <cit.>.

 
In addition, the relaxion is constrained by cosmological bounds on its late decays.  Significant
regions of the parameter space can be probed by big bang nucleosynthesis (BBN), the cosmic microwave
background (CMB) and extragalactic background light (EBL) as discussed in
Refs. <cit.>.  We do not present these bounds here because they
depend on various details like the reheating temperature, resulting in a strong model dependence
that may even allow to circumvent them <cit.>.  Instead, in Fig. <ref> we
show only lifetime contours to highlight where cosmological probes could play an important role, in
particular 1 s as an indication of the beginning of BBN, 10^17 s as the lifetime of the
universe, and 10^26 s to indicate where the relaxion is unbounded by EBL <cit.>.

In the 10 eV–1 MeV region, only astrophysical probes are sensitive enough to probe relevant
regions, i.e. below the line from Eq. (<ref>).  Together with the bounds
from invisible kaon decays shown in Fig. <ref>, neutron scattering experiments and
atomic precision measurements can in principle probe Higgs portal models, but they are not yet
sensitive enough to set competitive bounds <cit.>. Above 1 MeV, rare meson decays, beam dump and collider
experiments have sensitivity to the physical parameter region below the line.

Future experiments are being planned to improve on these existing limits.  For instance both
SHiP <cit.> and the beam dump run of NA62 <cit.> will extend the
coverage in the MeV–GeV mass range.  Moreover, the MATHUSLA surface detector will be able to
constrain the couplings of light, long-lived scalars below ≃ 5 <cit.>. In the following we will analyse how a relaxion in the 5 GeV–35 GeV mass range
can be probed by the HL-LHC and future lepton colliders.  Heavier relaxions with a mass of
35 GeV–62.5 GeV are already ruled out by the present limits on untagged Higgs decays as shown in
Ref. <cit.>.  Other scalars that mix with the Higgs could still exist in this mass
range.




 §.§ CP-odd couplings



Pseudoscalar couplings of the relaxion originate from the backreaction sector, and thus they depend
on the details of the specific model as we briefly discussed in Sect. <ref>.  The
perspectives to probe the relaxion via these -odd couplings are subject to the relative size of
the -odd and -even couplings.
In the following we focus on the mass region above 5 GeV and we refer to
Refs. <cit.> for a detailed discussion of the
phenomenology of axion-like particles with sinθ =0.  As it was shown in
Ref. <cit.>, heavy-ion collisions at the LHC can provide the best limits on -odd
photon couplings in the 5 GeV–100 GeV mass range. Considering Pb-Pb collisions with a luminosity
of 1 nb^-1 yields

    f/c̃_γγ≳ 500  .
 
LEP mono-photon searches provide a slightly weaker bound,
f/c̃_γγ≳ 300 <cit.>.

The ϕγ Z dual coupling is instead constrained by rare Z decays. For Higgs-like decays of
the relaxion, the strongest limit of

    f/c̃_Z γ≳ 1    TeV

comes from ALEPH <cit.> as we derive in
Sect. <ref>. This decay pattern does not only appear when the Higgs-relaxion
mixing dominates, but also in scenarios with -odd couplings to SM fermions as in
Ref. <cit.>. The case of  was studied in
Ref. <cit.>.
Similarly, the analogous ϕ Z Z coupling gives rise to the rare Z decay
Z →ϕ Z^*(f f̅), and to a weak bound via the additional contribution to the total
Z-width, Γ_Z^NP <cit.>.

Furthermore, since the relaxion has both -odd and -even couplings, they possibly generate a
contribution to electric dipole moments (EDMs). The strongest bound comes from the electron EDM,
d_e. In particular the present upper bound of
d_e/e ∼ 8 × 10^-29 cm <cit.> implies for the product of the
-odd coupling to photons and -even one to electrons <cit.>,

    αc̃_γγ/πm_e/vsinθ/f≲  2
      ×  10^-13^-1 ⇒ f/c̃_γγ >  2.5× 10^4 sinθ ,
 
where we omitted the logarithmic mass dependence (and, for simplicity, we omit the model dependent
-odd couplings to fermions here).  In the coming years improvements of one order of magnitude
are expected <cit.>.  Consequently, as we will see in the following, a relaxion signal
could be expected at future EDM experiments and high energy colliders, such as the HL-LHC or future
electron colliders.



§ PROSPECTS FOR  COLLIDERS



After reviewing the constraints on the relaxion mass and mixing angle in the previous section, we
will now identify and discuss different sensitive channels at hadron and electron colliders that are
able to improve on the existing limits. The region of relaxion masses between 5 and 35
is so far only poorly constrained, but can be probed in various ways by runs with higher luminosity
at the LHC as well as at future colliders. Higgs searches at LEP and constraints from untagged Higgs
decays at the LHC set the so-far strongest limits in this region.

Before going to the details of the collider channels, we notice that both very low f≪ M_ Pl
and ≲𝒪(1) are required in order to obtain a relaxion mass in the GeV
range <cit.>, which implies that only the little hierarchy can be
addressed by the relaxion.  For >5, Eq. (<ref>) implies the following relations
between the underlying model parameters:

    f ≥{ 1, 10}   ≳{  0.3, 0.9  }
    
      s_0 ≥{ 1/√(2) , 0.1, 0.01 }    f ≲{ 1.5, 15, 150}
    ≤{ 1, 1.5 }    f ≲{ 10, 30} .

With f and Λ_ br≡ r_ br v potentially within the range of future colliders,
there is a chance to encounter the richer phenomenology of the UV completion of a relaxion model,
beyond a discovery of the relaxion itself.  Moreover, the heavier the relaxion, the looser the upper
bound on the mixing angle sinθ, see Eq. (<ref>).

In the following we will use the relaxion as a benchmark model. However, our results are more
universal and apply also to a general Higgs portal.  The only difference is the model-dependence of
the quartic Higgs coupling λ given in Eq. (<ref>) and the triple-scalar coupling
c_hϕϕ given in Eq. (<ref>).  Hence our interpretation of the limit on the
h→ϕϕ decay is relaxion-specific, while all the other presented bounds are general.




 §.§ Precision probes


In this section we study possible constraints on the relaxion parameter space via precision
measurements of Higgs and Z properties, both at the HL-LHC and future lepton colliders.



  §.§.§ Untagged Higgs decays



A powerful way to probe the relaxion at colliders is the exotic Higgs decay into a pair of
relaxions, h→ϕϕ <cit.>, exploiting the sensitivity to the triple coupling
c_ϕϕ h given in Eq. (<ref>). There are two complementary ways of deriving
bounds from this decay. One possibility are direct searches for the decay products of the relaxions;
the other one is indirect by virtue of the NP contribution to the Higgs width in a global fit of
Higgs couplings.  In this Section we focus on the latter constraint, whereas we review the various
direct searches in Sect. <ref>.  In Fig. <ref> both approaches are
compared, showing their strong potential in particular at lepton colliders.

Global Higgs coupling fits are sensitive to additional exotic decay channels for the Higgs that
remain untagged in the corresponding searches.  The measured Higgs rates allow for a global fit of
the Higgs coupling modifiers κ_i and the branching ratio into NP (h→) as an
additional parameter.  Because the measured Higgs rates consist of the product of production and
decay, the fit is model dependent.  For the case of the relaxion mixed with the Higgs, a
two-parameter fit is required, namely a universal modifier of the Higgs couplings to SM particles
that can be identified as κ≡cosθ (thus automatically κ_V≤ 1, V=W,Z),
and (h→) that is realized by the h→ϕϕ decay channel.  The total Higgs width is
given by

    Γ_h = κ^2 Γ_h^ SM + Γ_h^ .

In general, the NP contribution to the Higgs width consists of

    Γ_h^ = Γ_h^ + Γ_h^ ,

where Γ_h^ denotes the partial width into invisible particles and
Γ_h^ denotes the partial width into untagged final states that are not
necessarily undetectable, but were not accounted for in the data set included in the fit, see
e.g. Refs. <cit.>.  In the relaxion case, we are interested in constraining
Γ_h^=Γ_h^=Γ(h→ϕϕ), i.e. not the invisible width.  With
masses in the GeV range, the relaxion is short-lived and decays inside the detector even for small
sinθ <cit.>. For the case of only one universal coupling modifier κ,
Ref. <cit.> obtains BR(h→ NP)≤ 20% from Run-1 data of the LHC
and Ref. <cit.> reports a similar bound. This bound applies directly to the
relaxion setup, as long as the modifications of the Higgs phenomenology are sufficiently described
by sinθ and the impact of the backreaction sector remains small.

Current and future runs at the LHC as well as Higgs precision measurements at future
lepton colliders will tighten the bound on untagged Higgs decays.  In Tab. <ref>
we collect the projections and estimates of upper bounds on BR(h→ untagged) at the
LHC, ILC, CEPC, CLIC and FCCee running at different energies and luminosities, and
compare them to the bounds on BR(h→ invisible) from the literature.







  
Estimate of (h→)
While projections of the sensitivity to invisible Higgs decays have been published for all
considered lepton colliders <cit.>, to our knowledge such systematic
studies are still missing for Higgs decays to untagged final states. Therefore we estimate the upper
bounds on (h→) in the following way.

Due to the NP contribution to the total Higgs width in Eq. (<ref>), each branching ratio
into SM final states F is diluted as

    (h→ F) = ^(h→ F) · [1-(h→)] .


A bound can be set on (h→) via the precision δ_κ of the experimental
determination of κ <cit.>[We thank Tim Stefaniak for a helpful
  discussion regarding the approximate bounds on untagged Higgs decays at lepton
  colliders.]. Assuming that the observed signal rates of production times decay are SM-like within
n times the uncertainty, n·δ_κ, the modification of the production and decay must
result in the product as

    (1 - n·δ_κ)^2 ≤κ^2 · [1-(h→)] .

Hence, the branching ratio into NP can be bounded as

    (h→) ≤ 1 - (1-n·δ_κ/κ)^2 .


Coupling fits for a model with one universal modifier have not been performed for all of the
colliders considered here. Instead, several coupling modifiers for different particle species have
usually been included. However, the by far most precisely determined coupling at lepton colliders is
generally κ_Z, i.e. the coupling of hZZ due to the absolute measurement of the
Zh associated production cross section. Therefore this parameter dominates in the fit and is
expected to yield a similar result as in a one-parameter fit. A fit of only one κ and
(h→) would result in a stronger bound; thus the described approximation can be regarded as
a conservative estimate.

In the relaxion case, (h→)≡(h→ϕϕ) depends on  and
κ≡√(1-sin^2θ). Thereby, as shown in Fig. <ref>, we set
constraints on (,sin^2θ) approximately at the 95% CL. We obtain these by setting
n=2 in Eq. (<ref>), and using the projected precision of κ_Z from the
various lepton colliders: ^ CLIC = 0.8% <cit.>,
^ CEPC = 0.49% <cit.>, ^ ILC = 0.38% <cit.>,
^ FCCee = 0.16% <cit.>.

In the limit κ=1, a conservative bound on (h→) can be derived.  We summarize these
bounds for future lepton colliders in Tab. <ref> (denoted by "est") derived from
the precision goals of κ_Z.  The resulting exclusion contours based on these
κ-independent BR bounds are almost identical to those from constraining the product of
κ and the BR in Eq. (<ref>) because κ is very close to 1 for the
sin^2θ of interest and the κ-dependence of the product is dominated by the
κ-dependence of (h→ϕϕ).

For the FCCee, the Higgs BR into NP has been fitted in a multi-κ fit without requiring
κ_V≤ 1, resulting in (h →) ≤ 0.48% at the 1σ
level <cit.>. The order of magnitude agrees well with our 2σ estimate
of 0.64%, but due to these two restrictions we apply our result for a coherent comparability
among the lepton colliders.

For the LHC, (h→) has been worked out in a global Higgs coupling fit for different
luminosities in Ref. <cit.>. The bound based on the Run-1 data set is explicitly
calculated in the model with one universal coupling modifier and a NP branching ratio, and is also
presented at the 95% CL. In contrast, the projections for 300 and 3000 result from a
multi-dimensional κ fit and are reported only at the 68% CL. The two modifications with
respect to the Run-1 bound have opposite effects and partially compensate each other. A dedicated
global fit of Higgs couplings to the single-κ case is necessary in order to obtain 95% CL
bounds precisely applicable to the relaxion case.





  
Implication of the untagged Higgs decays on the relaxion parameter space
The fact that the coupling does not vanish in the limit of small sinθ—owing to the
cos^3θ term in Eq. (<ref>)—gives rise to a sinθ-independent bound
on the mass for small sin^2θ≲ 10^-3.  Fig. <ref> shows in blue
the projections arising from the BR(h→ untagged) bounds at future colliders. While
the present bound of the LHC Run 1 excludes ≳ 32, the HL-LHC is expected to exclude
≳ 20. These bounds may be significantly improved by future lepton colliders. In
particular, the potential of the FCCee with 10 allows for an exclusion of
≳ 12. Since the hϕϕ coupling is s_0-dependent, also the exact exclusion
contours depend on the endpoint of the rolling. However, by comparing numerically s_0 of
𝒪(1) and 𝒪(10^-2), we notice that the impact on the
BR(h→ϕϕ) is only very mild. Consequently, these strong bounds can be
regarded as quite robust within the relaxion framework.

In contrast, for general Higgs portal models, the hϕϕ coupling has a different dependence on
the mixing angle and accordingly, the region of >35 is not necessarily ruled
out. Consequently, searches for light scalars with such masses remain highly relevant.




  §.§.§ Higgs self-coupling


The knowledge of the the Higgs self-coupling λ is crucial for the understanding of the Higgs
potential.  As shown in Eq. (<ref>), the Higgs-relaxion mixing alters λ. However,
the sensitivity of present and future colliders to λ via Higgs pair production is
limited. While relative deviations Δλ/λ = λ/λ_ SM-1 from the SM
value of 50% are estimated to be in the reach of the HL-LHC <cit.>, the ILC might
reach 10% <cit.>[However, this estimate of 10% at the ILC is not
  presented with a detailed quantitative analysis.] whereas Ref. <cit.> obtains -34%
to +42%.  According to Ref. <cit.>, CLIC is expected to become sensitive to a
deviation of Δλ/λ=19% for an electron polarisation of P^-=-80% and combining
runs at 1.4 and 3 with luminosities of 1.5 and 2, respectively, whereas
Ref. <cit.> reports -18% to +28% at the 68% CL.  Combining the anticipated
measurements at the HL-LHC and the FCCee in a 13-parameter fit, Ref. <cit.> concludes
a sensitivity of 40%.

While the estimates of the different references may not be directly comparable, it becomes clear
that relatively sizeable deviations in the order of several 10% of the Higgs self-coupling from
its SM prediction will remain in accordance with uncertainties at future colliders. Yet, the
relaxion with a mixing of sin^2θ<0.1 induces deviations from λ_ SM of only less
than 10%. Consequently, the considered future colliders will not be able to constrain the relaxion
parameter space via this indirect probe.  A precision of 1% (0.1%) would be needed to
constrain sin^2θ<10^-2 (10^-3).



  §.§.§ Total Z width


Relaxion-mediated decays of the Z boson give rise to an additional contribution
Γ_Z^ NP to the total Z width. Since the total width measurement at LEP1 is in
agreement with the theoretical SM prediction <cit.>, a bound on the
NP contribution can be derived by limiting it to the combined experimental and theoretical
uncertainty. Currently the uncertainty is dominated by the experimental one of
δΓ_Z^LEP1=2.3 <cit.> while the theoretical one is
δΓ_Z^th=0.5 <cit.>. At TeraZ, the experimental uncertainty
is expected to shrink below δΓ_Z^TeraZ = 0.1 <cit.>,
necessitating further theoretical improvement. When including the so far missing 3-loop
contributions, the theoretical uncertainty is estimated to reach
δΓ_Z^th,3-loop=0.2 <cit.>. Hence, unless the improvement of
the theory uncertainty goes beyond this estimate, the theoretical and experimental uncertainties
will be comparable, and we will use the combined projection of
Γ_Z^ NP< 2×δΓ_Z^ th+exp≃ 2× 0.22 to set a bound at
the 95% CL. Saturating this bound by the relaxion-mediated 3-body decay Z→ϕ ff̅ for
vanishing -odd couplings, we obtain the constraint shown as a red line in
Figs. <ref> and <ref>.

This partial width can be expressed in terms of the -even and -odd couplings of the relaxion,
giving rise to the integrands ℐ, ℐ, respectively,

    Γ_Z^ NP   =
      Γ(Z→ϕ ff̅)
      = 1/m_Z^3 (2π)^3∑_fN_f^c ∫_4m_f^2^(m_Z-m_ϕ)^2dm_12^2 (ℐ + ℐ) .

The sum is over the fermion final states f of interest, and m_f and N_f^c are the mass and
number of colors of the respective fermion. We obtain the following -even and -odd
integrands

    ℐ   = ^2/18 α π/ √(1
                    -4m_f^2/m_12^2) √(λ_ K(m_Z, m_ϕ,
                    m_12))/(m_12^2-m_Z^2)^2
                    +m_Z^2Γ_Z^2×
         [(m_12^2+2m_f^2)(λ_ K(m_Z,m_ϕ,m_12)/m_12^2+12 
        m_Z^2)(a_f^2 + v_f^2) -72  m_f^2  a_f^2] ,
        
    ℐ   =α π/36(4π f)^2 √(1-4m_f^2/m_12^2) λ_ K^3/2(,,m_12)
                                [c̃_γ Z^2 Q_f^2 m_12^2+2m_f^2/m_12^4.
      
         . +c̃_ZZ^2/(m_12^2 +
                    2m_f^2)(a_f^2+v_f^2) - 6 m_f^2 a_f^2/(m_12^2-)^2+Γ_Z^2.
         .-2 c̃_γ Z c̃_ZZ Q_f  v_f/ (m_12^2 +
                    2m_f^2)(-m_12^2)/m_12^2[(m_12^2-)^2
                    + Γ_Z^2]] ,

respectively, with a_f, v_f being the axial- and vector-coupling of the fermion, Q its
electric charge, θ_W the Weinberg angle, and the Källén function is given by 

    λ_ K(m_1,m_2,m_3) = -(m_1+m_2+m_3)(m_1+m_2-m_3)(m_1-m_2+m_3)(-m_1+m_2+m_3) .

In the calculation of ℐ we neglected the effect of the loop-induced ϕ Zγ
interaction. We find agreement of the -even contribution and the result reported in
Ref. <cit.>. For =0 and summing over all kinematically allowed SM fermions we
obtain

    Γ/MeV   ≈ 23  sin^2 θ + 2.1×10^-5(TeV/f)^2 (1.1 c̃_ZZ^2 - c̃_ZZc̃_γ Z + 47 c̃_γ Z^2) ,

where the approximation is valid for ≲ 0.5 for the -even and
≲ 10 for the -odd contribution. However, the latter contribution is
completely negligible in the range of interest. This justifies showing only the sin^2θ
contribution in Figs. <ref> and <ref>. In principle, also other
electroweak precision observables could play a role.  Yet, Ref. <cit.> showed that
for the considered masses they do not constrain small enough mixing angles in the Higgs portal.




 §.§ Direct probes


In a complementary way to the indirect bounds discussed in the previous section, the relaxion
parameter space can also be constrained by direct searches in various production modes as discussed
below.



  §.§.§ Pair production in Higgs decay



While Higgs coupling fits are sensitive to the BR of h→ϕϕ irrespective of the decays of
the relaxions (see Sect. <ref>), one can also look directly for the relaxion decay
products.  Each relaxion, pair-produced in the Higgs decay, further decays into a pair of fermions
f, photons γ or gluons g resulting in a four-particle final state F. ATLAS and CMS
search for such signatures and report -dependent bounds on
σ_h/σ_h^ SM× BR(h→ϕϕ→ F), which can be
compared to the prediction in the relaxion framework and thereby be translated into a bound on
sinθ and .  In Tab. <ref> we summarize the status of these exotic Higgs
decay searches performed by ATLAS and CMS during Run 1 and Run 2 of the LHC, listing the final
states, the considered data set and mass range of . We conclude that none of the current
searches is sensitive enough to probe parts of the relaxion parameter space displayed in
Figs. <ref> and <ref>, i.e. 5≤≤ 35
and 10^-5≤sin^2θ≤ 10^-1.  Moreover, we estimate the potential reach of the
HL-LHC with 3000 by rescaling the current limits by the ratio of luminosities and, in the case
of Run 1 limits at √(s)=8, additionally by the ratio of Higgs production cross sections at
8 and 13 in the dominant channels <cit.>. In Tab. <ref> we state the
relaxion mass corresponding to the vertical asymptote of the exclusion contour in the
(, sin^2θ) plane.

The strongest direct bound at the HL-LHC is expected in the bbττ channel excluding >26 which is presented
in Fig. <ref> (orange, dashed). The projections for bbμμ and 4b with
Vh, V=W,Z production are similar, but somewhat weaker (see ^ HL in Tab. <ref>) and therefore not
shown in the overview plot.  Neither the ττττ and ττμμ nor the recent
γγ g g final states are expected to constrain the displayed parameter plane at the
HL-LHC, based on this extrapolation from the current LHC data.

The 4μ final state covers low masses of a few GeV. However, sinθ is in part of this mass
range already strongly constrained by rare B-decays. Furthermore, the 4μ bounds at
13 <cit.> are reported for a certain model only, namely dark
gauge bosons, a 2-Higgs doublet model with an additional singlet (2HDM+S) or the NMSSM for a fixed
value of tanβ, making the translation to other models less straightforward.  A
model-independent presentation of the updated bounds in all channels would be helpful.

As a conclusion, the only direct channels evaluated here that have the potential to constrain the
displayed parameter space at the HL-LHC, are bbll, l=τ, μ and 4b in Vh production due
to the sizeable branching fractions. However, the strongest vertical asymptote is at >26
whereas Higgs coupling fits at the HL-LHC can exclude >20. Hence at the HL-LHC, there is
no direct search for h→ϕϕ→ F that is stronger than the indirect HL-LHC bound.

A different picture emerges at lepton colliders running at √(s)=240 for which projections
of the upper bound on exotic Higgs branching ratios in particular in hadronic final states have been
worked out in Ref. <cit.>. Such signatures are very hard to distinguish from background
at hadron colliders, but can be promising at lepton colliders.  Taking CEPC with 5 as a
benchmark, a strong bound of (h→ϕϕ→4b)<3×10^-3 for
10≤≤ 60 has been obtained. The resulting exclusion contour in the
(, sin^2θ) plane is shown in Fig. <ref> (orange, dotted), excluding
>11.  Hence, this direct search is expected to constrain the relaxion parameter space
more strongly than the bound derived from the anticipated precision of the Higgs couplings, see the
corresponding CEPC contour in Fig. <ref> that excludes >16.  Likewise,
for the FCCee and the ILC at the same energy, a similar bound on h→4b is expected according to
their luminosities.  In contrast, it is not easily transferable to CLIC due to its higher energy far
above the Zh threshold. A dedicated analysis would be necessary and useful to determine the reach
of CLIC via these hadronic search channels of exotic Higgs decays.




  §.§.§ Production at the LHC



Similarly to the Higgs, the dominant production modes for the relaxion at the LHC are


  * gluon fusion: p  p →ϕ 

  * relaxion strahlung: p  p → Z ϕ, Wϕ 

  * { tt̅, bb̅}-associated production:
  p p →{ tt̅ϕ, bb̅ϕ}

  * vector boson fusion (VBF): p p →ϕ j j

In the left panel of Fig. <ref> we present (as solid lines) the production cross
sections σ(pp→ X) depending on  in different channels X=ϕ (gluon fusion),
Wϕ, Zϕ, tt̅ϕ, bb̅ϕ and ϕ jj (VBF) for sin^2θ=1. In all
processes (beside the resonant ϕ production), we require a minimal transverse momentum of the
relaxion of p_T(ϕ)>20 and evaluate the cross sections at √(s)=13.  The gluon
fusion cross section was calculated at N^3LO and with resummation up to N^3LL using
 <cit.>.  All other
cross sections for processes at the HL-LHC were obtained with
 <cit.> at NLO.  Compared to the production cross section of
the Higgs with a mass of m_h=125, which is of the order of pb, there is at most an enhancement
of two orders of magnitude (reached only for Wϕ at the lowest considered ) for
sin^2θ=1 and setting all -odd couplings c̃_i=0. The LEP search, however,
already constrains sin^2θ≲ 10^-2 in the considered mass range. As a consequence,
the ϕ searches need to be targeted at cross sections smaller than the SM values of the
125 Higgs.

Regarding the decay modes, the right panel of Fig. <ref> shows the branching ratios
of ϕ for only -even couplings. The leptonic channels, i.e. ϕ→ττ,
μμ have a similar branching ratio for >2m_b as for the 125 GeV Higgs and therefore
the expected rates are small. The other clean channel, ϕ→γγ, is suppressed by a
factor of 100 to 2 for  between 5 and 100 GeV, respectively, compared to the BR at a mass of
125 GeV due to the ^3-dependence.  Hence, an observation of this final state might only
become feasible if this partial width is enhanced by the -odd coupling
c̃_γγ.  While below the bb̅ threshold the decay into cc̅ is
dominant, BR(ϕ→ bb̅) approaches 1 where kinematically allowed. However, hadronic
final states of low-mass resonances pose a severe challenge. CMS has already performed a vector
resonance search with hadronic final states, yet only for >50. Further experimental
efforts focused on lower masses would be helpful. Regarding the status at the LHC, CMS searched for
pseudoscalars in the 25–80 mass range produced in association with b b̅ in the
19.7 fb^-1 data set from the 8 TeV run (see Ref. <cit.> for
τ^+ τ^- and Ref. <cit.> for μ^+ μ^-decay modes). Assuming a
similar efficiency for scalars the present bound is weaker than the one from LEP.  Moreover, a naive
rescaling by the increased luminosity at HL-LHC does not improve on the LEP bound
either. Nonetheless, an improvement beyond the increased luminosity might provide valuable input.

Ref. <cit.> investigates in a phenomenological study the discovery prospects for a
new pseudoscalar with a mass between 20 and 100 GeV produced via t t̅ϕ with ϕ
decaying into b b̅.  Ref. <cit.> instead studied the LHC phenomenology of a
light scalar in the same mass range mixed with the Higgs and concluded that the LHC has the
potential to slightly improve the LEP constraints for this channel for >80.  However, the
cross sections to be probed by this analysis are too large to be realized in a relaxion framework.
Recently, the sensitivity of the LHC to light axion-like particles produced in gluon fusion has been
studied for decays into γγ <cit.> and
τ^+τ^- <cit.>.  Their projected sensitivity could set competitive
bounds for -odd couplings, but for -even couplings the cross sections within reach
correspond to a mixing angle already probed by LEP.  Potentially more promising production modes are
the associate production with a gauge boson, i.e. Zϕ or Wϕ, which yields larger cross
sections, as shown in Fig. <ref>.

So far we have considered the case where the dominant couplings are the ones given by the
Higgs-relaxion mixing.  As we previously discussed, the relaxion also has -odd couplings with
the SU(2) gauge bosons which cannot only change the decay pattern, but also the production. In
particular the production p p → W(Z) ϕ can receive an additional contribution such as

    σ( p p →ϕ Z )_5    = 2.7× 10^4 ^2  + 2.2 (1/f)^2 (1.9 c̃_ZZ^2
                                                            + 1.1 c̃_ZZc̃_γ Z + c̃_γ Z^2)
    σ( p p →ϕ Z )_35   = 1.0× 10^4  ^2  + 2.0 (1/f)^2 (1.8 c̃_ZZ^2
                                                             + 1.0 c̃_ZZc̃_γ Z +
                                                             c̃_γ Z^2)
    σ(p p →ϕ W)_5    = 6.6× 10^4 ^2  + 2.7
                                                          (1/f)^2  c̃_WW^2 
    σ(p p →ϕ W)_35   = 2.3× 10^4 ^2  + 2.4
                                                           (1/f)^2  c̃_WW^2  ,

where the subscripts 5 and 35 denote  in GeV.
The partonic cross sections were calculated analytically with  version
9.2.0 <cit.> and folded with the LO parton distribution from
 version 3.1 <cit.>, interfaced with  version
2 <cit.>. As a cross check of our calculation, we compared the -even part with
the cross section calculated with  <cit.> and found agreement.






  §.§.§ Production at electron colliders



Due to the clean environment, lepton colliders are not only able to explore the relaxion parameter
space via precision measurements as discussed in the previous section, but also via direct relaxion
production.

When considering the -even coupling via mixing with the Higgs, there are two main production
channels for relaxions at lepton colliders, ϕ Z and ϕν_e ν̅_e via W-fusion. The
cross section for the third production mode, ϕ e^+ e^- via Z-boson fusion, is about one order
of magnitude smaller than that of W-fusion and will therefore only play a negligible role. For
ϕ Z associated production, the polarized cross section is given by

    σ_ϕ Z   = πα^2  ^2/24 s𝒫^2_P^+P^-/[λ_ K(√(s),m_Z,m_ϕ) + 12m_Z^2 s]
                        √(λ_ K(√(s), m_Z,m_ϕ))/s[(s-m_Z^2)^2+m_Z^2Γ_Z^2] ,

where P^± denote the positron and electron polarization, respectively, and

    𝒫_RL   =-2/
    𝒫_LR   =1-2/
    𝒫_RR   =𝒫_LL=0 ,

and λ_ K is given in Eq. (<ref>).
In the limit s≫ m_Z ≫, the expression simplifies to

    σ_ϕ Z   →πα^2  ^2/24 s𝒫^2_P^+P^-/[1+11 m_Z^2/s + 𝒪(^2/s, m_Z^4/s^2)] .

For √(s)=240 and =0,
the following approximate numerical values are obtained from the full expression in Eq. (<ref>)

    σ_ϕ Z   ≈.
                                              
                                                1.6    for  P^+P^-=LR
    
                                                1.0    for  P^+P^-=RL
                                              } ^2 .


The existing limits on this process from LEP2 (ALEPH, DELPHI, L3, OPAL) <cit.> were
presented in Ref. <cit.>.  Our approach here is to estimate the reach of future
lepton colliders running at the Z-pole or above by rescaling LEP2 bound by the square root of the
ratio of luminosities.  For the FCCee we assume a luminosity of 10 at a center-of-mass energy
of √(s)=240.  When extrapolating the LEP2 limits, we neglect the difference in cross
section due to the different center-of-mass energies within the various LEP2 runs (192–202 GeV)
and between LEP2 and the FCCee.  The ϕ Z cross section at √(s)=240 is by a factor of
about 2 smaller than at √(s)=𝒪(200).  Yet, at the same time, also the cross
sections of the background processes decrease so that our rescaling by luminosity ratios only is
justified as a rough approximation.

For lighter ϕ-masses it is more promising to consider the three-body decay Z→ℓℓϕ
at the run at the Z pole. The TeraZ limit is rescaled from the LEP1 (L3)
limit <cit.> by a factor of
√(N_Z^L3/N_Z^TeraZ) where
N_Z^L3=4.4× 10^6 and N_Z^TeraZ = 10^12 are the number
of Z bosons at the respective collider/experiment.  These estimates are presented by the green
dashed lines in Fig. <ref>.This is a conservative extrapolation of cut-and-count analyses and will be certainly outperformed by
modern analyses. Compared to the LEP analyses they will benefit from improved flavour-tagging
techniques, better detector design and also the advances in computational methods,
i.e. usage of machine-learning tools. The final states might be different from those
considered in the LEP analyses when taking -odd couplings to photons or to leptons into
account.  In order to enable an estimate of the potential of these machines independent of the
relaxion decay modes, we show in Fig. <ref> as dashed lines the unpolarized cross
section σ(e^+e^-→ Zϕ) at √(s)=240 for sin^2θ=1. The lower cross
section with respect to those at the LHC will be compensated by the clean environment, the possible
enhancement by polarization as well as the large planned luminosity, making future lepton colliders
powerful machines for ϕ production.

The production of a light scalar ϕ with >20 in association with a Z was studied in
the context of the ILC in Ref. <cit.>. In particular, for a luminosity of
ℒ=2000 and with polarized beams, the LEP2 bounds can be significantly
improved. Their result is shown as the green dotted line in Fig. <ref>.

The second relevant production mode is W-fusion leading to a ϕν_eν̅_e final state. The
cross section at √(s)=240 with p_T(ϕ)>10 is calculated with
 yielding

    σ(e^+ e^- →ϕν_e  ν̅_e)_5   = 61 s_θ^2 
    σ(e^+ e^- →ϕν_e  ν̅_e)_35   = 45 s_θ^2  ,

where the subscripts 5 and 35 denote  in GeV. This process played a marginal role at LEP, but
it is important to be included for instance in searches targeting missing energy final states. For
comparison, the cross sections for Z-fusion leading to the ϕ e^+ e^- final state are

    σ(e^+ e^- →ϕ e^-  e^+)_5   = 4.8 s_θ^2 
    σ(e^+ e^- →ϕ e^-  e^+)_35   = 3.6 s_θ^2  .

These values are calculated with an additional cut on the transverse momentum of the two leptons of
10 and the pseudorapidity η<2.5.  They are about one order of magnitude smaller than the
corresponding cross sections for W-fusion.





  §.§.§ Probing the -odd couplings and the relaxion CPV nature




We have briefly discussed above the possible (CPV) contributions of the relaxion to the electron EDM
that is proportional both to its -even and -odd coupling.  Here we consider the potential to
probe the unique  properties of the relaxion from processes involving its production at lepton
colliders.  We first discuss the case of measurements involving ϕ Z final states, and then move
to examine those with ϕγ . As for the former channel, the production from mixing,
i.e. through the -even vertex, is dominated by the relevant coupling to Z while the
latter is loop-induced and is described by a dimension-five operator in the effective field theory.
This is in contrast with the pseudoscalar interactions that are always induced by dimension-five
operators.



  
Production via -odd couplings
As for the ϕ Z final state we have already discussed above the contribution from mixing and
thus show now the production cross section via -odd couplings:

    σ̃_ϕ Z   =α/48 (4π f)^2
          sλ_ K(√(s),m_Z, m_ϕ)^3/2×
         [c̃_ZZ^2/(s-m_Z^2)^2
            +m_Z^2Γ_Z^2𝒫_P^+P^-^2 
            + 4 c̃_γ Z^2/s^2
            + 4 c̃_γ Z c̃_ZZ (s-m_Z^2)/s [(s-m_Z^2)^2+m_Z^2Γ_Z^2]𝒫_P^+P^-],

 where the polarization factors are the ones given in Eq. (<ref>), and for
 s≫≫, it simplifies to

    σ̃_ϕ Z   →α (2 c̃_γ Z + c̃_ZZ𝒫_P^+P^-)^2/48 (4π f)^2[1- m_Z^2/s6 c̃_γ Z + c̃_ZZ𝒫_P^+P^-/2 c̃_γ Z + c̃_ZZ𝒫_P^+P^- +𝒪(^2/s,m_Z^4/s^2)] .

For √(s)=240 and =0, Eq. (<ref>) yields

    σ̃_ϕ Z   ≈ 1.0(TeV/f)^2
                                              
                                                0.60 c̃_ZZ^2 + 1.6 c̃_ZZc̃_γ Z + c̃_γ Z^2    for  P^+P^-=LR
    
                                                0.39 c̃_ZZ^2 - 1.3 c̃_ZZc̃_γ Z + c̃_γ Z^2    for  P^+P^-=RL
     .

These values are about two orders of magnitude smaller than the Zϕ cross section from the
-even interaction given in Eq. (<ref>), for c̃_i ∼. Therefore the -odd interaction only plays a significant role in this process when the
mixing angle is much smaller than the c̃_i for f=1 TeV.  Allowing for both -even and
-odd couplings to be present, a measurement of the total cross section can only constrain one
(quadratic) combination of these.  However, an angular analysis of the distributions of the Z
decay products involving an up-down imbalance is sensitive to the  asymmetry of the decay (see
e.g. Ref. <cit.>). Consequently, one can in principle probe both the
-even and -odd components of the cross section.

In addition to the Z ϕ production we also consider the γϕ production. We first note
that for the latter the mixing contributions only arise at the one-loop level and thus this channel
is much less sensitive to these when compared to the Z ϕ channel. For simplicity we therefore
focus on the -odd couplings with the corresponding cross section given by

    σ̃_ϕγ =α (s-m_ϕ^2)^3/48 s (4π
          f)^2[
        c̃_γ Z^2
          𝒫_P^+P^-^2/(s-m_Z^2)^2 + m_Z^2Γ_Z^2 +
        4 c̃_γγ^2/s^2 + 4 c̃_γ Zc̃_γγ(s-m_Z^2)𝒫_P^+P^-/s [(s-m_Z^2)^2+m_Z^2Γ_Z^2]] ,

where the polarization factors 𝒫_P^+P^- are given in Eq. (<ref>). In the
limit s≫ m_Z^2 ≫ m_ϕ^2 it simplifies to

    σ̃_ϕγ   →α (2 c̃_γγ + c̃_γ Z𝒫_P^+P^-)^2/48(4π f)^2[1+m_Z^2/s2 c̃_γ Z𝒫_P^+P^-/2 c̃_γγ + c̃_γ Z𝒫_P^+P^- +𝒪(^2/s,m_Z^4/s^2)] ,

which agrees with the result presented in Ref. <cit.> within their assumptions. For
√(s)=240 and =0, Eq. (<ref>) yields

    σ̃_ϕγ   ≈ 1.6(TeV/f)^2
                                  
                                    0.6 c̃_γ Z^2 + 1.6 c̃_γγc̃_γ Z + c̃_γγ^2    for  P^+P^-=LR
    
                                    0.39 c̃_γ Z^2 - 1.3 c̃_γγc̃_γ Z + c̃_γγ^2    for 
                                    P^+P^-=RL
                                   .


Finally, rare Z decays into ϕγ involve c̃_γ Z. For this partial width we
obtain

    Γ(Z→ϕγ ) = c̃^2_γ Zm_Z^3/96 π(4 π f)^2(1-m_ϕ^2/m_Z^2)^3

in agreement with e.g. Ref. <cit.>.



  
Constraining the -odd couplings


The above cross section prediction for σ̃_ϕ Z as well as the decay widths of
Z→ℓℓϕ and Z →ϕγ can be used to constrain the -odd couplings
c̃_ZZ and c̃_γ Z. We evaluate the implications of the LEP measurements at
and above the Z-pole and estimate the sensitivity at TeraZ and the FCCee.

Regarding Z→ℓℓϕ, we translate the L3 bound from LEP1 <cit.> on
sin^2θ shown in Fig. <ref> into bounds on the -odd c̃_ZZ
coupling. For this we require

    _Z→ℓℓϕ(sin^2θ_ max)
     !=_Z→ℓℓϕ((c̃_ZZ/f)^2_ max) ,

i.e. ignoring differences in the distributions. Here  refers to the decay through
mixing where all -odd couplings are neglected whereas  refers to the decay based
on the -odd coupling c̃_ZZ. For TeraZ, we take the rescaled LEP1 limit on
sin^2θ.

Likewise we proceed for Z→ϕγ where only c̃_γ Z plays a role. The
strongest bound on (Z→ϕγ)×(ϕ→ ff) is reported by
ALEPH <cit.> for f=τ. Here we assume the decay of Z into ϕγ via
c̃_γ Z and the decay of the relaxion into ττ via mixing with the BR as shown
in Fig. <ref>.  Since (ϕ→ ff) is independent of sin^2θ, this
approach is valid for sufficiently small sin^2θ. We rescale the ALEPH limit by the ratio of
Z bosons used in the analysis <cit.> to the 10^12 expected at TeraZ.

Concerning e^+ e^-→ϕ Z, we analogously require the number of Z bosons, N_Z, produced in
this channel via c̃_i to be equal to N_Z produced via mixing at the sin^2θ that is
maximally allowed by LEP2 or its rescaling to FCCee, respectively. In order to obtain N_Z at LEP2,
we sum over the products of luminosity times cross section evaluated at the energies of the
respective runs <cit.>.

Fig. <ref> shows the current limits and future projections on c̃_γ Z
and c̃_Z Z.  While the runs around 200 set stronger bounds on c̃_Z Z than
those at m_Z, c̃_γ Z is best constrained at the Z-pole.

Finally at very high energy lepton colliders like the high-energy stage of CLIC, W-boson fusion
becomes an important process to produce axion-like particles <cit.>. This could be
used to set bounds on c̃_WW.



 §.§ Comparison of direct and indirect probes

 


The comparison of the direct and indirect probes in the summary plot of Fig. <ref>
highlights the complementarity and perspectives in exploring the relaxion framework at colliders.
There are ample regions of parameter space where the relaxion can be discovered both via deviations
in precision measurements of Higgs properties and via direct production.

Focusing on the discovery potential of the future colliders, we notice that the searches for Higgs
decays into a pair of relaxions can probe relaxion masses rather independently of the mixing angle.
Such bounds can be complemented by direct searches in other channels and confirm or rule out a
potential indirect evidence.  Moreover, direct searches are the only way to explore relaxion masses
below 20 GeV at the (HL-)LHC, providing a strong motivation to try to push the sensitivity below
this mass.  In this context, τ final states are particularly motivated due to the large
branching ratio in Higgs portal models and a reasonable reconstruction efficiency. At the price of a
lower branching ratio, but with the benefit of a cleaner signature, also other leptonic final states
should be explored, see Fig. <ref>.

Despite the precision at TeraZ and the expected improvement of the theoretical uncertainty, the NP
contribution to the total Z-width will not set new limits beyond the already excluded
range. Nevertheless, direct searches at TeraZ for the decay of Z→ϕ ff̅ do play an
important role in constraining a light ϕ.

The `collider region' of the relaxion is bounded from below by the LHCb bound of up to 5 GeV. While
for general scalars any higher masses are relevant, the relaxion window ends at 32 GeV on account
of the bound set by untagged Higgs decays.

The maximal relaxion-Higgs mixing angle as a function of the mass is limited by
Eq. (<ref>) indicated by the black line. Hence the current bounds hardly probe the
viable region whereas the projections for all of the considered future colliders cut deeply into
this physical parameter region.



§ CONCLUSION



In this paper we evaluated the potential sensitivity of future colliders to light scalars, motivated
by the relaxion framework.  The relaxion phenomenology is driven by two kinds of couplings: The
first kind is in their simplest form identical to those of a -even Higgs portal and originates
from the relaxion mixing with the Higgs; the second kind is similar to those of -odd
axion-like-particles and typically arises due to an anomalous backreaction sector. We studied
indirect effects and direct production modes both at the HL-LHC and at future lepton colliders,
namely the ILC, CLIC, the CEPC and the FCCee running at the Z-pole and above the Zh production
threshold. Our results are applicable to a large class of models beyond the relaxion framework.
Light (pseudo-)scalars can arise in a variety of models, for example Higgs portals, 2HDM with an
additional singlet, and supersymmetric versions thereof, as well as in models with axion-like
particles.


The fact that the relaxion acquires both -even and -odd couplings makes it an interesting model
to study, in particular given the powerful capabilities of future colliders.  While focusing mainly
on the implications of the -even interactions, we also point out in which channels the
-odd couplings can influence the collider phenomenology and derive constraints on them.

Among the indirect probes, we evaluated the NP contribution to the total Z-width, in this
case from the decay Z→ ff̅ϕ, whose current bound from LEP1 will be significantly
improved at TeraZ. Furthermore, we studied the sensitivity to the hϕϕ coupling via untagged
h→ϕϕ decays, resulting in strong, though relaxion-specific bounds.

Regarding the direct probes, we considered gluon fusion, associated production as well as
explicit searches for h→ϕϕ. We provided analytic expressions for polarized cross sections
at lepton colliders and semi-analytic functions for the processes at the HL-LHC, involving both the
-even and -odd couplings.

Our main findings are:


  * We chart the physical parameter region of relaxion models, expressed via the relaxion-Higgs
  mixing angle as a function of the mass. We find that only mixing angles smaller than twice the
  ratio of mass to the Higgs vacuum expectation value can describe actual relaxion models. It is
  quite interesting that the HL-LHC and future colliders are capable to significantly probe this
  "physical relaxion" region, albeit only for a very massive relaxion.

  * The HL-LHC can probe relaxion masses roughly above 20 by untagged Higgs decays. However,
  its discovery prospects in direct channels via its mixing with the Higgs are limited due to the
  large background for low-mass resonances and suppressed branching ratios for the clean final
  states.

  * Future lepton colliders have the potential to significantly improve on existing limits and
  HL-LHC projections on the parameter space via direct and indirect channels.  In particular masses
  roughly above 10, and mixing angles above few times 10^-3 will be probed by exotic Higgs
  and exclusive Z decays respectively.

  * The -nature of ϕ can be deciphered by the interplay of -even and -odd
  driven signals. In particular the interplay of collider observables (such as ϕ Z and possibly
  ϕγ production as well as the angular distribution of decay products) and EDMs is crucial.





§ ACKNOWLEDGMENTS
 
We would like to thank Liron Barak, Roberto Franceschini, Rick S. Gupta, Marumi Kado, Yevgeny Kats, Yosef Nir, Diego Redigolo, Francesco Spano and Tim Stefaniak for useful discussions.
The work of GP is supported by grants from the BSF, ERC,
ISF,  Minerva Foundation,  and  the  Weizmann-UK  Making  Connections Program.
EF is supported by the Minerva Foundation.




JHEP


