Varying gauge couplings and collider phenomenology

In this paper we investigate a natural extension of the Standard Model that involves varying coupling constants. This is a general expectation in any fundamental theory such as string theory, and there are good reasons for why new physics could appear at reachable energy scales. We investigate the collider phenomenology of models with varying gauge couplings where the variations are associated with real singlet scalar fields. We introduce three different heavy scalar fields that are responsible for the variations of the three gauge couplings of the Standard Model. This gives rise to many interesting collider signatures that we explore, resulting in exclusion limits based on the most recent LHC data, and predictions of the future discovery potential at the high-luminosity LHC.

In this paper we investigate a natural extension of the Standard Model that involves varying coupling constants. This is a general expectation in any fundamental theory such as string theory, and there are good reasons for why new physics could appear at reachable energy scales. We investigate the collider phenomenology of models with varying gauge couplings where the variations are associated with real singlet scalar fields. We introduce three different heavy scalar fields that are responsible for the variations of the three gauge couplings of the Standard Model. This gives rise to many interesting collider signatures that we explore, resulting in exclusion limits based on the most recent LHC data, and predictions of the future discovery potential at the high-luminosity LHC.

I. INTRODUCTION
In a fundamental theory, there should be no free parameters. For example, in string theory, all coupling constants are given by vacuum expectation values (vevs) of scalar fields known as moduli fields. In top-down constructions of low-energy effective theories from string theory, the fundamental constants should in principle be derived from such vevs at a high scale, but there are many obstacles involved in constructing the Standard Model (SM), or extensions thereof, in this way. In this paper, we instead want to consider the idea of couplings set by scalar fields in bottom-up constructions, where the low-energy theory is constructed at lower scales as an effective field theory with a cutoff scale. We will here consider the Standard Model, where all dimensionless parameters could in principle be controlled by scalar fields.
One can make the case that this argument in favor of varying coupling constants is as strong as the theoretical arguments in favor of supersymmetry. Just as with supersymmetry, it is based on a very simple and reasonable assumption about the properties of the underlying fundamental theory. Estimating the relevant energy scales for these extensions is much more difficult. In the case of supersymmetry, it has been argued through naturalness that supersymmetry at low energies could solve the observed hierarchy of scales. In the case of varying couplings, the observed existence of a low energy scalar such as the Higgs boson can be taken as an argument suggesting that the scalars we study in this paper could have reasonably low masses as well. The Higgs boson is special in the sense that it controls the breaking of a gauge symmetry, but in addition it controls the masses of the standard model particles. Normalizing with the help of the Planck mass-which is an obvious thing to do in a fundamental theory-it means that the Higgs boson is better thought of as controlling a set of dimensionless parameters in the Standard model. This is similar to the scalars we have in mind. We cannot see any strong reason for why these * ulf.danielsson@physics.uu.se † rikard.enberg@physics.uu.se ‡ gunnar.ingelman@physics.uu.se § tanumoy.mandal@physics.uu.se other scalars must have masses way beyond the Higgs mass. Our conclusion is that there are good arguments in favor of conducting searches for such particles at the LHC and beyond.
There have been many studies of models where the coupling constants change over time [1]. In particular, models with a varying α EM have been widely discussed; in these models the associated scalar field is light and α EM has not been constant over the history of the Universe. Today the constraints on this kind of coupling variation of α EM are very strong, suggesting the mass of any associated scalar field to be very high. This is not the kind of variation we have in mind in this paper.
Here we are rather interested in the particle physics implications of dynamical couplings in terms of new scalar particles associated with the scalar fields.
Even if the effects of varying couplings can only be seen at high energies, they can still be of crucial importance for cosmology. This includes varying Yukawa couplings leading to a strong electroweak phase transition [2,3] or recently a varying strong coupling α S leading to QCD confinement in the early Universe and possible scenarios for baryogenesis [4]. The latter work introduced a nontrivial potential for the scalar field governing the strength of the strong interaction. In this way, it was possible in [4] to obtain a temperaturedependent vev for the scalar field. This, in turn, allowed for QCD to remain confined at very high temperatures, which might help explain baryogenesis. In [3], it was natural to pick the energy scales of the new physics to be of order TeV. Our general analysis outlines how data from LHC may be used to put meaningful constraints on such models.
An early, consistent varying coupling model was proposed by Jacob Bekenstein in 1982 [5,6]. This was a model for electromagnetism where the electromagnetic coupling e is given by a light scalar field. The model is gauge invariant and invariant under rescaling of the scalar field, and reduces to ordinary electromagnetism if the scalar field reduces to a constant value. We will refer to this model as the Bekenstein model (BM) and will discuss it further in Sec. II.
The Bekenstein model was proposed in the context of cosmology and has so far only been considered with this purpose. Recently, however, we introduced the idea [7] that the scalar field in the Bekenstein model could be heavy, with a mass on the TeV scale, such that it could potentially be discovered at the LHC. In [7], we concentrated on the electromagnetic coupling and the corresponding scalar field that couples to photons.
In this paper, we generalize this model. In the SM, the three gauge couplings g 1 , g 2 , g 3 of the gauge group G SM = U(1) Y ×SU(2) L ×SU(3) C as well as the Yukawa couplings of the fermions, the CKM matrix, and the Higgs quartic self-coupling λ H should be controlled by scalar fields. We will here as a first step restrict ourselves to the electroweak theory and will write down the model where the U(1) Y × SU(2) L gauge couplings g 1 and g 2 and the strong coupling g 3 are set by scalar fields. We will briefly discuss the further generalization to the Higgs and Yukawa sector, but will save the detailed implementation for future work. Ultimately, one would like to replace all dimensionless constants of the SM with values set by their corresponding scalar fields.
The paper is organized as follows. In Sec. II, we briefly discuss the construction of varying finestructure constant theory. Following the same method, we discuss in detail the varying electroweak gauge couplings and varying strong coupling models. In these models, the variations of three gauge couplings of the SM are promoted to three different heavy scalar fields. In Sec. III, we show the branching ratios (BRs) and production cross sections of these scalars at the LHC. We derive exclusion limits of the free parameters of these models in Sec. IV by using the latest LHC data. In Sec. V, we estimate the projected reach at the highluminosity LHC (HL-LHC). Finally, we summarize and conclude our results in Sec. VI.

II. THE MODEL
In Ref. [7], we considered the model of varying electromagnetic (EM) coupling as our starting point. In this model, the EM coupling is not a fundamental constant but a function of space-time, i.e., e(x) = eε(x), where ε(x) is a dimensionless scalar field and where we have, for later convenience, extracted the vev e(x) = e such that ε = 1 at the minimum of the scalar potential. We use the following notation for the rest of the paper: a space-time varying constant κ(x) will be denoted by κ ≡ κ(x) and its vev by κ = κ. For instance, in this notation e ≡ e(x) and e is the constant value that we measure in vacuum. We also omit the argument "x" from a field for ease of notation. The field ε has noncanonical kinetic term of the form where Λ is an energy scale. This form of the kinetic term respects rescaling invariance, and is typical of kinetic terms for moduli fields in string theory. For particle physics applications it is convenient to deal with a scalar field φ that has physical mass dimension one and a canonical kinetic term (1/2)(∂ µ φ) 2 , instead of the noncanonically normalized ε. This can be done by first writing ε = exp (ϕ), with ϕ = ln ( e/e), and then rescaling by an energy scale Λ such that ϕ = φ/Λ. The field ε is expanded to the first order as ε 1 + ϕ = 1 + φ/Λ. Bekenstein showed in [5] that the theory defined using the field strength tensor accompanied by a modified gauge transformation of the field A µ given by is gauge and rescaling invariant. In particular, this means that in the ordinary QED theory, eA µ is always replaced by eεA µ . We define the ordinary QED field strength tensor as In [7], we derived the form of the modified QED Lagrangian suitable for particle physics analyses with the varying coupling e replaced by the constant coupling e and with couplings to the real scalar field φ of mass m φ introduced. To illustrate the main construction, let us consider a simple U(1) EM theory consisting of photon and a fermion f with EM charge eQ f . This derivation is briefly described in Appendix A. The resulting new terms for interactions of the photon with the fermion are (there are also interactions of the scalar with W ± , but we do not list these terms here; see [7]) Notice that the kinetic term for the photon field and the interaction term of the scalar with photons are not individually gauge invariant, but their sum is. We may use integration by parts and the equations of motion for F µν to remove the last term and rewrite the nextto-last term so that there is no direct coupling of the type φf f γ. The Lagrangian then instead takes the more familiar form (see Appendix B for details) with the well-known dimension-five operator coupling a singlet scalar to a gauge field. This form of the Lagrangian is related to (4) by a field redefinition, and the two forms are equivalent in the sense that they predict the same S-matrix elements, although the explicit amplitudes will appear different. The φF µν F µν term in Eq. (5) gives the new interactions in the model. The obtained theory is therefore regular QED, plus a real scalar field that couples to photons through a dimension-five operator, where interestingly the coupling is given by the inverse of the scale Λ and no other parameters.

A. Varying electroweak couplings
Above we discussed the model with varying EM coupling e. Now, we generalize this picture to varying gauge couplings in the electroweak (EW) sector. This generalization has been previously studied in Refs. [8,9] in the context of cosmology. Here, our primary aim is to study the phenomenology of this model in the context of the LHC. We let the two gauge couplings g 1 and g 2 of the EW sector (where g 1 and g 2 are the U(1) Y and SU(2) L gauge couplings, respectively) be specified by scalar fields. We consider two cases: one scalar theory (1ST), where both couplings are set by the same scalar field, and two scalar theory (2ST), where each coupling is set by a different field.
In 1ST, a single scalar field S 0 is responsible for the variation of both the couplings g 1 and g 2 . It is more natural to introduce two different scalar fields responsible for the variations of the two independent gauge couplings-this is the 2ST model which we are mainly interested in. However, for pedagogical reasons, we first discuss the 1ST model, due to its simple nature with only two free parameters.

One scalar theory
In the 1ST, the space-time varying gauge couplings are given by g i = g i ε, for i = 1, 2, where ε is the new scalar field that we write as ε = exp (S 0 /Λ 0 ) 1 + S 0 /Λ 0 by introducing a canonically normalized scalar field S 0 and the new physics scale Λ 0 . This is an obvious generalization of the formalism described before for varying EM coupling. Similar to Eq. (5), the interaction Lagrangian before EWSB can be written as where B 2 = B µν B µν and W k2 = W k µν W kµν (with k = 1, 2, 3) and B µν and W k µν are the field-strength tensors for the U(1) Y and SU(2) L gauge bosons, respectively. These interactions are nonrenormalizable effective dimension-five operators suppressed by the scale Λ 0 . The scalar potential of this theory consists of the SM Higgs doublet Φ and the new scalar S 0 , and the most general form is given by where µ, λ, λ n , and λ n are the free couplings of the potential. In this theory, S 0 will not get a vev 1 , but the Higgs doublet Φ gets a vev and breaks the EW symmetry as usual in the SM. After EWSB, an offdiagonal term S 0 h (h is the Higgs field) in the scalar mass matrix is generated from the λ 1 S Φ † Φ interaction. This leads to the mixing between the Higgs h and the new scalar S 0 . This mixing is severely constrained from LHC measurements [10][11][12]. Therefore, for simplicity, we do not consider the mixing of the new scalar with the Higgs and assume S 0 is the physical mass eigenstate. Note that the Weinberg angle θ w = arctan( g 1 / g 2 ) is not varying, since the variations of g 1 and g 2 are the same and cancel in the ratio. However, the gauge boson masses are varying, with variations given by These expressions for the gauge boson masses must be inserted in the gauge boson-Higgs Lagrangian coming from the covariant derivative in the kinetic term for the Higgs field, which leads to vertices of the type S 0 V V , S 0 hV V , and S 0 hhV V .
To rewrite the interactions after EWSB, we insert (6), and replace m Z and m W in Eq. (9) by the corresponding expressions from Eq. (8). After this replacement, we obtain the following interactions of S 0 with the physical gauge bosons: where The dots stand for terms with four-and five-point interactions that arise from the non-Abelian part of the W k µν field-strength tensor of Eq. (6). We do not explicitly show these interactions here, but they are included and are important to preserve the gauge invariance of the theory. The interactions of S 0 with the Higgs associated with two vector bosons can be readily obtained from Eq. (9).

Two scalar theory
In 2ST, the space-time variations of the two gauge couplings in the EW sector are given by two new scalars, Here, we introduce two new, canonically normalized fields S i and two corresponding scales Λ i . As before, the scalars S i will not get vevs, and mixing between S i and the Higgs can arise from the general scalar potential. Therefore, the S i fields are not, in general, physical mass eigenstates, but interaction eigenstates. Although we neglect the mixing of S i with the Higgs to satisfy the Higgs measurements, we do consider mixing between the two new scalars as there is no a priori reason to neglect S 1 ↔ S 2 mixing. This mixing can lead to interesting collider signatures which we will see later.
There will be two sources of interactions of S i with the SM fields. The first source is the interaction terms similar to Eq. (6) as follows: The second source of interactions comes from the dynamical nature of the Weinberg angle and the gauge boson masses. The tangent of the Weinberg angle θ w (denoted by t w ) is given by such that the sine ( s w ) and cosine ( c w ) are given by The squares of the gauge boson masses are then expressed as It is interesting to note that in this model the Fermi coupling G F is not dynamical, This is because the parameters in the Higgs potential do not vary in this model. The mass parameter µ 2 is dimensionful and is not allowed to vary, and since in this paper we focus on the gauge couplings we have not yet implemented a varying Higgs coupling λ H , but this should be considered in the future and would lead to interesting phenomenology.
To rewrite the interactions after EWSB, we replace, in Eq. (11) The expressions in Eq. (13) must also be inserted in the above linear combinations. This yields vertices with one or two new scalars and multiple gauge bosons. The vertices with more than one scalar are suppressed by higher orders of 1/Λ i , and we neglect these terms in the rest of the paper. The fields S i are, in general, not mass eigenstates, which could be linear combinations of them. We, therefore, define the mass eigenstates φ i by the following linear combinations: where α is a mixing angle determined from the free parameters of the scalar potential. As before, we replace m Z and m W in Eq. (9) by the corresponding expressions from Eq. (14). Finally, we collect the leading interactions of φ i with two gauge bosons in the following Lagrangian: Table I. A representative Feynman diagram for the φ i V 1 V 2 interaction originating from Eq. (17) is shown in Fig. 1.  Table I.

B. Varying strong coupling
Moving on to the strong interaction, we promote the strong coupling g 3 to a field the same way as done for the other interactions above, so that g 3 = g 3 ε 3 , and we have a new scalar S 3 defined through ε 3 = exp(S 3 /Λ 3 ) 1 + S 3 /Λ 3 . The interaction of the scalar with gluons is described by the dimension-five interaction where G a µν with a = 1, . . . , 8 is the usual gluon field strength tensor. This new interaction gives rise to S 3 gg, S 3 g 3 , and S 3 g 4 vertices at order Λ −1 3 . The scalar will interact with quarks in the same way as the Bekenstein scalar in the BM model interacts with electrons, as described in Sec. II, i.e., in the original formulation there will be a four-point vertex S 3q qg, but after the field redefinition is performed, there will not be any direct coupling to fermions. Instead, it only couples via internal particles in amplitudes. In general, S 3 can also mix with other scalars in the theory through the interactions in the general scalar potential and lead to many interesting signatures at the colliders. However, a systematic study will be extremely complicated due to large number of free parameters in that situation. Therefore, for simplicity, we discuss the varying strong coupling (VSC) theory with no mixing with the S 1 and S 2 scalars in the EW sector.

III. BRANCHING RATIOS AND CROSS SECTIONS
For the phenomenological analysis, we implement the dimension-five effective interactions shown in Branching ratios of S0 in the one scalar theory, as functions of MS 0 .
Eqs. (9), (11), and (18) in addition to the SM Lagrangian in FeynRules [13], to obtain the Universal FeynRules Output [14] model files suitable for the MadGraph [15] event generator. The fields S 1 and S 2 are not the physical ones, and they are replaced by the linear combinations of the physical states φ 1 and φ 2 as given in Eq. (16) to obtain the interactions in the mass basis. All possible new interactions that originate from these replacements are implemented in FeynRules, although not all of them are shown in Eq. (17). For photon initiated processes, we use the photon parton distribution function (PDF) obtained using the approach in [16,17] while for other partons, the MMHT14LO [18] PDF set is used. We use the default MadGraph dynamical factorization and renormalization scales in our analysis. Events are showered and hadronized using Pythia8 [19]. Detector environment is simulated using Delphes [20], which uses the FastJet [21] package for jet clustering. The antik T algorithm [22] is used for jet clustering with radius parameter R = 0.4.

A. Branching ratios
We use the generic symbol S = {S 0 , φ 1 , φ 2 } for the new scalars (mass eigenstates) in the 1ST and 2ST models. These scalars dominantly decay to two SM gauge bosons, i.e., S → V 1 V 2 . Interestingly, there is no S 0 γZ coupling in 1ST, and therefore S 0 → γZ decay is not possible. In 2ST, however, all four two-body decay modes are possible. Below, we show the expressions for the partial widths of the scalars φ i in the 2ST model, and similar expressions can be readily obtained for S 0 in the 1ST by changing the corresponding decay couplings. According to Eq. (17), the partial widths of the two-body decays of φ i are given by where V = {Z, W } and f Z = 1 and f W = 2. There are also three-body decay modes possible and among them φ i → γf f (where f is a charged fermion) mediated by an off-shell photon is substantial. For instance, Γ i γf f /Γ i γγ ≈ 0.18 for M φ = 1 TeV and increases to 0.21 for M φ = 3 TeV. We compute Γ i γf f analytically using the expressions given in Appendix C and also numerically using MadGraph [15]. In our analysis, the additional contributions of the three-body partial widths have been properly included in the total widths of the scalars.
The 1ST model is very economical and predictive as it has only the two free parameters M S0 and Λ 0 . In Note that the S 0 → γf f mode has a BR of about 4%. There are in total five free parameters in the 2ST, M φ1 , M φ2 , Λ 1 , Λ 2 , and α. In general, the Λ 1 and Λ 2 scales can be different, and we define a parameter ρ as Λ 2 = ρΛ 1 . In Fig. 3, we show the BRs of φ i as functions of c α ≡ cos α for five cases, ρ = 0.5, 1, 2, Λ 2 Λ 1 , and Λ 2 Λ 1 . For both φ 1 and φ 2 , the largest BR to γγ is about 55% for all c α in the limit Λ 2 Λ 1 , which we call the photophilic limit. Here, the W W mode disappears. In the other extreme limit, Λ 2 Λ 1 , the γγ mode disappears and the W W BR becomes the largest, about 70% for all c α . We refer to this limit as the photophobic limit. The other two modes γZ and ZZ are always present except for specific c α values depending on the ρ values. φ2 (e-h) into γγ, γZ, ZZ and W W , as functions of cos α for ρ = 0.5, 1, 2, Λ2 Λ1 (photophilic limit) and Λ2 Λ1 (photophobic limit). Here, φ1 and φ2 are the mass eigenstates in the two scalar theory, assuming M φ 1 ,φ 2 = 1 TeV for these plots.

B. Production channels at the LHC
The scalars in the 1ST and 2ST models can be produced from V 1 V 2 fusion processes, while the scalar S 3 associated with the strong coupling variation can be produced from gg fusion. In Fig. 4, we show the total production cross section contours of S 0 in the 1ST model in the M S0 − Λ 0 plane at the 13 TeV LHC. The total cross section includes the contributions from the γγ, ZZ and W W fusion production processes. We include photons as effective "partons" inside the proton. Therefore, the γγ initiated processes can be directly generated in MadGraph (using the syntax generate a a > S). However, in case of ZZ or W W fusion, the initial W and Z come from quarks. Therefore, in this case, the scalars are produced in association with at least two additional non-QCD jets. Similarly, in the 2ST model, the φ i are produced in association with at least one non-QCD jet in case of γZ initiated production.
The γγ, γZ, and ZZ production modes interfere among themselves. One has to take care to properly include the effects of these interferences, which should be regarded as components of signal, in the total production cross sections of S. We calculate the total cross sections by combining the following processes, avoiding any double counting: Process 1: γγ → S Process 2: γγ + γZ → Sj Process 3: γγ + γZ + ZZ + W W → Sjj.
In process 1, both photons come from the proton, as explained above. In process 2, one of the photons or the Z is radiated from a quark rather than coming from the proton, so there is one additional jet. In process 3, both incoming bosons are radiated from quarks, so this is a vector boson fusion process. The γγ → Sj channel in process 2 can also arise from the γγ → S with one initial or final state jet radiated, but the interference of γγ and γZ initiated processes can only be captured in process 2. Therefore, we take γZ initiated and (γγ + γZ)-interference contributions from process 2, while only the γγ initiated contribution is taken from process 1. Similarly, (γγ +ZZ) and (γZ +ZZ) interferences are captured in process 3 in addition to the pure ZZ and W W initiated We present cross sections for ρ = 0.5, 1, and 2. In the photophilic limit, Λ2 Λ1 we set Λ1 = 10 TeV and Λ2 → ∞, whereas in the photophobic limit Λ2 Λ1, we set Λ2 = 10 TeV and Λ1 → ∞. contributions. Note that W W initiated processes do not interfere with the others. See our previous paper [7] for details of the calculation. In Fig. 5, we show the total production cross sections of φ i in the 2ST model, computed using the method described above. In Fig. 6, we show the total production cross section contours of S 3 in the VSC model in the M S3 − Λ 3 plane. These cross sections include NLO QCD K-factors. Not all K-factors for the processes in Eq. (21) are available in the literature. For the V 1 V 2 fusion production modes, we assume a constant K-factor of 1.3 [23]. The actual K-factors can differ slightly from this value and can also vary for different masses of the scalars, but we have used a constant since it has very little effect on our results. For the gg → S 3 production, we compute the NNLO QCD K-factor using Higlu package [24] as a function of M S3 .

IV. EXCLUSION LIMITS
Depending on the parameter values of the model, the scalars can be substantially produced in the different V 1 V 2 → S fusion channels, followed by decays S → V 1 V 2 . Therefore, γγ [25,26], γZ [27,28], and V V [29,30] (ZZ and W W with different decays) resonance searches at the LHC can be used to constrain the parameter space for the scalars associated with the EW sector, and dijet resonance search data [31,32] can be used to constrain the parameter space of the VSC theory. These searches set upper limits on σ × BR for s-channel resonances of different spin hypotheses. By recasting these upper limits, we can derive exclusion limits on our model parameters. Various production modes of a heavy scalar and its decay to different final states have been discussed in Refs. [33,34].
In general, these searches are optimized for a resonance produced through gg fusion. For a specific final state, the selection cut efficiencies can be different for different production mechanisms. Moreover, these efficiencies can vary for different spins of the resonances and also for different widths of the resonances. Therefore, when deriving exclusion limits on the parameters by recasting upper limits on cross sections from various experiments, one must take care of the different cut efficiencies for different signal processes contributing to the same final state. For this, we use the following relation from Refs. [7,34]: (22) where the upper limit on the number of signal events is represented by N s . This can be expressed as the product of the signal cross section (for a particular production mechanism), branching ratio (for a particular decay channel considered in the analysis), the corresponding selection cut efficiency s , and the integrated luminosity L. In the case of different production mechanisms contributing to an experiment, N s = i (σ · BR) i × i × L where i runs over all the contributing processes.
We have checked that the latest ATLAS and CMS results produce very similar exclusion limits and, in this paper, in the interest of simplicity, we use only the ATLAS results in our analysis. Below, we briefly discuss each of the searches used.
ATLAS γγ resonance: [25] The ATLAS Collaboration has performed searches for spin-0 and spin-2 γγ resonances produced in gg fusion at the 13 TeV LHC with 37 fb −1 integrated luminosity. We recast the fiducial σ × BR upper limit for a spin-0 resonance with a 4 MeV width hypothesis (in the narrow width approximation). The definition of the fiducial signal region can be found in [25]. We extract the observed and expected σ × BR upper limits (along with 1σ and 2σ uncertainty bands) from HEPdata [35].
In Fig. 7, we show how the selection cut efficiencies vary as functions of the γγ resonance mass for different production mechanisms. To validate our analysis code, we compare our fiducial selection cut efficiency for the gg → S → γγ process with the corresponding ATLAS numbers and obtain close agreement within about ∼ 3%. We find that the cut efficiencies for different production modes do not vary much ( 15%). To simplify our analysis, we neglect this small variation in cut efficiencies while recasting this search.
ATLAS γZ resonance: [28] ATLAS has also performed searches for spin-0 and spin-2 γZ resonances in gg fusion with hadronic decays of the Z, with 36 fb −1 integrated luminosity. We recast the σ × BR upper limit derived for the spin-0 selection and, here too, we neglect the small variations in cut efficiencies. ATLAS V V resonance: [29] The final searches relevant for the EW scalars are diboson (both ZZ and W W ) resonance searches from gg and V V fusion with L = 36 fb −1 . As discussed above, the scalars S in our model can be produced in V V fusion (or in the kinematically similar γγ or γZ fusion processes) (21). In [29], ATLAS presented combined upper limits on σ × BR for different bosonic and leptonic decay modes of ZZ and W W , and we recast these combined limits for spin-0 resonance produced from V V fusion.
ATLAS jj resonance: [31] This is the excited quark q * search that uses the qg → q * → qg process, with 139 fb −1 integrated luminosity. In our model, we have a gg resonance produced from the gg → S 3 → gg process. We have checked by employing the selection cuts used in [31] that the cut efficiencies for these two processes differ only by a few percent. This is expected as the spin information of the resonance is mostly lost in the very boosted dijet system and they appear as two highly collimated jets in the detector. For simplicity, we neglect the effect of cut efficiency in the derivation of exclusion limits.
In Fig. 8, we show the lower limits on Λ 0 as functions of M S0 in the 1ST model from the γγ, ZZ, and W W resonance searches. Although the γγ BR is about 18% in 1ST, the strongest limit on Λ 0 comes from the γγ channel. For example, we obtain Λ 0 18 TeV for M S0 ∼ 1 TeV from γγ data whereas V V data can only exclude Λ 0 up to ∼ 5 TeV for M S0 ∼ 1 TeV.
In Fig. 9, we show the lower limits on Λ i in the 2ST model, as functions of M Si , from the V 1 V 2 resonance searches. In the photophilic limit (Λ 2 Λ 1 ), the γγ data can set strongest limits on Λ 1 , e.g., Λ 1 20 TeV for M S1 ∼ 1 TeV. Here, we only show limits for the case when S 1 and S 2 are mass eigenstates, i.e., with no mixing between them. Turning on S 1 ↔ S 2 mixing can relax the limits on Λ 1 since the total production cross section σ φ1 is proportional to c 2 α in the photophilic limit. In the photophobic limit (Λ 2 Λ 1 ) with no mixing, the S 1 decouples from the theory and the γγ data can only constrain Λ 2 up to a few TeV. Overall, the limits on Λ 2 are not very strong from the latest LHC data, but Λ 1 can be constrained strongly from the γγ data.
In Fig. 10, we show the lower limit on Λ 3 as a function of M S3 in the VSC. The 139 fb −1 dijet resonance search data from ATLAS [31] is very powerful in constraining Λ 3 . For example, we estimate the lower limit on Λ 3 to be about 16 TeV for M S3 ∼ 2 TeV.

V. REACH AT THE HL-LHC
In this section, we present the model parameter space that can be probed with more than 2σ confidence level (CL) at the HL-LHC, while at the same time satisfying the current data. We estimate the projected 2σ CL upper limits of the expected σ × BR by using a simple √ L scaling of the current expected val-   ues, i.e., where (σ × BR) obs and (σ × BR) exp are the observed and expected upper limits with present luminosity L p , respectively, and the future luminosity is denoted by L f . The theoretical (σ × BR) th is a function of model parameters and can be calculated in a given theory. When the background is very small compared to the signal, the scaling goes as L. In reality, the actual scaling lies in between √ L and L, so using √ L scaling is more conservative and does not exaggerate the prospects. At L = 3000 fb −1 , the expected (σ × BR) exp goes down by roughly 1 order of magnitude ( 36/3000 ≈ 0.11). When translated to the projected lower limit on the new physics scale Λ, one finds that it goes up by roughly a factor of 3. In case of the 1ST model, the obtained limits on Λ 0 as shown in Fig. 8 go up by a factor of 3 for L = 3000 fb −1 . Similarly, for varying strong coupling theory, the limit on Λ 3 in Fig. 10 goes up by a factor 2 for that luminosity.
The situation gets more complicated in the case of the 2ST model, where we have five free parameters. To derive the simplified projected limits, we perform a random scan over the three parameters Λ 1 , Λ 2 , and c α in the regions and fix the masses M φ1 , M φ2 = 1 TeV. In Figs. 11 and 12, we show the outcomes of the random scan for φ 1 and φ 2 , respectively. We present scatterplots in the planes of c α − Λ 1 (red dots), c α − Λ 2 (blue dots), and Λ 1 − Λ 2 (green dots) for the four types of V 1 V 2 resonances. These colored regions of the parameter space can be probed with at least 2σ CL at the HL-LHC with 3000 fb −1 integrated luminosity. They are obtained by using the condition in Eq. (23) for a specific V 1 V 2 resonance while satisfying others by the latest data, i.e., (σ × BR) th < (σ × BR) obs only. The white regions in all figures represent the presently excluded regions of the parameter space, and the gray background dots represent the regions which are allowed by the current LHC data. For φ 1 , the γγ resonance can very effectively probe the Λ 1 direction (for c α away from zero) but the Λ 2 direction is mostly insensitive. On the other hand, the W W resonance can be used to probe the Λ 2 direction but not the Λ 1 direction. For φ 2 , the γγ resonance can probe large Λ 1 for c α ∼ 0. Since the γγ resonance is a very powerful channel to probe most of the parameter space of our models, it is expected to observe these scalars in γγ resonance first before observing them in other channels.

VI. SUMMARY AND CONCLUSIONS
In this paper, we have investigated the phenomenology of heavy scalars associated with the variation of the gauge couplings of the SM. In our previous work [7], we studied for the first time the collider phenomenology of a theory with a varying fine-structure constant (or in other words, a varying EM coupling e) originally proposed by Bekenstein [5]. We introduced a TeV-scale mass of the scalar associated with the variation of e, which therefore can be searched for at the LHC. In this work, we extend this idea and explore for the first time the phenomenology of the three heavy scalars associated with the variation of the three gauge couplings of the SM. Since these scalars are heavy, one has to go to high enough energy to excite them. In the electroweak sector of the SM, the EM coupling e appears only after EWSB. Therefore, it is more natural to consider the simultaneous variations of the two gauge couplings g 1 and g 2 of the EW sector than to consider the variation of e only.
Above, we first discussed briefly the construction of the simplest and most economical model of a varying e, with only two free parameters M φ and Λ. Following the same method, we then presented the construction of the model with varying g 1 and g 2 . We discussed two different models in this context. The one scalar theory has a single scalar S 0 , which is responsible for the variation of both g 1 and g 2 . This model is very economical with just two free parameters, M S0 and Λ 0 , and hence is very predictive in terms of collider signatures. The two scalar theory, on the other hand, has two different scalar fields, S 1 and S 2 , which are responsible for the variations of g 1 and g 2 , respectively. To make our analysis general, we consider the mixing of these two scalars among themselves and, therefore, introduce a mixing angle α as one of the free parameters of the theory. This model is less predictive than the one scalar theory as it contains in total five free parameters. However, this general setup is theoretically more pleasing and phenomenologically very rich, and one would expect many interesting collider signatures, which we have discussed in detail. In addition, we also consider the model with varying strong coupling of the SM and the phenomenology of the associated scalar S 3 . In order not to make our analysis more complex, we do not consider any mixing of S 3 with the other scalars.  The scalars S 1 and S 2 can be produced in γγ, γZ, ZZ, or W W fusion at the LHC and can decay to γγ, γZ, ZZ, or W W . This can lead to various interesting signatures that we systematically discussed above. Depending on the free parameters of the theory, the branching ratios and production cross sections of these scalars show interesting patterns. We discussed these aspects by considering five benchmark scenarios, including two extreme limits which we called the photophilic and photophobic limits. In the photophilic limit, the scalars have a strong affinity to couple to photons. In the photophobic limit, on the other hand, the scalars mostly couple to W bosons. The strengths of their couplings also depend on the mixing angle α. We present the branching ratios and production cross sections of these scalars as functions of cos α.
We have derived exclusion limits of the model parameters by making use of the latest LHC γγ, γZ, ZZ, W W , and jj resonance search results. We found that the γγ search in particular is very powerful in constraining most of parameter region of the varying electroweak theory. For example, the γγ search sets lower limit on the new physics scale roughly around 20 TeV for a scalar mass around 1 TeV. Using dijet resonance search data, the scale Λ 3 can be excluded up to 15 TeV for the M S3 ∼ 2 TeV. We also presented a simplified estimation of the sensitivity to our model parameter space at high luminosity LHC.
In this paper, we limit ourselves to the varying gauge couplings of the SM. As discussed in the Introduction, this paper is related to several theoretical developments that may now be investigated with respect to observable signatures. In the same spirit, one can extend our work by promoting other dimensionless constants of the SM to scalar fields. For example, one may consider the variation of Yukawa couplings and the Higgs quartic coupling. The multiple scalars present in our framework might provide a strong first order electroweak phase transition, as required for baryogenesis. A quintessence potential, as introduced in [36] to address the vacuum stability of the Higgs potential, is of a similar structure to the potential and interactions that arise in our framework if we let the Higgs quartic coupling vary and promote the variation to a scalar field. These examples, and the first study presented in this paper, illustrate the future observational potential of the fundamental theoretical issue of whether what we now consider to be constants of nature have a deeper physical origin, as suggested in the considered extensions of the Standard Model.  In this appendix, we will for completeness derive the Lagrangian for a varying EM coupling shown in Sec. II, including the Bekenstein scalar and the constant coupling e. This is generalized to the other gauge groups in Sec. II A. We define the regular QED field strength tensor without the scalar field as F µν = ∂ µ A ν − ∂ ν A µ , while the varying field strength tensor is given by F µν = ε −1 [∂ µ (εA ν ) − ∂ ν (εA µ )]. We will be assuming small fields, and write ε 1 + ϕ, where ϕ is defined in Sec. II. In the following, we keep only terms linear in ϕ or its derivatives.
We begin with the coupling of ϕ to photons. Consider the kinetic term for the EM field, Expanding and keeping terms to lowest order in ε and ϕ give Because we start with a gauge invariant theory, the above Lagrangian is also gauge invariant. However, neither the kinetic term −(1/4)F µν F µν nor the interaction term ∂ µ ϕA ν F µν is gauge invariant by themselves. Indeed, the field-strength tensor transforms under the gauge transformation (3) as F µν → F µν − ε −2 [∂ µ ε ∂ ν α − ∂ ν ε ∂ µ α] , which for constant ε = 1 reduces to the QED result. The corresponding transformations for the kinetic and interaction terms, ignoring all higher order terms in ε, are given by and the sum of these terms, and therefore the Lagrangian, is gauge invariant, but not the individual terms. This is somewhat unusual, but shows that the interaction of the photon with the scalar is intimately connected to the kinetic term. We now switch from the dimensionless field ϕ to the canonically normalized field φ = Λϕ. The new scalar φ also couples to charged particles. The electric charge times photon field eA µ is to be everywhere replaced by eεA µ , so terms coupling φ to charged particles will be generated from the covariant derivatives for charged fermions.
The covariant derivative acting on charged fermions is given by D µ = ∂ µ − i eQA µ , where Q is the charge of the corresponding fermion. Using e = eε gives If we now define D µ = ∂ µ − ieQA µ as the more familiar covariant derivative of electromagnetism, the gauge invariant kinetic term for fermions is given by where the last term is a four-point interaction term for the vertex φγψψ. This is the only direct coupling of the scalar field to fermions. We finally allow the scalar φ to be massive. We then have the Lagrangian where ψ is a generic fermion field with charge Q. Note that heuristically, the new interaction vertices of this model are obtained by inserting a scalar into each possible QED vertex and propagator.