Singlet-assisted electroweak phase transition at two loops

We investigate the electroweak phase transition in the real-singlet extension of the Standard Model at two-loop level, building upon existing one-loop studies. We calculate the effective potential in the high-temperature approximation and detail the required resummations at two-loop order. In typical strong-transition scenarios, we find deviations of order $20\% - 50\%$ from one-loop results in transition strength and critical temperature for both one- and two-step phase transitions. For extremely strong transitions, the discrepancy with one-loop predictions is even larger, presumably due to sizable scalar couplings in the tree-level potential. Along the way, we obtain a dimensionally-reduced effective theory applicable for non-perturbative lattice studies of the model.


I. INTRODUCTION
Phase transitions occurring in various extensions of the Standard Model (SM) have attracted a lot of attention in recent years.Many candidate theories for physics beyond the SM (BSM) involve scalar fields in addition to the SM Higgs, and it is interesting to ask if there can be phase transitions associated with the scalar potential.A cosmological first-order phase transition around temperatures comparable to the electroweak (EW) scale could produce gravitational waves with observational prospects at future detectors [1,2], and may provide the necessary conditions for EW baryogenesis [3,4] (see [5] for a review).
In the SM, the phase structure of the EW sector is well established: there is no phase transition associated with the Higgs mechanism [6,7].Instead, the early universe smoothly interpolates between the hightemperature "symmetric" and "broken" Higgs regimes, and a first-order discontinuity between the two would require the Higgs particle to be unphysically light (mass < ∼ 70 GeV).For any BSM scenario to change this conclusion, drastic modifications to the Higgs potential are necessary.First-order electroweak phase transition (EWPT) can be achieved by coupling the Higgs to new scalars at sufficient strength, or the phase structure itself can be more complicated with the theory undergoing multiple phase transitions at temperatures comparable to the EW scale.In both cases the new particles cannot be much heavier than the Higgs, making the EWPT a potential target at present and planned collider experiments [8].For cosmology, one is mainly interested in strong transitions with large latent heat available for the production of gravitational waves, or with significant discontinuity in the Higgs VEV (v/T > ∼ 1) as required for efficient baryogenesis.
Unfortunately, studying the EWPT with perturbation theory is complicated because (i) first-order transitions in simple scalar extensions typically require O(1) couplings among the scalars, and (ii) finite-temperature perturbation theory is sensitive to infrared (IR) physics and has a different structure than at zero temperature, often leading to slow convergence.Because of the second point, perturbation theory breaks down near second-order transitions where the Higgs correlation length diverges, making it impossible to determine the order of the EWPT with perturbative methods alone.Hence, when applying perturbation theory to first-order transitions one simply assumes that the transition remains first order also at the non-perturbative level, and that predictions are not significantly altered by higher order effects.Because of points (i) and (ii), it is far from obvious whether such assumptions are justified in typical BSM settings.
The situation can be contrasted to that of the SU(2)-Higgs theory relevant for the SM, where the dominant contributions to the free energy come from gauge loops.These are under perturbative control for v/T > ∼ 1, although there are still quantitatively important corrections at two-loop order [9].In many BSM scenarios, corrections from the scalar sector dominate over the gauge loops, and predictions may well be sensitive to higherorder scalar diagrams even when the transition is strong.
The standard tool for studying the EWPT in BSM settings is the one-loop effective potential with resummation of thermal "daisy" diagrams.The purpose of this paper is to extend this calculation to two-loop level, and scrutinize the effects of the two-loop corrections on basic characteristics of the EWPT (critical temperature, latent heat, v/T ).We work in the real-singlet extension of the Standard Model which we shall refer to as xSM, following [10].Several one-loop analyses of the EWPT exist in this minimalistic yet realistic BSM model, see [10][11][12][13][14][15][16][17][18][19][20][21][22][23] and references therein, but we are not aware of existing two-loop calculations in the finite-T context.Two-loop studies in other scalar extensions and in the MSSM have reported on considerable deviations from the one-loop behavior [24][25][26][27][28][29], motivating a detailed two-loop investigation of the xSM as well.
Our calculation proceeds in two stages.First, we apply effective field theory (EFT) methodology to integrate out hard thermal loops (loops involving thermal excitations at the scale πT ) in the static limit, leaving a three-dimensional (3d) EFT describing the long-distance physics of the phase transition [30,31].This dimensional reduction automatically incorporates the required resummations while also summing a subset of higher-order corrections using the renormalization group (RG).We then calculate the effective potential in the resulting EFT.The 3d EFT constructed in this paper can also serve as a starting point for future non-perturbative lattice studies [32] of the xSM phase structure.
Applying the two-loop effective potential, we compute the strength and critical temperature of the EWPT in different scenarios, including two-step and purely radiatively generated transitions.The analysis is concentrated on a handful of benchmark points to better illustrate the impact of two-loop corrections on physical predictions.The main takeaway is that the corrections are significant even if the transition occurs through a tree-level barrier, and should be taken into account if reliable predictions of e.g. the associated gravitational-wave background are sought.
The outline of this paper is as follows.We begin by introducing the model in Section II.Section III reviews some formal aspects of thermal resummation and high-T dimensional reduction that form the backbone for understanding the remainder of the paper.In Section IV we introduce a formal power counting scheme for various parameters in the theory and discuss subtleties specific to the high-T behavior of the xSM.Sections V and VI contain details of the EFT construction and collect the relevant 4d → 3d matching relations.In Section VII we calculate the resummed effective potential using the EFT approach.Finally, our results concerning the importance of two-loop corrections are presented in section VIII.Section IX concludes.There are two appendices: Appendix A details our renormalization prescription and loop-corrected relations to physical observables.Appendix B calculates the two-loop correction to the effective potential.

II. MODEL
We work on a model consisting of the usual SM field content and an additional real scalar S, that is a singlet under the SM gauge group.The Lagrangian in Euclidean signature is where L YM is the Yang-Mills part for the SM gauge fields (including hypercharge), L fer contains the fermion kinetic terms and y t is the top Yukawa coupling.Yukawa couplings of the light fermions are orders of magnitude smaller than other parameters in the theory and are neglected.Covariant derivative for the Higgs doublet φ is , where g and g ′ are the SU(2) L and U (1) Y gauge couplings respectively.The scalar potential is As usual, the zero-temperature theory can be studied by assuming m 2 φ < 0 and specifying components of the doublet as where v 2 = −m 2 φ /λ is the gauge-fixed Higgs vacuum expectation value (VEV) at tree level, and G ± , G are the would-be Goldstone modes.The singlet VEV and hence the linear coupling b 1 can be absorbed into a redefinition of the tree-level parameters, i.e. the theory contains five free parameters instead of six.At tree level, the condition S = 0 fixes b 1 = −a 1 v 2 /4.In the limit b 1 , a 1 , b 3 → 0 the theory becomes symmetric under S → −S; we shall call this the Z 2 symmetric case.Outside the Z 2 symmetric limit the physical scalar excitations h 1 , h 2 are mixtures of h and S: The angle θ is chosen to diagonalize the mass matrix at tree level and we will identify the lighter h 1 as the SM-Higgs-like excitation with mass 125 GeV.
The input parameters we use for our analysis are b 3 , b 4 , a 2 , sin θ and the masses of h 1 , h 2 (the VEV v can be fixed using the Fermi constant).Inverting the defining tree-level equations for these leads to the following expressions for the parameters in (2): Without loss of generality, sin θ can be restricted to values | sin θ| < 1/ √ 2 [11].
These relations obtain corrections at higher orders of perturbation theory.At one-loop, the corrections are parametrically of the same order as two-loop corrections to the thermal effective potential.Hence, a two-loop investigation of the finite-temperature behavior should be accompanied by a one-loop renormalization of the zerotemperature theory.This is a separate calculation from the main focus of the paper, and we perform it in Appendix A using the MS scheme.This fixes our renormalized parameters at the Z boson pole.The main effect of the one-loop improvement is to shift the singlet mass parameter m 2 S , by as much as 50% relative to the treelevel result in some of our benchmark points below.The parameters b 1 and λ also get considerable corrections.
The input parameters are constrained experimentally by electroweak precision measurements and direct searches for new physics [33,34].We will study the theory near its Z 2 symmetric limit, see the discussion in Section IV below, taking | sin θ| < 0.1 and assuming the new excitation to be heavier than the SM Higgs.The remaining parameter space is then relatively unconstrained by current experiments [14,20,23,35].We will not consider experimental constraints further in this paper as our focus is on the finite-T behavior of the theory.

III. REVIEW OF THERMAL RESUMMATIONS AND HIGH-T EFFECTIVE THEORIES
Moving on to the finite-temperature theory, let us briefly review the premise of thermal corrections that are important for phase transitions.Much of the discussion here is collected from refs.[9,25,36,37].Readers familiar with thermal resummations and dimensionally-reduced theories may want to skip this section.
Equilibrium quantum field theory is equivalent to a Euclidean field theory with a compact time direction, where the periodicity of the imaginary time coordinate τ is identified with 1/T .The fields have momentum-space representations of the form where the Matsubara frequency is ω n = 2πnT for bosons and ω n = (2n + 1)πT for fermions.The bosonic zero modes with ω n = 0 carry no momentum in the time direction, so their dynamics is effectively three-dimensional.At length scales ≫ 1/T the zero-mode contributions dominate, leading to the well-known IR problems in perturbation theory because the leading interaction terms in 3d have dimensionful couplings [38].In contrast, modes with ω n = 0 (non-zero Matsubara modes) are safe in the IR due to additional factors of (πT ) 2 in the propagators.Qualitatively, the EWPT can occur because of a Tdependent correction to the Higgs mass parameter from thermal loops, schematically m 2 φ → m 2 φ (T ) = m 2 φ + g 2 T 2 .At large enough T the VEV is relaxed to zero, and at the critical temperature the "symmetric" and "broken" phases have equal free energies.In the absence of a singlet VEV, this happens, at the mean field level, when m 2 φ (T ) = 0 and requires a cancellation between the vacuum and thermal masses. 1 Such large corrections to the tree-level behavior call for resummation of specific thermal corrections.For the one-loop effective potential this is achieved by the ring ("daisy") resummation where the propagators are replaced with thermally-corrected ones [39].
At two-loop level, a consistent prescription has been given in ref. [9]: thermal masses for the zero modes are generated by integrating out the heavier modes with ω 2 n ≥ (πT ) 2 , and resummation of the non-zero modes is not needed at all.Such prescription is possible in the formal high-T limit where the IR degrees of freedom are lighter than πT .Outside this limit, there is no simple way of choosing thermal masses for resummation and the procedure becomes ambiguous beyond one-loop level [27].In many cases the high-T assumption is justified because the (thermally corrected) Higgs mass is parametrically small compared to T near the critical temperature, and possible BSM excitations cannot be arbitrarily heavy if they are to produce a first-order EWPT [8].
To understand the necessary resummations beyond one-loop order, it is convenient to reformulate the problem in the language of effective field theory.The finite-T theory contains several mass scales in addition to those present at zero temperature.The most obvious one is the scale πT , characterizing the non-zero Matsubara modes.There are also the scales associated with Debye screening of electric and magnetic fields, parametrically gT and g 2 T , respectively.The latter appears in non-Abelian gauge theory and is inherently non-perturbative [40].Conventional terminology for the scales πT , gT , and g 2 T is "hard", "soft" and "ultrasoft" (or sometimes "superheavy", "heavy" and "light"; we will stick to the former).In a weakly coupled theory, these are all distinct scales.For scalars, a similar hierarchy of scales between the zeroand non-zero Matsubara modes is present in a window around T c where the thermal mass is comparable to the vacuum one.
The strategy is now as follows.Making use of the thermal scale hierarchy, we construct an effective theory for the Matsubara zero-modes by integrating out all field modes with ω n = 0, which includes all fermionic fields.This leaves a 3d theory for the zero modes, describing the static properties of the finite-T theory at length scales > ∼ (πT ) −1 .The effects of the hard scale are incorporated in tree-level parameters of the EFT, and the process of thermal resummation at any order of perturbation theory is then equivalent to calculating the EFT parameters to corresponding accuracy.Consequently, the effective potential calculated within the dimensionally-reduced EFT reproduces the resummed high-T effective potential [36].
The form of the EFT is governed by 3d gauge invariance.In particular, under static gauge transformations the gauge-field time components transform like adjoint scalar fields.They generate a thermal mass due to Debye screening of electric fields in the thermal medium.Therefore, the dimensionally-reduced theory corresponding to (1) has an action of the form2 + (operators of dimension 5 and higher) (11) and is valid in the high-T limit.Here L YM is the Yang-Mills part in three spatial dimensions and A 0 is the adjoint scalar corresponding to the time component of the SU(2) L gauge field; the mass m D is the Debye screening mass.The corresponding scalars for the hypercharge field and QCD gluons are not shown for simplicity.The potential V (φ, S) is of the same form as in (2) but with modified parameters.The overall factor 1/T comes from a trivial integration over the imaginary time and can be absorbed in a rescaling of the fields and couplings.At leading order one matches masses at one loop and couplings at tree level, which is precisely the standard procedure for daisy resummation.
The higher-dimensional operators in (11) are suppressed by powers of the ultraviolet (UV) scale πT , but Ref. [42] pointed out that their effects do not fully decouple even in the T → ∞ limit.In many cases it is still a good approximation to neglect them, provided that all mass scales appearing in the EFT (including cubic couplings) are small compared to πT [32].It should be emphasized that the 3d approach is not suited for describing heavy BSM fields if their zero-mode masses are comparable to πT .However, thermal effects from such fields are Boltzmann suppressed and can be incorporated by integrating out the field together with the hard Matsubara modes.It is also difficult to obtain a first-order EWPT in such a scenario, unless there are non-negligible effects from the higher-dimensional operators arising from integration over the heavy field [35,43].
If the higher-dimensional operators are dropped by truncation, the remaining EFT is super-renormalizable.In this case a non-perturbative approach to the phase transition is possible using relatively simple lattice simulations [32,44].Hence, our calculations in the following sections are a pre-requisite for future non-perturbative studies of the singlet-extended SM, but in the present paper we restrict ourselves to perturbation theory.
A further simplification can be achieved by explicitly integrating out the adjoint scalar A 0 , whose mass is of the order gT .The resulting theory is formally valid for momenta at the ultrasoft scale, k < ∼ g 2 T , and has the same form as (11) but without the A 0 field.This is often a good approximation even when parametric mass hierarchies are not strictly satisfied, because of small numerical factors involved in the integration over A 0 [37].
The reader may wonder if purely perturbative studies of the effective potential actually benefit from the EFT approach, which does seem like an extra step compared to the more direct calculation described in [9].The two approaches are, of course, equivalent, since the resummation of ref. [9] also utilizes approximate decoupling of hard Matsubara modes.In the perturbative context, the main advantage of the EFT lies in factorization: physics at the hard thermal scale is accounted for in the EFT matching and can, in principle, be calculated to any order in perturbation theory, while the IR sensitive contributions are obtained by separate calculations within the EFT.There are also additional resummations of coupling constants and those related to the A 0 field that are typically not included in the 4d calculations but are straightforward to implement using EFT methods.In the context of hot QCD, where dimensionally-reduced theories are widely used, including these higher-order effects is known to improve the convergence of perturbation theory [45,46].

IV. POWER COUNTING AND PARAMETRIC ESTIMATES
Let us now discuss the accuracy goal of our calculation.Near a phase transition, there is competition between the vacuum and thermal masses, so we assume that our tree-level masses are not considerably heavier than T .We express this parametrically as m φ , m S ∼ gT , where g is a formal expansion parameter that will be identified as the SU(2) L gauge coupling for concreteness.For the dimensionful cubic couplings, we consequently require a 1 , b 3 ∼ g 2 T ; otherwise ratios such as a 1 /m φ are unsuppressed and perturbativity is questionable.For other couplings, the renormalization structure of the theory suggests the parametric countings λ, a 2 , b 4 ∼ g 2 and y t , g ′ , g s ∼ g (here g s is the QCD coupling).Adopting this schematic counting facilitates the discussion of loop corrections where combinations of the various couplings appear.
The structure of perturbative expansion at finite T is modified by resummations and differs from that of its zero-temperature counterpart.Loops involving hard Matsubara modes contribute to the effective potential at orders g 2 , g 4 , . . .while the resummed zero modes provide corrections also at g 3 , g 5 , . . . .A one-loop calculation with daisy resummation gives the field-dependent terms in the effective potential at O(g 3 ) with partial O(g 4 ) contributions, but reaching g 4 accuracy in the quadratic terms requires a two-loop calculation.In particular the relative difference between a one-loop daisy resummed potential and a two-loop calculation is O(g), while at zero temperature the difference between one and two loops is O(g 2 ).Thus the two-loop improvement is relatively more important in the finite-T theory.We will perform a full O(g 4 ) computation applying the 3d EFT approach.
The cubic couplings appearing in the xSM have implications for the parametric accuracy of dimensional reduction.In the absence of chemical potentials, the leading higher-dimensional operators appearing in the dimensionally-reduced EFT of the SM are of dimension six (in 4d units) [47].This is also true for the Z 2symmetric case of xSM.However, in the more general case there can be Z 2 -breaking operators like c 4,1 (φ † φ) 2 S at dimension five.The coefficient is parametrically c 4,1 ∼ g 4 a 1 /T 2 .With this operator present in the 3d EFT, it contributes, among others, to the φ † φS 2 correlator via the diagram Since the typical mass scale of the 3d EFT is m ∼ gT , this diagram contributes at order a 4 1 /T 4 .But there is a one-loop contribution of the same parametric order from integrating out the hard thermal modes.We conclude that that as far as arbitrary Green's functions are concerned, there is no reason to keep terms proportional to a 4  1 in the EFT matching relations unless dimension five operators are included as well. 3 To keep the discussion simple, we will only include operators up to dimension four in the EFT.In terms of our power counting, the maximal accuracy that can be reached for field-dependent terms in the effective potential is then O(g 4 ), sufficient for our purposes.The EFT reproduces hard-mode contributions at this accuracy provided that we match 1-and 2-point Green's functions at two-loop order, and 3-and 4-point functions at one loop.Terms such as a 4  1 and g 2 a 2 1 will be dropped, following the power counting described above.An additional advantage of dropping the higher-dimensional operators is that the resulting EFT admits a simple continuum limit, facilitating non-perturbative lattice studies in the future.

V. MATCHING TO THE 3D EFT
Building on the discussion of section III, the dimensionally reduced EFT for the full xSM has the Euclidean action We use "natural" 3d units where all couplings are dimensionful and fields have mass dimension 1/2.Here σ a are the Pauli matrices and F ij , B ij are 3d field strengths for the SU(2) and U(1) Y gauge fields whose couplings are denoted by g 3 and g ′ 3 .The fields A 0 , B 0 , C 0 are scalars in adjoint representations of SU(2), U(1) Y , SU(3), respectively (i.e.B 0 is a singlet). 4Their masses are parametrically gT, g ′ T, g s T due to Debye screening; we will 3 This conclusion can also be reached via the general powercounting arguments given in [47]. 43d gauge invariance also allows for operators with an odd number of the adjoint fields, such as φ † A a 0 σaφ and φ † φB 0 .However, these can only appear in the presence of non-vanishing chemical potentials [47].
shortly integrate these fields out to considerably simplify the EFT.
Couplings between the singlet and the adjoint fields are generated by φ loops in the 4d theory and the φ † φC 0 C 0 term by quark loops.In principle the SU(3) gluons are also present in the EFT, but these couple to the Higgs sector only through C 0 and any such contribution to the V eff is of higher order than our accuracy goal.Interaction vertices among the adjoint scalars are O(g 4 ) and contribute only at higher orders [49].The scalar poten- FIG. 1. Singlet contributions to the φ † φ self energy.Oriented dashed lines refer to Higgs doublets, double lines to the singlet S, and the two types of wiggly lines describe SU(2) and U(1) Y gauge propagators.Black blobs represent counterterm insertions.
In our power counting, two-loop diagrams involving cubic couplings a1 or b3 contribute at higher orders and are not shown.The diagrams were drawn with Axodraw [48].
tial for φ and S in the EFT is To avoid clutter, we will use the same notation for fields in the original theory and those in (13).
To relate the EFT parameters to those appearing in the full xSM in 4d, we identify the fields in (13) as the Matsubara zero-modes of the original theory and integrate out all loops involving a non-zero (hard) Matsubara mode.Effects of the hard modes approximately decouple in the high-T limit, apart from a T -dependent renormalization of parameters and fields in the zero-mode sector.It is convenient to make these corrections explicit in the EFT by writing, for instance m 2 φ,3 = m 2 φ + δ T m 2 φ , where δ T m 2 φ is the T -dependent correction from hard modes ("thermal counterterm").Here m 2 φ,3 is the renormalized parameter; a separate counterterm δm 2 φ,3 is introduced to absorb UV divergences in dimensional regularization.In the matching we treat both counterterms as perturbations.We then calculate 1PI Green's functions both in the EFT and in the full theory at soft external momenta, solving for the matching corrections order-by-order in perturbation theory by requiring that the Green's functions match (note that hard Matsubara modes cannot appear in reducible internal legs).
The required matching relations have been worked out previously in ref. [50], where terms such as g 2 a 2 1 and a 4 1 were included.As discussed above, one then has corrections of the same order from higher-dimensional operators, and the resulting EFT becomes complicated.Here we shall stick to the power counting of the previous section and omit such terms, whose effects are negligible unless the cubic couplings are unnaturally large compared to masses appearing in the quadratic terms.We have checked this numerically in the benchmark points discussed below.
For the matching, we can expand the fields around any constant background value that is parametrically smaller than the UV scale πT .Any such field shifts affect only the IR sector and have no effect on the matching.Therefore we choose to expand around the origin so that no background fields are explicitly present in the calculation.The resulting EFT correctly describes the phase structure of the full xSM, provided that VEVs after the EWPT are not considerably larger than πT ; otherwise there will be unsuppressed contributions from higherdimensional operators.We take the linear term b 1 S to be an additional interaction rather than part of the free particle action.This is allowed because it only creates oneparticle-reducible insertions that can be incorporated in the EFT by matching the singlet 1-point correlator at O(g 4 ).
Detailed discussions of matching calculations in the high-T context can be found for instance in [37,41,50,51].Here we cut down on details, showing explicitly only the matching of singlet contributions to the Higgs mass parameter m 2 φ,3 .The relevant diagrams are shown in figure 1.At two-loop level there are diagrams where both soft and hard momenta flow in loops, however these are absorbed in the EFT by one-loop diagrams containing a thermal counterterm insertion [37]. 5Hence, only diagrams with hard momenta in all loops are needed for the matching.The sum-integrals are calculated using high-T expansion, i.e. an expansion in m 2 /(πT ) 2 , and we assume m ∼ gT .For O(g 4 ) accuracy we need to expand the oneloop diagrams to next-to-leading order, while at two-loop level the leading order suffices.
We denote sum-integrals over non-zero modes as where P = (ω n , p), d = 3 − 2ǫ, μ is the MS scale and γ is the Euler-Mascheroni constant.At one loop, the relevant bosonic sum-integrals are [52] At vanishing external momentum, the one-loop diagrams in fig. 1 give Here the a 2 1 m 2 S term exceeds our desired accuracy.The first term is the O(g 2 ) singlet contribution to thermal Higgs mass, δ T m 2 φ = (δ T m 2 φ ) SM + T 2 24 a 2 + O(g 4 ).On the EFT side, we renormalize the field as φ → (1 + δ T Z φ ) 1/2 φ, where the counterterm arises from integration over the hard modes at soft external momenta in the full theory.Then, in the EFT Lagrangian we have , plus a mass counterterm needed to cancel 1/ǫ poles in 3d.For O(g 4 ) matching, the field renormalization counterterm needs to be known at O(g 2 ) accuracy.Momentum-dependent contributions from the singlet are ∼ k 2 a 2 1 /T 2 ∼ g 6 T 2 for soft k ∼ gT , therefore the SM result is sufficient (equation (B.37) in [53]): which holds in a general covariant gauge (ξ is the gaugefixing parameter).Here L b and L f are logarithms that arise frequently from one-loop bosonic and fermionic sum-integrals: Gauge dependence of ( 19) cancels against two-loop contributions in the final result for m 2 φ,3 .We then turn to the two-loop diagrams in fig. 1, which we evaluate at vanishing external momentum and leading order in the high-T expansion (i.e.propagators can be taken massless).Diagrams (c), (d), (e) are products of one-loop sum-integrals: The "sunset" diagram (f ) does not contribute to matching at O(g 4 ) because of the identity [54] P Q Diagrams (g 1 ), (g 2 ) contain gauge fields, but the resulting integrals are trivial because the gauge propagator D µν (P ) with loop momentum P gets contracted with P µ P ν , leaving only the longitudinal part.In a general covariant gauge, the result is Finally, there are one-loop diagrams with counterterm insertions.Our counterterms for the MS scheme are listed in Appendix A.Here we note that mass counterterms contribute finite parts at O(g 6 ) in our counting and that the singlet does not need field renormalization at one loop.Thus diagram (c 1 ) does not contribute while (c 2 ) is proportional to δa 2 .We also need an additional diagram proportional to δλ (not shown in fig.1), because the Higgs self-interaction gets renormalized by singlet loops.The required diagrams are where, for the present discussion, we only include singlet contributions to δλ.Summing all the diagrams, expanding in ǫ and discarding terms that go beyond O(g 4 ) gives Γ (S) which is the O(g 4 ) correction to the Higgs self energy due to hard singlet loops.There are also leftover 1/ǫ divergences that are canceled by contributions from the soft modes, but this is irrelevant for our present discussion.Here The self energy (26) is to be reproduced by tree-level terms in the EFT, so we match as from which the matching correction δ T m 2 φ can be solved order-by-order, i.e. the one-loop result is used in the δ T m 2 φ δ T Z φ term (the m 2 φ term actually contains no singlet contributions).Accounting for the full SM part, given in eq.(A28) of [55], our result for the effective Higgs mass parameter is where gauge dependence duly cancelled.This is a twoloop generalization of the thermal Higgs mass.
Eq. ( 29) can be improved by considering RG behavior of the 3d EFT.Applying RG evolution to the parameters in (29), one can check that apart from the ln(3T /μ)-term, the effective mass is RG invariant up to higher-order effects.This means that RG evolution down from the matching scale is governed by the ln(3T /μ)term.But we can also solve the running directly in the super-renormalizable EFT, using the renormalization scale μ3 , to obtain an exact6 RG equation at two-loop order [36,37].Replacing the ln(3T /μ)-term with the solution of this RG equation incorporates higher-order effects to the running of m 2 φ,3 , in this case we simply replace and similarly in the SM part (c.f.[55]).Written in this form, the effective mass is independent of the matching scale μ at order O(g 4 ).We apply this same improvement to parameters b 1,3 and m 2 S,3 below, while all couplings are RG invariant in 3d.This is an example of a higher-order improvement that is difficult to incorporate directly in the 4d effective potential.Numerically the effect is small, for instance in benchmark points BM2 and BM6 below the difference between improved and unimproved results is less than 1%.
The matching calculation outlined above is straightforward to generalize to other parameters in the effective theory (13).At O(g 4 ), the results are: Matching the interactions between adjoint fields and the singlet at one-loop level gives: x Couplings between S and the temporal gluon C 0 do not appear at one loop.
The remaining parameters m D , m ′ D , m ′′ D , h 3 , h ′ 3 , h ′′ 3 , ω 3 receive no contributions from the singlet at this order and the results can be read from eqs. (A10)-(A12) and (A18)-(A20) in ref. [55] (by setting N t → 0, a 2 → 0 in those expressions). 7For the Debye masses a one-loop matching is sufficient, as they contribute to the Higgs effective potential only through loops.
In practice, we have automated the matching calculation using in-house FORM [56] software based on a standard Laporta algorithm [57] for the reduction of sumintegrals.All required integrals can be read from [49], and [50] collects them in a reduced form.While installing a general covariant gauge, we have verified that all matching relations are manifestly gauge invariant and independent of the matching scale up to higher-order corrections.For the ensuing numerical analysis we fix the matching scale μ = 4πe −γ T ≈ 7T , for which the L b terms vanish identically and which should be optimal for minimizing higher-order logarithms [36].
For later convenience, we also match a set of dimension five and six operators (in 4d language).These are used merely to estimate the accuracy of our EFT and are excluded from the numerical analysis below.The effective potential is most sensitive to non-derivative operators constructed from φ and S. Up to dimension six, these are where the coefficients are to be matched at one-loop.In the SM, the operator c 6,0 (φ † φ) 3 is dominated by contributions from the top quark [37].Here we include the top quark loop and pure scalar contributions at O(g 6 ) in our power counting.Corrections from gauge fields are neglected, because we expect dominant effects to come running is exact up to terms of order ω 2 3 ∼ g8 .There are no further corrections to the running at higher orders. 7Eq.(A12) in [55] from the scalar sector. 8At leading order in high-T expansion, the matching reads Since the Yukawa and BSM contributions to c 6,0 come with opposite signs, the error from neglecting the (φ † φ) 3  operator can be smaller than in the minimal SM.

VI. INTEGRATING OUT THE DEBYE SCREENED GAUGE FIELDS
As a final step of our EFT construction we integrate out the adjoint fields A 0 , B 0 , C 0 , incorporating additional ring resummations at the soft scale.These are not typically included if calculating the potential directly in 4d.This procedure is possible because the scale gT of electric Debye screening is parametrically larger than that of magnetic screening, g 2 T .The resulting EFT is simply where the scalar potential is as in ( 14) but with modified parameters, which we denote by an additional overline (for instance m 2 φ,3 → m2 φ,3 ).Likewise, the gauge couplings are ḡ3 , ḡ′ 3 .Again, the notation is left unchanged for fields.
Matching proceeds in complete analogy to that of the dimensional reduction step.Because the singlet couples to the adjoint fields only at O(g 4 ), corrections to couplings ā1,3 , ā2,3 , b3,3 , b4,3 start at O(g 5 ) and can be dropped.For S and SS correlators at two-loop, we likewise neglect diagrams containing singlet-adjoint vertices.This leaves diagrams where the interaction with adjoint scalars occurs through the Higgs, but these are proportional to m D /m φ,3 and are absorbed by one-loop diagrams involving a counterterm in the IR theory (this is analogous to the cancellation of mixed hard and soft loops due to a thermal counterterm in dimensional reduction).Thus one-loop matching is sufficient for b1,3 , m2 S,3 : Even these are formally higher-order corrections, but we nevertheless retain them here.
For the Higgs mass parameter, we again neglect all two-loop contributions involving singlet-adjoint couplings.Then the result for m2 φ,3 is the same as in the SM: The two-loop part can be read from equation (B.97) of [53].In our power counting the SU(3) term ω 3 m ′′ D is O(g 5 ); the reason we include it here is because ω 3 ∼ g 2 s y 2 t is numerically large and affects the mass more than the hypercharge contribution.The remaining parameters ob-tain no BSM corrections at one-loop: Finally, corrections to the higher-dimensional operators in eqs.( 43)-( 49) are numerically small compared to the hard scale contributions [37] and will be neglected in our error estimates.This completes the construction of our EFT.

VII. THE EFFECTIVE POTENTIAL
We are now ready to construct the effective potential V eff using the EFT (50), with resummation of hard thermal loops implemented by the 4d → 3d matching relations ( 29)-( 41) and the Debye screened adjoint fields integrated out.The effective potential can be obtained by integrating over all non-zero momentum modes in the path integral around constant background fields.In 3d this is essentially a zero-temperature calculation: temperature dependence appears only in the definitions of the renormalized 3d parameters.
We shift the Higgs field in analogy to eq. ( 3), denoting the associated background field in 3d units by v, and the singlet is shifted as S → S + x.The background fields v and x are treated as free real parameters.As in the zero-temperature case, the neutral components of φ and S mix to form eigenstates h 1 , h 2 , but now the mixing angle depends on the background fields.In the diagonal basis, obtaining the one-loop correction to V eff is straightforward [36].Including the gauge fields, the 3d effective potential up to one-loop order is where d = 3 − 2ǫ, and we have used Landau gauge ξ = 0.
The integral appearing at one-loop is UV finite and given by Masses of the gauge and would-be Goldstone fields are while mh,1 and mh,2 are given in Appendix B. The two-loop correction to V eff requires more work but is straightforward, as there is no need for further resummations and the required integrals are known [29].The correction is obtained from two-loop vacuum diagrams in the shifted theory and is UV divergent, unlike the oneloop correction.Details of the calculation are relegated to Appendix B.
A qualitative difference to 4d effective potentials is that in 3d, RG evolution starts only at two-loop order and can be solved exactly, with no corrections at higher or-ders (ignoring the running of field-independent additive terms).The only parameters requiring renormalization are b1,3 , m2 φ,3 , m2 S,3 , whose values depend on the 3d renormalization scale μ3 , see section V. Accounting for their running in the tree-level part yields a manifestly RGinvariant result for the two-loop V eff [36].In practice, we use the RG-evolved parameters even inside the loop corrections.This results in residual μ3 -dependence that would be cancelled by logarithms at 3-and 4-loop orders.
The vacuum structure can be studied by minimizing the effective potential perturbatively around the treelevel minima, i.e. v = v0 + v1 + 2 v2 + . . .where is a formal loop-counting parameter (and similarly for x).This " expansion" gives gauge-invariant results for V eff in its minima, order-by-order in [58,59].Unfortunately, this construction can fail near phase transitions [44].For instance in the Z 2 symmetric xSM, the critical temperature for "symmetric → Higgs" type transitions is determined, at tree level, by the condition m2 φ,3 (T c ) = 0.In the expansion one needs the two-loop correction evaluated at the tree-level minimum, i.e. at vanishing m2 φ,3 .This leads to an explicit IR divergence in the two-loop potential, invalidating the method near T c [29,60].
Such divergences can be avoided by solving the equations ∂V eff /∂v = ∂V eff /∂ x = 0 "exactly", in practice by minimizing the potential numerically.This sums a subset of one-particle-reducible corrections that render the potential IR-finite at its minima, and the difference to a strict expansion is formally of higher order [44].However, the exact minimization generally introduces residual gauge dependence to the results that is only cancelled at higher loop orders.
In this paper, we will take the latter approach and minimize V eff numerically in Landau gauge.Our reasons for not implementing a gauge-invariant expansion are as follows: 1. We will study the xSM near the Z 2 symmetric limit, where the tree-level condition for EWPT is m2 φ,3 (T c ) ≈ 0 and performance of the expansion is questionable.

For small gauge parameters |ξ| <
∼ 1, one expects the remaining gauge dependence to be numerically comparable to the gauge-independent corrections from higher loop orders.Landau gauge ξ = 0 should minimize this effect, and earlier studies in other models support this assumption [61,62].

Many existing studies in the xSM apply Landau
gauge, and since our goal is to investigate the convergence of perturbation theory in a typical BSM setting it makes sense to present our results in the same gauge.
4. For strong transitions in the xSM, we expect that dominant uncertainties come from residual RG scale dependence in scalar loops and not from gauge loops where the gauge dependence appears.This expectation is supported by the one-loop study of ref. [21].
It is worth noting that the residual gauge dependence arises only in the effective potential, while our 4d → 3d matching relations are gauge invariant.Of course, while the gauge-dependent approach avoids explicit IR divergences at two-loop, it does not remove the fundamental IR problem of high-T perturbation theory, which is an intrinsic property of the Matsubara zero modes.There can still be divergences at higher orders due to massless gauge bosons in the symmetric phase, and slow convergence is to be expected if the scalars are light compared to the dimensionful couplings in the 3d scalar potential.

VIII. NUMERICAL RESULTS
To identify phase transitions, we look for temperatures where the effective potential has degenerate local minima.For a set of input parameters M h2 , sin θ, b 3 , b 4 , a 2 and the temperature T , the analysis consists of the following steps: 1. Solve the renormalized parameters at MS scale μ = M Z from the one-loop corrected relations (A27)-(A35) in Appendix A.
2. Run the parameters to the matching scale μ = 4πe −γ T using one-loop β-functions (c.f.Appendix A).This ensures that any leftover μ dependence in the EFT matching is of higher order than the calculation.
4. Obtain the two-loop effective potential from eq. ( 56) and Appendix B and minimize V eff (v, x) with respect to the background fields using a differential evolution algorithm.
Steps 2 through 4 are repeated for different temperatures until a critical temperature T c is found (we vary T in steps of 0.1 GeV).Outside the minimum V eff can develop an imaginary part, signaling an unstable field configuration [63].In such cases, we simply discard the imaginary part.
At T c , we calculate the jumps ∆v and ∆x in scalar VEVs across the phase transition, and the latent heat L = T ∆(∂V eff /∂T ).To avoid confusion, we use familiar 4d units when discussing the results: i.e.V eff and the fields in 3d get rescaled by approriate powers of T .We also compute the discontinuity in the quadratic Higgs condensate, ∆ φ † φ = ∆(∂V eff /∂ m2 φ,3 ), which provides a gauge-invariant definition of the Higgs VEV through v phys = 2∆ φ † φ [32] (although our results contain residual gauge dependence, In BM10 the 2-loop results (solid curves) also show abrupt restoration of the Z2 symmetry associated with the singlet around T ≈ 130 GeV, while at 1-loop the restoration occurs at much higher T ≈ 300 GeV, not shown here.
see above).The difference between ∆v and v phys is small in all of our benchmarks.We do not study metastability ranges of the phases nor nucleation processes, so our latent heat will be given at T c instead of the nucleation temperature.Therefore our results give a lower bound on the actual transition strength.
We focus mainly on strong transitions with ∆v/T c > ∼ 1.0 and examine how the transition changes as the twoloop correction is included.In searching for strong transitions we make use of the one-loop parameter scans appearing in [19,35].In this paper we concentrate our discussion on selected benchmark (BM) points that represent typical strong-EWPT scenarios, rather than fullfledged scans of the free parameter space.This approach suffices to quantify the effect of two-loop corrections on basic quantities like ∆v/T c , but it would be interesting to perform systematic two-loop scans over the free parameters in a future study.
Table I collects our results for direct transitions between the EW symmetric and broken minima.In the non-Z 2 symmetric xSM, the minima are separated by a potential barrier already at tree level due to the singlet VEV being non-zero (points BM1-BM5), while in the Z 2 symmetric limit the barrier is generated radiatively (BM6 and BM7).In the Z 2 symmetric case there are also twostep transitions where the EWPT is preceeded by spontaneous breakdown of the Z 2 symmetry, which is restored in the final EW minimum.Here we consider only the second, EW-breaking transition, for which our results in the two-step scenario are in table II.Because the Z 2 -breaking phase remains metastable even at zero temperature, a realistic two-step scenario necessitates rapid enough nucleation into the EW phase to avoid contradictions with standard cosmology.While this requirement can restrict the parameter space of allowed two-step EWPT considerably [19], accounting for it requires calculation of the bubble nucleation rate and will not be considered further in this paper.Evolution of the background fields in the different types of EWPT are illustrated in fig. 2.
For comparison, we show also the results from a oneloop calculation.Here "one loop" refers to the one-loop potential in (56) with the EFT parameters matched at O(g 2 ), in other words it is the daisy-resummed 4d potential to leading order in high-T expansion (with an additional resummation of A 0 loops), correct at O(g 3 )."Two loop" refers to the full O(g 4 ) calculation described in detail in the above sections and in appendix B. In both cases we use the same T = 0 renormalization prescription for input parameters.
With the exception of benchmark points 10 through 12, we find lower critical temperatures at two-loop order and the difference to the one-loop temperature is quite large in almost all points.We expect this shift to originate from the O(g 4 ) corrections to thermal masses.For one-step transitions with a tree-level barrier, we observe that transitions tend to be stronger at two loop, and the latent heat is increased by at least 20% relative to the one-loop result.In contrast, two-loop corrections make the transition significantly weaker in the Z 2 symmetric xSM.Here, the transition may appear to be strongly first order at one loop, but the radiatively generated potential barrier can be destroyed by higher-order effects if the perturbative expansion convergences slowly.In fact, scanning over M h2 ∈ [150, 450] GeV and a 2 ∈ [1, 5] with b 4 = 0.3 in the Z 2 symmetric case, we found no transitions that would remain strong (∆v/T c > 1.0) at the two-loop level.Indeed, one-loop results should be taken with a grain of salt when the coupling a 2 is this strong.These considerations call for more thorough investigations of higher-order effects in the Z 2 symmetric xSM in the future.Two-step transitions are quite sensitive to the singlet self-interaction and here we concentrate on b 4 = 0.3.Contrary to naive expectations, we find the transition strength to be very sensitive to two-loop corrections, despite the presence of a tree-level barrier.Particularly in BM12 the latent heat is an order of magnitude smaller  Points BM6 and BM7 are in the Z2-symmetric xSM.δEFT estimates the validity of our dimensionally-reduced theory as described around eq. ( 60).The numbers are obtained using fixed RG scale μ3 = T .
than the one-loop value, signaling definite breakdown of perturbation theory at least at low orders, probably due to the large value of a 2 .For points BM8 and BM9 where the coupling is smaller, two-loop effects are moderate but quantitatively important.We conclude that transitions involving a tree-level barrier in the potential are not protected from large loop corrections.The presence of such a barrier implies that loop effects from the Matsubara zero modes are relatively less important than in radiativelyinduced transitions.This argument does not apply to loop corrections involving hard Matsubara loops, whose effect is to modify the long-distance coefficients of the zero-mode EFT.Large loop corrections for two-step transitions were also reported in [29] in a triplet-extended model.
Our analysis was performed at a fixed renormalization scale μ3 = T and it is important to estimate how sensitive the results are to variations of this scale (as discussed above, the remaining scale dependence is of higher order than the calculation at hand).In our EFT approach the soft IR modes evolve according to the exact RG equations in 3d, so we may expect weaker scale dependence than what is observed in 4d effective potentials.In fig. 3 we plot the effective potential at different values of μ3 in three benchmark points, together with the one-loop potential for comparison.The scale dependence is seen to be extremely mild when the transition occurs through a tree-level barrier.For the radiatively generated transition in BM6 the scale variations can shift the critical temperature by a few GeV, but the transition itself is so weak that perturbation theory may not be accurate due to non-perturbative gauge loops.
What is not shown in fig. 3 is possible sensitivity to the matching scale at which the EFT is constructed, corresponding to residual scale dependence of the hard thermal corrections.By varying the matching scale in the interval [5T, 9T ] we have checked that this scale ambiguity affects the two-loop results even less than the residual dependence on the 3d scale μ3 .In contrast the one-loop is quite sensitive to the matching scale, which is not surprising since we already concluded that twoloop corrections are non-negligible.Dependence on the renormalization scale is further discussed in [64].
The validity of our high-T prescription is estimated as follows.After minimizing the two-loop potential at a given T , we add to tree-level V eff the dimension five and six operators from eq. ( 42) that capture the dominant higher-order corrections to our high-T EFT.We then locate the minimum once more using the modified potential and look at the relative change in the Higgs VEV, This quantity estimates the effect of higher-dimensional operators, that were neglected in the main analysis, on properties of the broken phase and particularly on the transition strength.The tables show values of δ EFT at T c , where δ v is evaluated on the v = 0 side of the transition.We also considered the shift in the singlet VEV, but in all of our points it was smaller than the respective change in v.In most of our benchmarks the relative error is small, less than one percent.Hence our 3d EFT should give a very accurate description of at least the broken phase.
The excellent performance of high-T approximations has also been demonstrated in [27] for a two-Higgs doublet model.The above accuracy considerations suggest that the uncertainty in our two-loop results should be at most of order 10% and considerably less in most points.One exception is BM12, where the a 2 coupling runs rapidly to non-perturbative values.Of course, the simple arguments presented here cannot conclusively rule out the possibility of sizable corrections at higher loop orders or at the non-perturbative level.Considering how much the two-loop terms affect T c and L/T 4 c in several points, we may even anticipate large higher-order effects in at least some region of the interesting parameter space.

IX. DISCUSSION
In the paper at hand, we carried out a two-loop study of the electroweak phase transition in the singletextended Standard model (xSM), quantifying the magnitude of two-loop corrections on equilibrium characteristics of the transition.The study involved radiativelyinduced transitions and transitions through a tree-level barrier, as well as two-step phase transitions.In all cases we found considerable deviations from one-loop predictions, in particular the critical temperature and latent heat typically change by at least 20%, but particularly for very strong two-step transitions the differences can be 100% or more.
Although we focused on a limited selection of points in the free parameter space of the xSM, the main conclusion should be clear: one-loop analyses are generally not very accurate in predicting the strength nor the interesting temperature range(s) of the EWPT.In our opinion, this result by itself is not very surprising, because strong transitions in the xSM are typically associated with at least O(1) coupling constants in the scalar sector, and because the perturbative expansion converges more slowly at high temperature than it does at T = 0. We consequently anticipate that two-loop corrections can significantly impact predictions of, for example, the bubble nucleation rate and other quantities relevant for gravitational waves and baryogenesis.Our results strongly motivate improved studies of the EWPT and its cosmological consequences at the two-loop level, not only in the xSM but also in other scalar extensions.
On the technical side, we wish to emphasize that the dimensionally-reduced approach to the EWPT described in this paper provides a systematic and very intuitive approach to thermal resummations beyond one-loop level.It is also quite straightforward to apply the method to other BSM models and in fact, dimensionally-reduced investigations have already been performed in several such models [25,28,29].The approach is limited to theories in the formal high-T limit, but tests of its validity indicate good accuracy even when there are relatively heavy excitations in the spectrum [27][28][29].In the present paper, we found dimensional reduction accurate also in the xSM.This observation motivates one more approach of studying phase transitions in the xSM, namely that of non-perturbative 3d simulations, which allow for a very accurate determination of the relevant thermodynamics at the cost of increased numerical effort.
Note added on March 2024.In previous versions of the paper there was an error in our numerical implementation of eq.(A28) which affected our benchmark results somewhat.This has been corrected in the present version with figures and tables updated accordingly.In addition, we have corrected misprints in equations (38) and (39).See the published Erratum [65] for further details.
The SU(2) L gauge coupling can be fixed by relating it to Fermi constant, which describes muon decay in the lowenergy theory obtained by integrating out weak interactions.Conventionally, the matching relation between the full EW theory and the Fermi EFT is written, in the on-shell scheme, as ∆r contains contributions from weak gauge bosons and scalars, but not from the QED sector which cancels when performing the matching.The on-shell (OS) coupling is related to our MS gauge coupling by g 2 (μ) = g 2 os (1+δg 2 os /g 2 os ), where δg 2 os is the regularized OS counterterm (this relation follows because the bare coupling must be equal in both schemes).
Sum is over all fermions.Finally, for the top quark we have  The integrals D and I 3 1 can be found in the Supplemental Material of ref. [29].In the above, c ψi...ψj denote vertex coefficients for generic fields ψ i and are defined as minus the respective coefficient in the Lagrangian, including the combinatorial factors arising from contractions.For momentum-dependent vertices the momentum is absorbed inside the integral definition.An exhaustive list of required vertex coefficients read

FIG. 2 .
FIG.2.Evolution of Higgs and singlet VEVs, v and x respectively, in different types of transitions: (a) symmetric → broken transition in the presence of a tree-level barrier, (b) radiatively generated symmetric → broken transition in the Z2 symmetric case, (c) two-step phase transition with spontaneously broken Z2 symmetry.The dashed lines are predictions from a one-loop calculation.In BM10 the 2-loop results (solid curves) also show abrupt restoration of the Z2 symmetry associated with the singlet around T ≈ 130 GeV, while at 1-loop the restoration occurs at much higher T ≈ 300 GeV, not shown here.

FIG. 3 .
FIG. 3.Real part of the effective potential, plotted in the Higgs direction at one and two loops at their respective critical temperatures.In (a) and (c) we have accounted for the changing singlet VEV.In (b), the potential barrier vanishes almost completely at two-loop level.The dashed and dotted lines correspond to variations of the RG scale μ3.The potential generates an imaginary part near the local maximum as a scalar eigenstate becomes tachyonic, and at specific values of v/T there is a massless mode.In (a) and (c), this results in "kinks" in the two-loop potential, since 3d perturbation theory is not well behaved at such points.This complication does not affect our numerical results, which depend only on properties of the potential in its local minima.

TABLE I .
. The potential generates an imaginary part near the local maximum as a scalar eigenstate becomes tachyonic, and at specific values of v/T there is a massless mode.In (a) and (c), this results in "kinks" in the two-loop potential, since 3d perturbation theory is not well behaved at such points.This complication does not affect our numerical results, which depend only on properties of the potential in its local minima.Benchmark scenarios for direct, one-step EWPT.

TABLE II .
Collected results for two-step EWPT in chosen benchmark points.We consider only transitions in the Higgs direction.