Searching for LFV Flavon decays at hadron colliders

The search for Flavons with a mass of $\mathcal{O}$(1) TeV at current and future colliders might probe low-scale flavor models. We are interested in the simplest model that invokes the Froggatt-Nielsen (FN) mechanism with an Abelian flavor symmetry, which includes a Higgs doublet and a FN complex singlet. Assuming a CP conserving scalar potential, there are a $CP$-even $H_F$ and a $CP$-odd $A_F$ Flavons with lepton flavor violating (LFV) couplings. The former can mix with the standard-model-like Higgs boson, thereby inducing tree-level LFV Higgs interactions that may be at the reach of the LHC. We study the constraints on the parameter space of the model from low-energy LFV processes, which are then used to evaluate the Flavon decay widths and the $gg\to \phi \to \tau\mu$ ($\phi=H_F,\,A_F$) production cross section at hadron colliders. After imposing several kinematic cuts to reduce the SM main background, we find that for $m_{H_F}$ about 200-350 GeV, the decay $H_F \to \tau \mu$ might be at the reach of the LHC for a luminosity in the range 1-3 ab$^{-1}$, however, a luminosity of the order of 10 ab$^{-1}$ would be required to detect the $A_F \to \tau \mu$ decay. On the other hand a future 100 TeV $pp$ collider could probe masses as high as $\mathcal{O}$(10) TeV if it reaches an integrated luminosity of at least 20 ab$^{-1}$. Therefore, the 100 TeV Collider could work as a Flavon factory.


I. INTRODUCTION
After the discovery of a Higgs-like particle with a mass m h = 125 − 126 GeV [1,2], the search for new physics (NP) has become the one of the next goals of the LHC. Although current measurements of the spin, parity, and couplings of the Higgs boson seem consistent with the standard model (SM) [3], its light mass seems troublesome, i.e. the hierarchy problem, and calls for new physics (NP). The SM has also other open issues, such as the flavor problem, unification, etc. [4,5], which also encourages the study of NP models.
The couplings of the Higgs particle to a pair of massive gauge bosons or fermions have strengths proportional to the masses of such particles. However, the LHC has tested only a few of such Higgs couplings, namely the ones to the gauge bosons and the heaviest fermions. Along these lines, many studies have been devoted to analyze the pattern of Higgs couplings derived from LHC data, for instance [6,7]. However, non-standard Higgs couplings, including the flavor violating (FV) ones, are predicted in many models of physics beyond the SM [8][9][10][11]. In particular, the observation of neutrino oscillations, which is associated with massive neutrinos, motivates the occurrence of lepton flavor violation (LFV) in nature [12]. Within the SM, LFV processes vanish at any order of perturbation theory, which motivates the study of SM extensions that predict sizeable LFV effects that could be at the reach of detection. Apart from decays such as l i → l j γ and l i → l jlk l k , particularly interesting is the decay h → τ µ, which was studied first in Refs. [13,14], with subsequent analyses on the detectability of the signal appearing soon after [15,16]. This motivated a plethora of calculations in the framework of several SM extensions, such as theories with massive neutrinos, supersymmetric theories, etc. [11,[17][18][19][20][21][22]. After the Higgs boson discovery, the decay h → τ µ offers a great opportunity to search for NP at the LHC. Although a slight excess of the h → τ µ branching ratio was reported at the LHC run I, with a significance of 2.4 standard deviations [23], a subsequent study [24] ruled out such an excess and put the limit BR(h →μτ ) < 1.2 × 10 −2 with 95% C.L. In the model we are interested in, LFV effects are induced at the tree-level in the scalar sector, so it is thus worth assessing their phenomenology. Another open issue in the SM is the flavor problem [25], which has long been the focus of interest, with several proposals meant to address it, such as textures, GUT-inspired relations, symmetries, radiative generation, etc. In particular, a flavor symmetry approach can be supplemented with the Froggatt-Nielsen (FN) mechanism, which assumes that above some scale M F , there is a symmetry (perhaps of Abelian type U (1) F ) that forbids the appearance of Yukawa couplings; SM fermions are charged under this symmetry. However, the Yukawa matrices can arise through non-renormalizable operators. The Higgs spectrum of these models could include a light Flavon H F , which could mix with the Higgs bosons when the flavor scale is of the order of the TeVs. Quite recently, the phenomenology of Higgs-Flavons at particle colliders has been the focus of attention [26][27][28][29][30].
Although the states we are interested in arise from the mixing of Higgs bosons and Flavons, we will still call them Flavons for short. Depending on the particular model, there could be several potentially detectable Flavon decays, which would be indistinguishable from the decays of a heavy Higgs boson, therefore, to search for a distinctive signature, we will focus on the one arising from the LFV decay H F → τ µ, with τ µ = τ − µ + + τ + µ − . The minimal model that introduces the FN mechanism with an Abelian flavor symmetry includes a scalar sector consisting of a Higgs doublet and a FN complex singlet. From now on we will refer to this model as the FN extension of the Standard Model (FNSM). Such a model predicts a CP -even Flavon H F and a CP -odd one A F . Also, the couplings of the light SM-like Higgs boson would deviate from the SM ones, such that two possible scenarios are possible: firstly, the mixing of the real part of the doublet with the real component of the FN singlet could induce sizable LFV Higgs couplings of the light physical Higgs boson, which might affect the light Higgs phenomenology; secondly, the CP -even Flavon could be very heavy, so its mixing with the light Higgs boson would be negligible and unconstrained by LHC Higgs data, though in such a case the CP -odd state A F could be the lighter one, thereby giving rise to a potentially detectable LFV signal.
In this paper we are interested in studying the possible detection of both CP -even and CP -odd Flavons at the LHC and a future 100 TeV pp collider via their LFV decays. It has been pointed out that a 100 TeV pp collider would allow for a detailed study of several topics of interest in particle physics, such as Higgs physics and the electroweak symmetry breaking mechanism [31]. It will also be useful to search for possible signals of dark matter, SUSY theories, and other extension models [32].
We will start our analysis by considering the constraints on the parameter space of the FNSM obtained from Higgs data at the LHC, low energy LFV, and the muon magnetic dipole moment. A set of benchmarks will then be used to estimate the Flavon decay modes, focusing on the LFV ones, as well as their production cross sections by gluon fusion at the LHC and a future 100 TeV pp collider. We will then explore the possibility that the CP -even and CP -odd Flavons could be detected via the τ µ decay channel, for which we will make a Monte Carlo analysis of the signal and the SM main background.
The organization of our work is as follows. In Sec. II we describe the realization of the FN mechanism within the simplest model, namely that with one Higgs doublet and one FN singlet. In particular, we present the Higgs potential and Yukawa Lagrangian, from which the Flavon couplings can be extracted. Section III is devoted to the constraints on the parameters space of the model and the benchmarks of parameter values we will be using in our analysis, whereas Sec. IV is focused on the analysis of the decay modes of both a CP -even and a CP -odd Flavon as well as their production cross sections via gluon fusion at the LHC and a future 100 TeV pp collider. We also present the Monte Carlo analysis of the µτ signal and its main background. The conclusions and outlook are presented in Sec. V.

II. THE SCALAR SECTOR OF THE MINIMAL FNSM
The scalar sector of the FNSM includes the usual SM Higgs doublet and a complex singlet where v denotes the SM vacuum expectation value (VEV) and u that of the FN singlet, whereas G + and G 0 will be identified with the pseudo-goldstone bosons that become the longitudinal modes of the W and Z gauge bosons.

A. The Higgs potential
We turn now to discuss the minimal CP -conserving Higgs potential with a softly-broken U (1) global symmetry, which is given as follows We are therefore left with the following U (1)-symmetric terms (m 2 1 , m 2 s1 , λ 1 , λ s , λ 11 ), and the U (1)-soft-breaking term m 2 s2 . The latter is required to avoid a massless Goldstone boson when S F = 0. An extensive analysis of this potential was presented in reference [33], where the parameter space that allows a viable model was identified.
After imposing the minimization conditions on the potential, the following relations are obtained:

B. The Scalar Mass Matrix
In a CP -invariant potential, the CP -even (real) and CP -odd (imaginary) components of the mass matrix do not mix. In this case the mass matrix for the real components in the basis (φ 0 , s) is given by: whereas the mass matrix for the imaginary components, in the basis (G Z , p), reads We notice that the mass scale for the CP -odd state arising from the FN singlet A F = p is different from the VEV u, which is the U (1)-breaking scale, and therefore it could be much lighter. As for the mixing of the real components of the doublet Φ and the singlet s, the mass eigenstates are obtained through the standard 2 × 2 rotation: In what follows we will identify the mass eigenstate h as the SM-like Higgs boson with m h = 125 GeV, while the mass eigenstates H F and A F will be assumed to be heavier. Although they arise from Flavon-Higgs mixing, in the present work we will still refer to H F and A F as Flavons for short. The properties of the CP -even Flavon will depend on the size of its mixing with the lightest state. On the other hand, the CP -odd state, which does not couple to gauge bosons, it does couple to the SM fermions, including both diagonal and non-diagonal interactions.
Our analysis of Flavon decays requires the knowledge of cubic interactions, such as the trilinear vertex H F hh, which is given in the minimal model by:

C. Yukawa sector and LFV interactions
The FN Lagrangian of the model includes the terms that become the Yukawa couplings once the U (1) F flavor symmetry is spontaneously broken. It is given by: where q f ij (f = u, d, l) denote the Abelian charges that reproduce the observed fermion masses, for each fermion type. The Flavon field S F is assumed to have flavor charge equal to -1, such that L Y is U (1) F -invariant. Then, the Yukawa couplings arise after the spontaneous breaking of the flavor symmetry, i.e. λ x = ( <S F > Λ F ) nx , where < S F > denotes the Flavon VEV, whereas Λ F denotes a heavy mass scale, which represents the mass of heavy fields that transmit such symmetry breaking to the quarks and leptons. For a detailed discussion of the viable structures for the Yukawa matrices for the FN multi-Higgs model see Ref. [27]. Here we shall discuss the generic features that can be identified by studying specific Abelian charges for the charged leptons. With this purpose we consider the charge assignment used by one of us in Ref. [30], where the Yukawa matrix Y l is of the form: which can be justified with the following Abelian charges Then, Y l ij ∼ λ |αi−βj | , which means that choosing Y l 23 ∼ λ 4 and Y l 33 ∼ λ 2 implies |α 2 − β 3 | = 4 and |α 3 − β 3 | = 2. Other examples of Abelian charges are presented in Ref. [27,34,35]. We then focus on the 2-3 submatrix, which dictates the φτ µ interaction. The corresponding squared mass matrix can be written as: with Y l 33 ∼ λ 2 = (m τ − m µ )/v. This mass matrix can be diagonalized by a 2 × 2 rotation with mixing angle θ ∼ (m µ /m τ ).
In the unitary gauge we set G ± → 0, G Z → 0, thus we can write the doublet Φ as follows: whereas the powers of the Flavon field can be expanded as Then, after substituting the Higgs boson and Flavon mass eigenstates, one gets finally the following interaction Lagrangian for the Higgs-fermion couplings where we use the usual short-hand notation s α ≡ sin α and c α ≡ cos α.
Here,M f is the diagonal mass matrix, whereas the information about the size of FV Higgs-Flavon couplings is contained in given in the flavor basis as which remains non-diagonal once the fermion mass matrices are diagonalized, thereby giving rise to FV scalar couplings. For the Yukawa matrix we are considering, the Z matrix is given as follows (for the 2-3 subsystem): Thus we can use the approximate diagonalization of the mass matrix to express Y l ij and Z l in terms of the mass ratios. To leading order we have and similary for up-type quarks, whose Higgs couplings play a fundamental role for collider phenomenology. With the assignment of Abelian charges chosen in [30], the Z u matrix for the 2-family case is given by: However, in the mass-eigenstates basis, the rotated matrix takes the following form to leading order: where s u = sin θ u , with θ u the mixing angle. It turns out that the φtt coupling vanishes in the weak basis, but it is non-vanishing as long as there is mixing between the top quark and the mass eigenstates basis. The diagonal and non-diagonal interactions of the h, H F , and A F scalar bosons with massive fermions are given by: where the Feynman rule for the A F f i f j vertex includes a γ 5 Dirac matrix and r s = v/( √ 2u). Besides the Yukawa couplings, we also need to specify the scalar-to-gauge-boson couplings. They can be readily extracted from the kinetic terms of the Higgs doublet and the singlet, which transforms trivially under the SM gauge group. Thus after substituting Eq. (8) in the kinetic term, we obtain that the h and H F couplings to gauge boson pairs are SM-like, with the coupling constants given by Thus the coupling constants are We are interested in the possible detection of CP -even and CP -odd Flavons with masses of the TeV order at both the LHC and a future pp collider with a center-of-mass energy of 100 TeV. Depending on the particular model, there could be several potentially detectable decays of such Flavons, but some of them would also arise from heavy Higgs bosons, for instance within multi-Higgs doublet models. Thus, in order to search for a distinctive Flavon signature, we shall focus on the one arising from the LFV decay φ → τ µ (φ = H F , A F ). In order to determine the detectability of this decay, we will present a Monte Carlo analysis of the signal and the most relevant SM backgrounds.

III. CONSTRAINTS ON THE PARAMETER SPACE OF THE FNSM
In order to evaluate the Flavon decays and production modes at a hadron collider we need to analyze the most up-to-date constraints on the model parameters. For the mixing angle α we can use the data obtained by the LHC collaborations on the Higgs boson properties, whereas the LFV couplings can be constrained via the experimental data on the muon anomalous magnetic dipole moment (AMDM) a µ , the LFV decays of the tau lepton τ →l ilj l j and τ → l i γ, as well as the experimental constraint on the h → µτ decay. All the necessary formulas to perform our analysis below are presented in Appendix A.

A. Mixing angle α
We shall use the universal Higgs fit of Ref. [7], which presents constraints on the parameters X , defined as (small) deviations of the Higgs couplings from the SM values, i.e. χ h X = g hXX /g SM hXX = 1 + X . For the W and Z gauge bosons, the corresponding constraints are W = −0.15 ± 0.14 and Z = −0.01 ± 0. 13. Regarding the fermion couplings, we notice that this universal Higgs fit is valid for the CP -conserving case, which we are considering here; therefore we can apply them to constrain the properties of the CP -even SM-like Higgs boson. Furthermore, the constraints derived from the gauge interactions provide the strongest constraints on the mixing angle α. Since the lightest scalar boson h couples with the SM gauge bosons with a strength that deviates from the SM couplings by the factor c α , to satisfy the bound on Z we need to have 0.86 < c α < 1. We will use a conservative approach and use the benchmark c α = 0.95 in our analysis below. B. DiagonalZ l 22 andZ l 33 matrix elements In this work we will use the 2-family approximation, neglecting FV with the leptons of the first fermion family. We will also assume that there is no CP -violating phase, which means that we have three free parameters:Z 22 ,Z 23 , andZ 33 . To constrain the diagonalZ 33 matrix element, we refer again to the universal Higgs fit of Ref. [7] and consider the constraint on the deviation of the SM hτ τ coupling, namely, τ = 0 ± 0. 18. We thus show in Fig. 1 the allowed area on the u −Z l 33 plane for two values of c α . We observe that in order to agree with the universal Higgs fit constraint when c α = 0.9,Z l 33 must be of the order of 10 −3 for u = 0.5 TeV and 10 −2 for u = 2 TeV, but when c α = 0.95 we must have values of the order 10 −2 in the complete u interval. Since we are considering c α = 0.95, we will useZ l 33 = 10 −2 as benchmark.  We will consider the current experimental bounds on the muon AMDM a µ [39], the tau decay τ → µγ [39], and the Higgs boson decay h → τ µ [23,40] to constrain the Z l ij matrix elements [41]. Notice that the current bound on the τ → 3µ decay width gives very weak constraints, so we will omit such a process in our analysis. Two scenarios arise when dealing with constraints on LFV couplings: • Scenario I: the FNSM is assumed to be responsible for the current discrepancy between the theoretical and experimental values of a µ , which requires a positive contribution from new physics. Along this line, the oneloop contribution from a CP -even scalar boson is positive, whereas that of the CP -odd scalar boson is negative. Therefore a suitable scenario would be that with a relatively light CP -even Flavon and a heavy CP -odd Flavon, which would suppress the negative contribution. When m H F m A F , the Flavon contributions would largely cancel one another out for c α 1 and the remaining (positive) contribution would arise from the corrections to the SM Higgs boson couplings, unless there were large extra positive contributions arising from the CP odd Flavon at the two-loop level: such contributions could arise as long as A F is very light and there was an enhanced one-loop-induced A F γγ coupling, as discussed in [42]. We do not expect such an enhancement in the FNSM, so large values of the mixing matrix element Z l 23 would be required to solve the a µ discrepancy. • Scenario II: the FNSM Flavons give a negligible contribution to a µ and fail to explain the discrepancy, though the model still can remain viable. It could happen for instance that extra positive contributions to a µ would arise from the ultraviolet completion of the FNSM, thereby solving the puzzle. Also, it is still possible that more precise determinations of the SM hadronic contribution and the experimental measurement would settle the discrepancy in the future without requiring any NP effects. In any case, by requiring that the Flavon contribution to a µ is negligible, it is enough to satisfy the LHC constrain from h → τ µ in order to have a viable parameter space in this scenario.
We will now assess the viability of both scenarios. To avoid large corrections to the diagonal lepton scalar couplings, we takeZ 22 10 −3 andZ 33 10 −2 . We then show in Fig. 2   constraints on the τ → µγ and h → τ µ decays for the indicated values of the Flavon masses and the mixing angle c α . The blue strip is the region where u andZ 23 must lie in order to alleviate the a µ discrepancy. We observe that the most stringent limits onZ 23 are obtained from the experimental bound on the h → µτ decay, which requiresZ 23 to be of the order of 10 −2 (10 −1 ) for u = 0.5 (1) TeV. However, in order to solve the a µ discrepancy, values ofZ 23 as large as 1 would be required for u = 1 TeV, which is due to the fact that we need large (positive) contributions arising from the CP -even Flavon and the new corrections to the SM Higgs boson couplings. We note that the blue strips shifts downwards as the mass of the CP -odd Flavon increases as in this case the (negative) contribution is smaller. Therefore, scenario I discussed above is not favored by the current experimental data and we assume that the scenario II is fulfilled, with the Flavon contribution to a µ being rather small and not responsible for the discrepancy between experimental and theoretical values. We thus take as benchmark Z 23 0.1, which is allowed by the experimental bound on the τ → µγ and h → µτ decays for u 1.5 TeV. According to the universal Higgs fit [7], the allowed value for the deviation of the SM htt coupling is t = −0.21±0. 22. We will use again the two-family approximation and take the valuesZ u 33 = 0.01 andZ u 23 = 0.1.

E. Summary of Benchmarks for the model parameters
In summary, in our study below we will use the following benchmarks: 1. Mixing angle α: As discussed above, to satisfy the fit on the 125 GeV Higgs couplings measured at the LHC, the following constraint must be obeyed 0.86 < c α < 1. We will thus use the benchmark c α = 0.95. 2. FN singlet VEV u: It appears in the Flavon couplings and the LFV SM Higgs coupling. We will consider the value 1.5 TeV.
3.Z f ij matrix: It determines the strength of the LFV scalar couplings. We will use the 2-family approximation and takeZ l 22 = 10 −4 ,Z l 33 = 10 −2 , andZ l 23 = 10 −1 , which are consistent with the constraint on the LFV decay h → µτ . 4. H F hh interaction: This vertex depends on a combination of parameters that appear in the Higgs potential.
However, these parameters could be traded by an effective coupling λ H F hh , which can take values of the order of O(1). We will thus fix λ H F hh 0.1u.
A summary of the benchmarks we are going to consider in our analysis is presented in Table I.

IV. SEARCH FOR LFV FLAVON DECAYS AT HADRON COLLIDERS
As stated above, the aim of this work is to analyze the detectability of the LFV signal arising from the Flavon decays, as predicted by the FNSM, at the LHC and a future 100 TeV pp collider. Below we will present an analysis   concentrating on the main Flavon production mechanism, i.e. gluon fusion, as well as the branching ratios of its dominant decay modes. Then, we will present the Monte Carlo analysis of the H F → τ µ and A F → τ µ decay signatures, including the study of the potential SM background. We will present a conservative analysis meant to find out whether it is possible to have evidence of our signal at the LHC and the future 100 TeV pp collider.

A. Production cross-sections of the CP -even and CP -odd Flavons
We now turn to analyze the main production mode of both H F and A F Flavons at hadronic colliders, namely, by gluon fusion. The High-Luminosity Large Hadron Collider project aims to increase potential discoveries contemplating a luminosity of up to L=3 ab −1 about 2025. As far as a 100 TeV pp collider is concerned the future circular collider (FCC) contemplates an integrated luminosity of until O(10 ab −1 ). We consider the integrated luminosities shown in tha table II. In Fig. 3 we show the pp → H F X production cross section of a CP -even Flavon as a function of its mass m H F at the LHC and a future 100 TeV pp collider. We also show the event numbers on the right axis of each plot. As far as the CP -odd Flavon is concerned, the respective production cross section and event numbers are presented in Fig. 4.  The dominant contribution to gluon fusion arises from loops carrying the top quark as shown in Fig 5. This explains the suppression of the production of the CP -odd Flavon as compared to that of the CP -even one, as observed in Figs. 3 and 4, which stems from the appearance of the coupling g A Ft t in the corresponding cross section. For instance, taking into account the parameter values of Table I, we have g 2 A Ft t ∼ 10 −6 , whereas for the CP -even Flavon g 2 H F tt ∼ 10 −2 .

Two-body decays
We now analyze the behavior of the branching ratios of the dominant decay modes, including the FV ones, of both the CP -even and CP -odd Flavons. Analytical expressions for the partial decay widths are presented in Appendix B. It is worth mentioning that a crosscheck was done by comparing the numerical results obtained via our own C language code implementing the analytic expressions of Appendix B, and those computed with the aid of the CalcHEP package [43], for which we used an implementation of the FNSM Feynman rules obtained with LanHEP [44]. In Fig.  6 we present the relevant branching ratios of the decays modes of a CP -even Flavon as functions of its mass for the benchmarks of Table I.  Table I.
We observe in Fig. 6  As far as the CP -odd Flavon A F is concerned, since it does not couple to gauge bosons at the tree-level, its main decay modes are into fermion pairs. The corresponding branching ratios are shown in Fig. 7. We can conclude that the decay modes A F →f f could have branching ratios up to two orders of magnitude larger than the analogue branching ratios of the CP -even Flavon . For the parameter values used here, the main decay channels of the CP -odd Flavon would be A F →tc, A F →bb, A F → τ µ, whereas the A F →tt decay would be very suppressed due to the A Ft t coupling. As for the one-loop decays A F → γγ and A F → γZ, they have very small branching ratios, of the order of 10 −8 . Along this line, in Fig. 7 we also have included the branching ratio for the one-loop induced decay A F → hZ, which proceeds via the fermion loops shown in Fig. 8. We have calculated the corresponding decay width, which is shown in Appendix B. It has pointed out recently [45] that this decay could have a relevant branching ratio in models with a CP -odd scalar boson arising from a complex singlet and it could be useful to look for evidences of CP violation. We observe in Fig. 7 that this decay has a somewhat suppressed branching ratio in the scenario of the FNSM we are interested in: BR(A F → hZ) 10 −4 for m A F in the interval between 200 GeV and 1000 GeV. This value is much larger than the branching ratios of the loop induced decays A F → γZ and A F → γγ but it is about the same size than BR(A F → gg). A rough estimate shows that other loop induced decays such as A F → ZZ and A F → W W have also suppressed branching ratios. This stems from the fact that not only they are of higher order in the coupling constants but are suppressed by the loop factor 1/(16π 2 ), which enters squared into the decay width. In general, the one-loop induced branching ratios BR(A F → V V ) (V = W, Z) and BR(A F → hZ) are suppressed with respect to those of the tree-level induced decays A F →f f by a factor of (3g 2 m t /(16π 2 m W )) 2 3.3 × 10 −4 , assuming that the main contribution arises from the top quark. Although these one-loop induced decays may be interesting by themselves, they have no impact on our study of the lepton flavor decays of the CP -odd Flavon.

Three-body decays HF →fifjh
At the tree-level the CP -even Flavon also couples with a fermion pair and a SM Higgs boson, thus it is worth analyzing the behavior of the three-body decay modes H F →f i f j h. The corresponding decay width is presented in Appendix B. In Fig. 9 we show the branching ratios for these three-body decays. We observe that, for a relatively light Flavon with mass around 300 GeV, the decay channel H F →bbh would have a branching ratio as large as 10 −1 . On the other hand, other kinematically allowed H F →f i f j h decays would reach branching ratios as high as 10 −2 . For a heavier Flavon with m H F > 600 GeV, the decay H F →tth would become open and could be the dominant three-body decay mode for m H F 900 GeV. Although these decay channels seem worth a more detailed study, we will content ourselves with obtaining the event numbers that could be achieved at the LHC and the FCC. We consider values for the luminosities of the table II. In the Fig. 10 we shown the number of events for the processes pp → H F → qqh with q = b, t.    √ s = 100 TeV for an integrated luminosity from L=3 ab −1 (lower limit) to 30 ab −1 (upper limit). The parameters of Table I were used. C. Search for LFV Flavon decays at the LHC and a future 100 TeV pp collider We are interested in the possible detection of the CP -even and CP -odd Flavons via their LFV decay into a τ µ pair at the LHC and the FCC. We thus show in Fig. 11 the event numbers for the processes pp → H F → τ µ and pp → A F → τ µ, for √ s = 14 (100) TeV and an integrated luminosity for the values displayed in the Table II. We also use the same set of parameter values of Table I. We note that for a CP -even Flavon with a mass about 1 TeV, there would be about O(10) µτ signal events at the LHC and O(10 4 ) events at the FCC. These event numbers would decrease by about two orders of magnitude for a CP -odd Flavon. We will now analyze the signature of the LFV Flavon decays H F → τ µ and A F → τ µ, with τ µ = τ − µ + +τ + µ − , and their potential SM background. We are inspired in the analysis carried out by the CMS collaboration in Refs. [23], ATLAS collaboration [46] and the work of the authors of [47]. The ATLAS and CMS collaborations considered the two following tau decay channels: electron decay τ e µ and the hadron decay τ h µ. For our analysis we will concentrate instead on the electron decay. As far as our computation scheme is concerned, we first use the LanHEP routines to obtain the FNSM Feynman rules for Madgraph [48]. In this way, the signal and background events can be generated by MadGraph5 and MadGraph5 aMC@NLO, respectively, interfaced with Pythia 6 [49] and Delphes 3 [50] for detector simulations. The background events were generated at NLO in QCD and the signal at LO, using the CT10 parton distribution functions [51]. The signal and main background events are as follows: • Signal: the signal is gg → φ → τ µ → eν e ν τ µ with φ = H F , A F . The electron channel must contain exactly two opposite-charge leptons, one an electron and the other a muon. Then, we search for the final state eµ+miss energy. We consider the specific case for which an integrated luminosity L of 1 ab −1 for the LHC is considered and in the interval 1 − 20 ab −1 for the FCC.
• Background: the main SM background arises from Z production via the Drell-Yan process, followed by the decay Z → τ τ as well as W + W − and ZZ pair production and jets. In this work we will only consider the main background to assess how our signal could be searched for.  Table III, where we also show the event numbers of the signal (S) and background (B) after the kinematic cuts are applied, along with the signal significance S/ √ S + B for m H F = 200 GeV. The effect of the cuts on the signal and background event numbers is best illustrated in Fig. 12, where we show how the efficiencies signal and background evolve after each cut is successively applied. One can observe that once the kinematic cuts are applied, the resulting signal efficiency is about 0.31, whereas that for the background is around 0.02. With a luminosity of 0.3 ab −1 , the signal significance is about 2σ. However, if we take into account an integrated luminosity of 3 ab −1 it increases up to ∼6.5σ. The net effect is shown in Fig. 13, where the signal significance is plotted as a function of the luminosity for the chosen values of m H F .  Table  III are successively applied. We use mH F = 200 GeV and the parameters of Table I.
Therefore, it seems troublesome that a CP -even Flavon with a mass greater than 300 GeV could be detected at the LHC via the H F → τ µ decay channel. In such a case, we can turn to the decay modes H F → W W and H F → ZZ, which are the dominant ones. These decay channels seem more promising for the detection of a heavy Flavon. Along this line, to assess the potentiality of the W W and ZZ decay channels for the Flavon detection, we show in Fig. 14 the corresponding event numbers that could be produced at the LHC and a 100 TeV pp collider. We note that about 10 5 W W events would be produced at the LHC for m H F = 1 TeV, whereas the number of ZZ events would be slightly smaller. This seems more promising for the signal detection, though a more detailed analysis of the background would be required to draw a definitive conclusion. As for the CP -odd Flavon, we notice that it has a smaller pp → A F → τ µ production rate, therefore, in order to have evidence of the A F → τ µ decay, higher luminosities, of the order of 10 ab −1 , would be required. This seems inaccessible for the LHC, however, we expect that the search for this decay could be possible at the FCC.

Analysis at the future 100 TeV pp collider
The building of a 100 TeV pp collider is under consideration [53]. The luminosity and center-of-mass energy are crucial factors to allow detectable LFV signatures of the FNSM. Luminosity goals for a 100 TeV pp collider are discussed by the authors of Ref. [53]. Again, we consider the values of the Table II, although up to 30 ab −1 might be reached. In our analysis we consider the conservative kinematic cuts in order to give a general overview assess how our signal could be searched for. The applied cuts to both the signal and background are shown in Table IV, whereas in Fig. 15 we show the signal significance as a function of the luminosity and the Flavon mass m φ , with φ = H F , A F . We observe that at a 100 TeV pp collider with integrated luminosity of 10 ab −1 , it would be possible to probe CP -even Flavon masses in the multi-TeV range, up to about 4 TeV. We also notice that CP -odd Flavon masses could be searched until the order ∼1 TeV.

V. CONCLUSIONS
In this work we have explored the possibility that the LFV τ µ decay channel of CP -even H F and CP -odd A F Higgs Flavons with a mass of a few hundreds of GeVs could be at the reach of detection at the LHC and a future 100 TeV pp collider, which would serve as a possible probe of low-scale flavor models. For the theoretical framework, we have considered the simplest Froggatt-Nielsen model, with an Abelian flavor symmetry and a CP -conserving Higgs sector that includes a Higgs doublet and a Froggatt-Nielsen complex singlet. In this model the CP -even Flavon can mix with the SM-like Higgs boson, thereby inducing tree-level LFV interactions mediated by the latter. In this work we concentrate instead on the LFV couplings of both the CP -even and CP -odd Higgs Flavons. After studying the constraints on the parameter space of the model from low-energy LFV processes, we choose a set of benchmarks and estimate the relevant decay modes and the production cross section of the Flavons via gluon fusion at the LHC and a future 100 TeV pp collider. We then consider a set of kinematic cuts for both the signal and the SM main background. It is found that the LHC has the potential to discover the LFV decay H F → τ µ for m H F between 200 and 350 GeVs provided that luminosities, of the order of 1-3 ab −1 , are achieved. In such a case other decay channels would be more appropriate to search for the signal of a Flavon at the LHC. As far as a future 100 TeV pp collider is concerned, it would be able to probe the LFV τ µ decay channel for Flavon masses as heavy as 10 TeVs, as long as an integrated luminosity of at least 20 ab −1 was available, which has been deemed viable in the literature regarding the possible construction of such a collider. Therefore, besides other physics goals, a 100 TeV Collider might also work as a Flavon factory.
Appendix A: Flavon contributions to LFV h and τ decays and the muon anomaly In this Appendix we present the analytical expressions necessary to obtain the constraints on the LFV Flavons couplings shown in Fig. 2. Although these results were meant for CP -even and CP -odd scalar bosons, they are also valid for the Flavons. In the FNSM, the LFV decay h → τ µ proceeds at the tree-level. The decay width can be obtained from Eq. (B1) of Appendix B in the m h m τ m µ limit. The result is given by The CMS collaboration reported a bound on the respective branching ratio: BR(h →μτ ) < 1.2 × 10 −2 [24]. As far as the τ → µγ decay is concerned, it arises at the one-loop level and receives contributions of the SM Higgs boson and the Flavons via the Feynman diagram of Fig. 16(a). The respective decay width is where the C S and C P coefficients stand for the contribution of CP -even and CP -odd scalar bosons, respectively, which in the limit of g φτ τ g φµµ g φee (φ = h, H F , A F ) and m τ m µ m e , can be approximated as [54] C Two-loop contributions can be relevant and the respective expressions are reported in [54]. The current experimental limit on the branching ratio is BR(τ → µγ) < 4.4 × 10 −8 [39]. As for the τ → µμµ decay, it receives contributions from the exchange of a SM-Higgs boson and the Flavons as depicted in the Feynman diagram of Fig. 16(b). The tree-level decay width can be approximated as where S φ = g φµµ g φµτ . It has been pointed in Ref. [54], however, that the one-loop contribution is dominant. We refrain from presenting the corresponding expression as this process, for which the experimental limit on the respective branching ratio is BR(τ → µμµ) < 2.1 × 10 −8 [39], gives very weak constraints on the FNSM parameters. Finally, the muon AMDM also receives contributions from the SM Higgs boson and the Flavons, which are induced by a triangle diagram similar to the diagram of Fig. 16(a) but with two external muons. The corresponding contribution can be approximated for m φ m l as [54] δa µ ∼ m µ where one must take into account the NP corrections to the g hµµ coupling only. If the Flavons are too heavy, the dominant NP contribution would arise from the SM Higgs boson.

The decay AF → hZ
Although the coupling of the Z gauge boson to a pair of CP -even and CP -odd scalar bosons AhZ appears at the tree level in multiple Higgs doublet models, in the FNSM the scalar singlet does not couple to SM gauge bosons, so the A F → hZ decay proceeds up to the one-loop level via fermion exchange through the Feynman diagrams of Fig. 8. It has been pointed out recently that such a decay [45] can be relevant to find evidences of CP violation. We have calculated the corresponding amplitude via the Passarino-Veltman reduction scheme. Using the nomenclature of Fig.  8 we obtain where w a = m 2 a /m 2 A F , whereas the one-loop induced coupling g A F hZ (w h , w Z ) is given as with g f A = 1 2 (− 1