Theoretical constraints on models with vector-like fermions

We illustrate the impact of VLF on the running of the Higgs quartic coupling $\lambda$ and EW vacuum stability given by (10) assuming $n=1$ (as we checked, behaviour of models with larger numbers of VLF families is qualitatively very similar, see also discussion in Section IV.5).

Plots in Fig. 1 show how varying VLF Yukawa couplings influences the running of $\lambda$, the gauge couplings $g_{1}$ and $g_{2}$ and the top Yukawa coupling $y_{t}$. Two main regimes can be distinguished depending on the magnitude of $y$:

They have positive impact on the stability of the EW vacuum, as compared to the SM, eventually even leading to the stability up to the Planck scale. This can be understood at the 1-loop level as follows. For $y=0$, VLF contributions affect directly the running of $g_{1}$ and $g_{2}$ gauge couplings (see Eq. (33)), which in turn leads to bigger values of $\beta_{\lambda}$ (see Eq. (31)). Adding VLQ multiplets leads to larger effects than adding VLL ones due to larger number of degrees of freedom, affecting RGE running through a bigger colour factor $N_{c}^{\prime}$ in Eq. (33). This effect in the case of $g_{1}$ is more model-dependent due to differences in hypercharges of quarks and leptons.

Increasing the VLF Yukawa couplings typically has negative impact on the stability of the EW vacuum. For these couplings bigger than some critical value the EW vacuum becomes unstable. This follows again from the 1-loop RGEs: the leading direct VLF contribution to the running of $\lambda$ is proportional to $(2\lambda-y_{F}^{2})y_{F}^{2}$ (see the first equation in (33)) which is positive for small values of $y_{F}^{2}$ but becomes negative for $y_{F}^{2}>2\lambda$. The negative contribution grows with increasing $y_{F}^{2}$ and very soon overcomes the indirect positive contribution mentioned in point 1 above. This effect is further strengthened by the positive contribution of VLF to $\beta_{y_{t}}$ which indirectly leads to bigger negative term in $\beta_{\lambda}$ proportional to $y_{t}^{4}$ (see the first equation in (31)).

The maximal value of the Yukawa coupling consistent with the stability condition (10) depends not only on the cut-off scale $\Lambda$ but also on the masses of VLF. Any effects from VLF are stronger when their masses are smaller because lighter particles modify RGEs starting from lower energy scales. This applies to both types of effects: those which improve and those which worsen stability of the EW vacuum. In Fig. 2 we illustrate this feature for two values of $y$ in scenario I with VLQ.

As we argued above, VLF with Yukawa couplings larger than some critical values have negative impact on the EW vacuum stability. The bigger are such couplings the lower is the energy scale at which $\lambda$ becomes negative. Thus, demanding $\lambda$ to be positive up to some chosen cut-off scale $\Lambda$ results in an upper limit on the Yukawa coupling, $y_{MAX}$, in a given model²²2We have assumed here that all VLF’s Yukawa couplings are equal. In more general models the corresponding upper bound applies to some combination of different Yukawa couplings. (we remind that by couplings without explicit scale-dependence we denote these couplings renormalized at the scale $\mu=M_{t}$).

Negative contributions to $\beta_{\lambda}$ from the VLF increase with the corresponding Yukawa couplings. Thus, the maximal value of the Yukawa coupling in a given model, $y_{MAX}$, is a decreasing functions of the cut-off scale $\Lambda$. Such dependence for scales above the $10^{9\div 10}$ GeV (close to the scale at which $\lambda$ becomes negative in the SM) is very different for VLL and VLQ models. In the case of VLL models $y_{MAX}$ drops very quickly to zero, while in VLQ models $y_{MAX}$ decreases very slowly. The reason for this difference follows from the impact of VLF with small Yukawa couplings on the RGE evolution of $\lambda$ discussed in the point 1 above (and illustrated in the upper row of plots in Fig. 1). Namely, VLL may improve the EW vacuum stability only sightly and only for very small Yukawa couplings, while addition of VLQ with not too large Yukawa couplings may stabilise the EW vacuum even up to the Planck scale.

A few examples of $y_{MAX}(\Lambda)$ for VLL and VLQ models are shown in Fig. 3. The biggest possible value of $y_{MAX}$, corresponding to the lowest cut-off scale considered in this paper, $\Lambda=100$ TeV, in scenarios I and II for VLL and VLQ models and for two values of $M_{F}$, are collected in Table 2.

Results presented in this Section show that in all considered scenarios the maximal values of VLF Yukawa couplings allowed by the EW stability condition are always relatively small, $y_{MAX}<1$, which, as we shall see later on, suggests limited range of phenomenological consequences of such BSM extensions which could be experimentally tested in near future. Many studies considering phenomenology of VLF models (e.g. [27, 30, 31, 32, 33, 34, 35, 36]) often require large values of VLF Yukawa couplings to produce significant observable effects. The present analysis shows that obtaining those effects may be not possible in perturbative models with stable EW vacuum.

The discussed above impact of VLF on the EW vacuum stability is in line with similar study [21], and we were able to reproduce their results. However, in our work we studied a wider set of VLF model scenarios, and focused more on their possible observable effects on phenomenology (which we discuss in Section V).

Our analysis obviously does not exhaust all possible VLF scenarios, but shows that the simple EW vacuum stability requirement strongly constraints allowed parameter space of such models. Testing VLF models with more complex pattern of different couplings would require case by case study, which is beyond the scope of this work.

In order to explore the freedom of construction of VLF models and enlarge their allowed parameter space, we consider also SM extension containing also the real scalar singlet field, $S$, with a tree-level potential given in Eq. (7). We start with reviewing the theoretical constraints on SM extended by the real scalar singlet only.

The main effect of adding $S$ comes from the coupling between such scalar and the SM Higgs scalar, $\lambda_{HS}$, which gives a positive contribution to $\beta_{\lambda}$ (see the first equation in Eq. (34)). The singlet self-coupling, $\lambda_{S}$, also (indirectly) leads to bigger values of $\lambda$ via its positive contribution to $\beta_{\lambda_{HS}}$ (second equation in Eq. (34)). Thus, increasing scalar couplings, $\lambda_{HS}$ and $\lambda_{S}$, and/or decreasing $M_{S}$ result in larger values of $\lambda$ and so has positive impact on the stability of the EW vacuum. For example one gets its absolute stability up to the Planck scale for $M_{S}=1$ TeV, $\lambda_{S}=0$ and $\lambda_{HS}\gtrsim 0.3$.

However, increasing the values of $\lambda_{S}$ and $\lambda_{HS}$ too much may lead to the loss of perturbativity, i.e. condition(s) (11) and/or (12) may be violated below the chosen cut-off scale $\Lambda$. Fig. 4 shows regions on the $\lambda_{S}$-$\lambda_{HS}$ plane compatible with conditions (11)–(12) for chosen values of $\Lambda$ and the singlet mass $M_{S}$. As expected, increasing value of $M_{S}$ weakens the upper bound on $\lambda_{HS}$ while increasing $\lambda_{S}$ or $\Lambda$ makes this bound stronger.

Our approach provides new and independent way of studying allowed values of the parameters in the real scalar model and can be seen as complementary to the previous papers (e.g. [72, 73, 74]) by further constraining allowed parameter space of this model. It emphasises the importance of going beyond the naive study of perturbativity as given by the condition (11) by taking into account the differences between 1- and 2-loop contributions to the RGEs .

Since, as discussed in the previous subsection, adding the scalar singlet to the model improves stability of the Higgs potential, it should also help to alleviate constraints on VLF Yukawa couplings summarised in Table 2. However, relaxing the constraints on $y_{MAX}$ in the presence of extra singlet scalar field should also take into account potential violation of perturbativity conditions (11) and (12) when increasing the values of new scalar couplings.

In a model with the SM extended by the singlet scalar only (SM+S), all the conditions (10)-(12) are fulfilled in the SM+S model for values of $\lambda_{S}$ and $\lambda_{HS}$ which are not too large in order to avoid problems with perturbativity. As it turns out, adding VLF with Yukawa couplings big enough to have meaningful phenomenological effects further strengthens the constraints from the conditions (10)–(12), as VLF have both negative impact on the vacuum stability (see Section IV.1) and lead to stronger problems with perturbativity. The reason of the latter effect is obvious from Eq. (35): additional fermions give positive contribution to $\beta_{\lambda_{HS}}$ proportional to the sum of squares of their Yukawa coupling.

Impact of VLF on the RGE running of $\lambda$ is illustrated with examples in Fig. 5. The Higgs self-coupling $\lambda$ in the models with $S$ and VLF (green curves) is lower at low scales but bigger at high scales as compared to model without VLF (red curves). Both effects result in stronger constraints on singlet scalar couplings $\lambda_{S}$ and especially $\lambda_{HS}$. Stability condition (10) leads to lower bound on $\lambda_{HS}$ as function of $\lambda_{S}$, while perturbativity conditions (11) and (12) lead to corresponding upper bound. The impact of VLF on the allowed region for $\lambda_{HS}$ and $\lambda_{S}$ couplings is shown, depending on value of their Yukawa couplings, in Fig. 6 (for comparison, see corresponding plot without VLF shown in Fig. 4) and in Table 3

The allowed region in the $\lambda_{S}$-$\lambda_{HS}$ parameter space decreases with increasing VLF Yukawa couplings, eventually shrinking to the point defining its maximal allowed value $y_{MAX}$. Adding singlet $S$ to the theory may allow for maximal Yukawa couplings bigger by at most up to about 50% comparing to pure VLF models, as can be seen by comparing Tables 2 and 3.³³3We accept points in the parameter space which don’t satisfy condition (12) for $\beta_{\lambda}$ with $\delta=1$, when it’s clear that it’s the consequence of $\beta_{\lambda}^{(1)}(\mu)\approx 0$ and all other conditions are satisfied.

So far we assumed that there are no tree-level interactions between VLF and the singlet scalar $S$ or the SM fermions. As we discuss below, adding such couplings in general results in further shrinking of the allowed parameter space.

1. VLF and real scalar couplings.

The simplest way of including the VLF-$S$ couplings can be realised by adding the following terms to the Lagrangian:

The presence of non-vanishing VLF-$S$ Yukawa couplings leads to even stronger bounds from perturbativity condition (11). The reason is as follows: $\beta_{\lambda_{HS}}$ has a positive contributions proportional to squares of Yukawa couplings present in Eq. (LABEL:eq:LSMSVLF:Coupled:1), increasing $\lambda_{HS}$ during evolution. This indirectly, via term proportional to $\lambda_{HS}^{2}$ in $\beta_{\lambda}$ (see eq. (34)), leads to a faster growth of $\lambda$. Such effect is illustrated in the left plot of Fig. 7 where for simplicity, we assumed all new Yukawa couplings to be equal - $y_{XS}\equiv y_{S}$, $X=Q,U,D,L,N,E$. Fig. 7 shows how non-vanishing $y_{S}$ pushes $\lambda$ to larger values, implying more severe problems with perturbativity than in the case of $y_{S}=0$ and in consequence leading to smaller allowed maximal values of VLF-Higgs Yukawa coupling.

Similarly, we can write down Lagrangian terms corresponding to the interactions between VLF and SM fermions. Assuming for simplicity that the only non-vanishing couplings are between the third family of SM quarks ($q_{L}$, $t_{R}$ and $b_{R}$) and VLQ in the $n_{Q}=n_{U}=n_{D}=1$ scenario we get the following extra terms in the Lagrangian:

Non-vanishing VLQ-SM quarks Yukawa couplings have similar impact of the allowed parameter space as ordinary VLF Yukawa couplings discussed earlier in this Section. They give additional negative contributions to $\beta_{\lambda}$ resulting in problems with vacuum stability more severe than in the SM or when VLF and SM fermions are not interacting. Example of effects of non-vanishing VLQ-SM quarks Yukawa couplings for a simple case where $y_{Qt}=y_{Uq}=y_{Qb}=y_{Dq}=y_{Qq}$ is presented in the right panel of Fig. 7. One can see that even relatively small values of these additional Yukawa couplings may destabilise the EW vacuum.

Fulfilling the vacuum stability and perturbativity conditions (10)–(12) is obligatory for the theoretical consistency of any BSM model. However, one can also impose other, more optional, requirements on the considered model field content and coupling structure, in order to achieve additional advantages in its predictions. One such possible extra condition may concern gauge coupling unification.

Let us define scales $\mu_{ij}$ at which pairs of gauge couplings, $g_{i}$ and $g_{j}$, have equal values, i.e. when $\Delta g_{ij}\equiv g_{i}(\mu_{ij})-g_{j}(\mu_{ij})=0$. In the SM these three scales are quite different with $\mu_{23}$ bigger than $\mu_{12}$ by more than 3 orders of magnitude, thus there is no common point of unification of gauge couplings.

The idea of grand unification favours models for which all three scales $\mu_{ij}$ are close to each other. In addition, higher unification scales, at least of order $10^{15}$ GeV, are preferred in order to avoid too fast proton decay. It is interesting to check whether unification of the gauge couplings may be achieved in models with extra VLF multiplets modifying their RGE evolution (we do not consider models with scalar singlet here).

In Fig. 8 we show how the energy scales at which the gauge couplings may unify depend on the number of VLF multiplet families for $M_{F}=3$ TeV (we point out that VLF Yukawa $y$ has no impact on the running of gauge couplings at 1-loop). Closer look at the plots on those figures reveals the following (compare Section II for model classification):

• •

  SM + VLL models (crosses in Fig. 8):

  • –

    scenario I (top panel) and scenario III (bottom panel): gauge couplings convergence points become closer to each other with increasing $n$ but the intersection points move to lower energy scales;

  • –

    scenario II (middle panel): gauge couplings converge closer to each other for $n=1,2$ but already for $n=3$ they start diverging;

• •

  SM + VLQ models (circles in Fig. 8):

  • –

    scenario I (top panel): the $\Delta g_{23}=0$ point is above the Planck scale already for $n=1$;

  • –

    scenario II (middle panel): gauge couplings converge closer to each other for $n=1,2$ but already for $n=3$ they start diverging, with all intersection points moving towards higher energy scales;

  • –

    scenario III (bottom panel): for $n=1$ couplings converge closer to each other than in the SM, but already for $n=2$ the $\Delta g_{12}=0$ point is above the Planck scale.

As we checked, varying VLF mass $M_{F}$ between 1 TeV and 10 TeV does not make significant difference, general picture and conclusions remain the same.

Based on the discussion above, one can see that VLQ scenario II with $n=2$ (to a lesser extent also $n=1,3$) and scenario III with $n=1$ are the most promising ones for the idea of grand unification of gauge couplings in VLF models. All other cases are disfavoured due to different reasons. All VLL scenarios I, II and III with not too large $n$ have a positive impact on unification of couplings when compared to the SM. However, the corresponding unification scales are too low. The VLQ scenario I has a negative impact on unification and already for $n=2$, $g_{2}$ and $g_{3}$ couplings meet above the Planck scale. The situation is similar for the VLQ scenario III with $n\geq 2$.

Models with large cut-off scales and big number of VLF multiplets are disfavoured not only from the point of view of possible unification of gauge couplings but also by perturbativity conditions (11) and (12). New fermions give positive contributions to the gauge coupling $\beta$-functions. Too many fermions lead to too fast grow of (at least some of) gauge couplings. Demanding the perturbativity up to the Planck scale gives upper bounds on the number of VLF multiplets presented in Table 4.

An additional remark concerning models with VLL is in order here. As we discussed in Sec. IV.1, EW vacuum in such models becomes unstable at scales above about $10^{9\div 10}$ GeV. Therefore, it is reasonable to consider much bigger cut-off scales only if the singlet scalar $S$ with suitable couplings is also added to the model, improving its stability. The singlet does not change the RGEs for the gauge couplings so the results presented in this Section are valid for SM+S+VLL models.

The enhancement of the double Higgs production cross section in models with new VLF multiplets and heavy scalar $S$ can occur by modification of the triple Higgs coupling or VLQ contributions to the box and triangle gluon fusion diagrams of Fig. 9⁴⁴4Results in this Section were obtained using FeynRules [95] and MadGraph5_aMC@NLO [96] packages.. We denote differences between cross sections for the gluon fusion processes of single and double Higgs boson production in the SM and its extensions as:

where $X$ is either $S$ or VLF, and “final” state can be $H$ for single and $HH$ for double Higgs production. From the experimental point of view we are interested in scenarios which predict significant enhancement of $K^{X}_{HH}$ while keeping $K^{X}_{H}$ within present experimental bounds.

Interactions of the Higgs scalar with new particles can modify the triple Higgs coupling $\lambda_{3}$ which, if sizeable, may directly impact double Higgs boson production. This coupling is defined as:

with tree level potential given by Eq. (7) and its Coleman-Weinberg part by Eq. (45). The current allowed $95\%$ confidence level intervals for the trilinear Higgs coupling modifier $\kappa_{\lambda_{3}}=\lambda_{3}/\lambda_{3}^{SM}$ still leave a lot of room for the BSM physics and read:

Corrections to the leading tree level contribution $\lambda_{3,\text{tree}}^{SM}=3M_{H}^{2}/v$ predicted by the SM can be classified as follows.

• •

  Leading SM one-loop top-quark induced contribution:

  $$\Delta\lambda_{3}^{SM}\approx-\frac{3M_{t}^{4}}{v^{3}\pi^{2}}.$$

  (18)

• •

  Leading one-loop real scalar singlet induced contribution:

  $$\Delta\lambda_{3}^{S}\approx\frac{\lambda_{HS}^{3}v^{3}}{32\pi^{2}M_{S}^{2}}.$$

  (19)

• •

  Loop contributions from heavy VLF multiplets:

  $\displaystyle\Delta\lambda_{3}^{VLF}\approx n_{F_{1}}\frac{N_{c}^{\prime}v^{3}y_{F_{1}}^{6}}{8\pi^{2}M_{F_{1}}^{2}}+n_{F_{2}}\frac{N_{c}^{\prime}v^{3}y_{F_{2}}^{6}}{8\pi^{2}M_{F_{2}}^{2}},$

  (20)

  where $N_{c}^{\prime}$ is the number of colours of VLF multiplets and one should sum over $n_{F_{1}}\in\{n_{U},n_{N}\}$, $n_{F_{2}}\in\{n_{D},n_{E}\}$ and $y_{F_{1}}\in\{y_{U},y_{N}\}$, $y_{F_{2}}\in\{y_{D},y_{E}\}$.

Effects of the scalar singlet are presented in Fig. 10 which shows the corresponding relative enhancement of the triple Higgs coupling, $\kappa_{\lambda_{3}}^{SM+S}=\lambda_{3}^{SM+S}/\lambda_{3}^{SM}$. The contours of maximal allowed values of $\lambda_{HS}$ for chosen cut-off scales $\Lambda$ are obtained using perturbativity conditions (10)– (12).

For the cut-off scale of 100 TeV or lower, the enhancement of $\lambda_{3}$ from scalar singlet contributions is at most $\sim 10\%$ for $M_{S}\geq 200$ GeV, leading to a moderate, at most $K^{S}_{HH}\approx-8\%$, decrease of the double Higgs production cross section (which is due to destructive interference between the box and triangle diagrams of Fig. 9).

Addition of VLF also has some impact on $\lambda_{3}$, however, this effect is very small. We can estimate its maximal possible value by using numbers from Table 3 and plugging them into the Eq. (20). As a result, we get negligible contributions with $|\Delta\lambda_{3}^{VLF}/\lambda_{3}^{SM}|\ll 1\%$ for all considered cases, leaving no chance for significant triple Higgs coupling modification by VLF fields only.

Apart from modifying the triple Higgs boson coupling $\lambda_{3}$, heavy VLQ may impact double Higgs boson production through the box and triangle loop gluon fusion diagrams presented in Fig. 9. Considering such corrections, one should take into account that loop contributions to the double and single Higgs boson production are closely related, the Feynman diagrams for the latter can be obtained from Fig. 9 by replacing one of the external Higgs fields by the VEV insertion. As a consequence, all relevant amplitudes are proportional to the same combination of VLQ parameters:

Therefore, apart from theoretical constraints considered already in this work, another important source of limits on VLQ model parameters comes from the current limits on the single Higgs production. Available experimental data [89] puts strong constraint, at around $10\%$ at $68\%$ CL and around $18\%$ at $95\%$ CL, on the deviations of the single Higgs boson production rate via gluon fusion from its SM value. This condition may (depending on the choice of $\Delta$ in (12)) lead to stronger constraints on the maximal values of VLF Yukawa couplings then listed in Table 3, as illustrated in Table 5.

Table 5 shows that the enhancement of the double Higgs production rate due to the VLQ loop contributions, $K_{HH}^{VLQ}$, within the parameter space obtained in this work, is always small (at most $\sim 15\%$). Even such minor enhancement should be considered as optimistic scenario. As discussed before, in models with scalar singlet this positive contribution may be in addition compensated by scalar contribution to the triple Higgs coupling. Increasing $M_{F}$ can also only further suppress the VLF effects. Finally, effects of increasing number of VLF multiplets can be deduced from the form of the amplitudes in (21) - for example, assuming identical masses for all multiplets, case with $n=1$ and some value of $Y$ is equivalent to $n=N$ with $Y\to Y/\sqrt{N}$, so that the corresponding bounds on $K_{HH}^{VLQ}$ are not much affected by varying $n$ (as can be seen in the 4th column of Table 5).

Comparing numbers in Table 5 with still rather lax current experimental constraints on the Higgs boson pair production at $95\%$ confidence level:

one can conclude that adding vector-like multiplets to the SM likely cannot affect the double Higgs production in a way which could be experimentally observed at present or in the foreseeable future.

Additional constraints on the parameter space of models considered in this work may arise from the electroweak precision observables. In our analysis, we focus on $\mathbb{S}$, $\mathbb{T}$ and $\mathbb{U}$ oblique parameters, which parametrize one-loop contributions to the electroweak gauge bosons self-energies [97, 98]. In the context of VLF (the singlet scalar does not contribute to $\mathbb{S}$, $\mathbb{T}$ and $\mathbb{U}$ parameters due to $\mathbb{Z}_{2}$ symmetry of its potential), there have been a number of extensive studies addressing this issue, see e.g. [55, 56, 22, 57]. Current experimental constraints, obtained with $\mathbb{U}$ fixed to zero (motivated by the fact that $\mathbb{U}$ is suppressed by additional factor of $M_{i}^{2}/M_{Z}^{2}$ comparing with $\mathbb{S}$, and $\mathbb{T}$, see Ref. [99]), read [92]:

Formulae for the corrections to $\mathbb{T}$ and $\mathbb{S}$ parameters in the presence of VLF have been worked out in [55] and are collected in Appendix B. Those parameters by definition come only from the BSM sector, so that in the absence of mixing between VLF and the SM fermions one has:

Comparison of constraints on VLF models based on theoretical considerations of model consistency and/or of single Higgs production presented in the Section V.1 (see Table 5), with those based on the precision $\mathbb{T}$ and $\mathbb{S}$ electroweak observables is illustrated in Table 6 and Figure 11.

Several comments are in place. In the model scenarios considered in this work, we assume the same number of left- and right-handed VLF multiplets. This significantly simplifies formulae from [55]. Moreover, we assume identical VLF masses $M_{VLF}$ and Yukawa couplings $y_{VLF}$, which limits the effect of VLF on $\mathbb{T}$ and $\mathbb{S}$ parameters. More specifically, for “Scenario I” in which we have the same number of VLF doublets and “up-type” and “down-type” singlets one has $\mathbb{T}=0$ as such setup preserves custodial symmetry. For other scenarios, contribution to $\mathbb{T}$ is significantly larger than to $\mathbb{S}$, but still way below the experimental constraints. Finally, VLL scenarios are far less constrained than VLQ due to the lack of enhancement by the $N_{c}$ colour factor (see formulae in Appendix B). Since increasing $M_{VLF}$ only leads to further decrease in the values of $\mathbb{T}$ and $\mathbb{S}$ parameters (even after accounting for the increase of allowed maximal value of VLF Yukawa couplings), we can conclude that for the model scenarios studied in this paper, the theoretical stability constraints are much stronger than those arising from the bounds on the oblique parameters $\mathbb{T}$ and $\mathbb{S}$.

The one-loop effective potential in the SM with additional singlet scalar field (for review see e.g. [100]) reads:

with $V_{0}$ a tree-level potential (7), $V_{CW}$ the one-loop correction known also as the Coleman-Weinberg potential [101] and $V_{T}$ the finite temperature contribution (for further details see Appendix C).

Models are considered to be interesting in the context of baryogenesis if they allow for a strong enough first order phase transition. In general, a 1st order electroweak phase transition may occur in the presence of a barrier in the effective Higgs potential separating two degenerate minima at $h=0$ and $h=v_{C}$ (here, $v_{C}$ is the VEV of the Higgs field at the critical temperature, $T=T_{C}$). The phase transition is considered strong if the following condition is satisfied:

In the real scalar model with the tree-level potential given by Eq. (7), depending on the values of the scalar sector couplings and masses, we can get either 1- or 2-step EWPT. The former takes place when the system moves from the high-temperature minimum with $\left<H\right>=\left<S\right>=0$ to the standard EW vacuum with $\left<H\right>\neq 0$, $\left<S\right>=0$. The latter is realised if for some range of temperatures the minimum of the effective potential is located at $\left<H\right>=0$, $\left<S\right>\neq 0$.

• •

  1-step EWPT is realised if $\mu_{S}^{2}>0$. This region is not affected by the value of $\lambda_{S}$, therefore we choose $\lambda_{S}=0$ in order to maximise the allowed value of $\lambda_{HS}$ coupling. The allowed parameter space for the 1-step EWPT is presented on the left panel in Fig. 12. The blue line indicates the boundary between 1- and 2-step regions. The red and black lines indicate maximal allowed values of $\lambda_{HS}(M_{S})$ obtained with constraints (11) and (12) assuming $\Lambda=100$ TeV for two values of the VLF Yukawa couplings. The orange dashed lines indicate values of $\lambda_{HS}(M_{S})$ corresponding to the strength of the EWPT, $\xi$ equal 0.6 and 1.0. As one can note, the region for which $\xi\geq 0.6$ in the singlet scalar extension is basically excluded in comparison e.g. with corresponding area presented in [72]. This follows from our more careful treatment of perturbativity conditions.

• •

  2-step EWPT may be realised if $\mu_{S}^{2}<0$ with further condition guaranteeing that at zero temperature the EW vacuum corresponds to the global minimum of the potential:

  $$V_{0}(v,0)<V_{0}(0,w)\,,$$

  (27)

  which can be translated to the bounds on $\lambda_{HS}$. A 2-step EWPT is possible if:

  $$2\frac{M_{S}^{2}}{v^{2}}\leq\lambda_{HS}\leq 2\frac{M_{S}^{2}}{v^{2}}+\sqrt{2}\frac{M_{H}\sqrt{\lambda_{S}}}{v}\,.$$

  (28)

  The allowed 2-step region is indicated by the red colour on the right plot in Fig. 12. The lower and upper limits from the Eq. (28) for $\lambda_{S}=0.5$ are indicated by the blue lines. The choice of this particular value of $\lambda_{S}$ allows for reasonably large $\lambda_{HS}$ while still leaving space for the 2-step EWPT to occur. The red line, as in the case of 1-step transition, show the maximal possible values of $\lambda_{HS}$.

In principle VLF can also affect the mechanism of the EWPT through thermal and loop corrections to the effective potential (see Appendix C). Fermions do not contribute much to inducing potential barrier for the Higgs field, contrary to the case of bosons which in general lead to cubic terms in the potential in the high temperature expansion of thermal effects. However, in the case in which the critical temperature is much lower than fermion masses, the small temperature expansion is more suitable. In the latter case, which fits very well the scope of this work in which VLF are relatively heavy, the leading contribution to the thermal potential is similar for fermions and bosons. This could potentially lead to non-trivial impact of VLF on EWPT and thermal history of the Universe (see e.g. [27, 27]).

However, because of strong constraints on the values of the VLF Yukawa couplings, discussed in Section IV, the actual impact of these particles on the EWPT is very limited. For the maximal allowed values of $y_{F}$ in SM+VLF models, presented in Table 3, we obtained $\xi_{SM+VLF}$ quite close to $\xi_{SM}$. Addition of the real singlet scalar allows for slightly larger maximal values of $y_{F}$, see Table 3, but these larger values are not sufficient to have any significant impact on the EWPT, so $\xi_{SM+S+VLF}\approx\xi_{SM+S}$. We checked that changes of $\xi$ due to VLF are at most of the order $\mathcal{O}(10\%)$.