Three-dimensional effective theories for the two Higgs doublet model at high temperature

Due to the infrared problem of high-temperature field theory, a robust study of the electroweak phase transition (EWPT) requires use of non-perturbative methods. We apply the method of high-temperature dimensional reduction to the two Higgs doublet model (2HDM) to obtain three-dimensional effective theories that can be used for non-perturbative simulations. A detailed derivation of the mapping between the full four-dimensional and the effective three-dimensional theories is presented. The results will be used in future lattice studies of the 2HDM. In the limit of large mass mixing between the doublets, existing lattice results can be recycled. The results of such a study are presented in a companion paper.


Contents
1 Introduction

Electroweak phase transition, gravitational waves and baryogenesis
As the search for new particles at collider experiments continues, the full structure of the scalar sector remains an active subject of theoretical study. Beyond-the-Standard-Model (BSM) theories assuming a richer Higgs sector are motivated by unanswered phenomenological questions in the Standard Model (SM) and also by cosmological observations suggesting, e.g., the existence of an unknown dark matter particle. One particularly interesting BSM theory is the Two-Higgs-Doublet Model (2HDM), which augments the SM with an additional scalar doublet and predicts new bosons that could have observable signatures at present particle accelerators [1][2][3]. In particular, perturbative studies of the 2HDM at high temperatures suggest that it may be possible to explain the observed matter/antimatter asymmetry by means of electroweak baryogenesis (EWBG) [4,5]. EWBG is a mechanism for generating an excess of baryonic matter during the electroweak phase transition (EWPT) via non-perturbative sphaleron processes near the bubble walls that form during a first-order phase transition [6]. However, it is widely known from lattice simulations performed in the 1990's that the EWPT in the SM with a physical Higgs mass of 125 GeV is a cross-over transition instead of first order, ruling out EWBG in the minimal SM [7][8][9][10]. It has also been demonstrated that an another necessary ingredient for baryogenesis, CP violation, is too weak in the SM [11][12][13]. However, models with multiple scalar doublets provide a mechanism for CP violation beyond that of the CKM matrix via a mixing term between the doublet fields [4,14], making the 2HDM a viable candidate for the realisation of EWBG. A thorough study of the nature of the EWPT in the 2HDM could thus provide insight on both the phenomenology of the model and the cosmological question of matter/antimatter asymmetry.
First-order phase transitions at the electroweak scale are also a source of gravitational waves, peaked at a characteristic frequency given by the bubble radius, which would be in the mHz range today [15,16]. This is within the sensitivity regions of LISA, so if the phase transition were strong enough, its existence and properties could be probed through the gravitational-wave power spectrum it left behind [17]. Studies of gravitational waves from the 2HDM have been carried out in the past in [18]. With the results of this paper, we will aim to improve the precision of these investigations.

Non-perturbative analysis
Frequently, properties of the EWPT are studied in terms of a finite-T effective potential; that is, perturbatively in the full four-dimensional theory [19][20][21][22][23][24][25]. This approach contains a serious disadvantage: It is well-known that in perturbation theory the symmetric phase is accociated with disastrous infrared (IR) problems [26,27]. Yet, in order to to find the critical temperature in perturbation theory-from the condition that the symmetric and broken minima are degenerate-information about the value of the potential at the broken minimum, as well as the value of the potential at the origin, is required. While the former quantity can be determined, as at sufficiently large field values perturbation theory is applicable, the latter quantity cannot be computed due to the non-perturbative nature of the symmetric phase. This means that an accurate determination of the critical temperature of the phase transition-as well as some other thermodynamic quantities-is beyond the scope of perturbation theory. The determination of these quantities in perturbation theory is always inaccurate without information about the behaviour of the potential near the origin. For this reason, reliable determination of these quantities requires use of non-perturbative methods: in practice, lattice Monte Carlo simulations. In the non-perturbative studies of the SM in [7], considerable deviations from perturbative computations of the effective potential were found at small field values.
In the paper at hand, we take a different approach from earlier perturbative studies and seek to apply existing non-perturbative results. We first carry out a procedure known as finite-temperature dimensional reduction, explained in detail in [28], to construct an effective three-dimensional theory at high temperature. This theory is then readily studied non-perturbatively using previously-obtained lattice results [7,29]. By taking the effective theory approach, one is able to construct a lower-dimensional field theory that reproduces the long-distance physics of the full theory without suffering from IR problems [28,30,31]. Physically, this dimensional reduction is made possible by the fact that in thermal equilibrium, the four-dimensional fields can be expressed in terms of three-dimensional Matsubara modes that generate thermal masses proportional to πT , a scale which we shall refer to as "superheavy". This causes all non-zero modes to decouple from physics at high temperatures. The remaining fields in the effective theory are just the bosonic zero modes.
In practice, DR is performed by matching parameters of the three-dimensional theory to those of the full theory so that the long-distance Green's functions match. This requires perturbative calculations of correlation functions in the four-dimensional theory to a given accuracy. We will perform the analysis at O(g 4 ) accuracy, requiring one-loop accuracy in the four-point functions and two-loop accuracy in the scalar two-point functions. This matching procedure can, in principle, be performed at arbitrary accuracy, but it requires one to include higher-dimensional operators beyond order O(g 4 ).
Temporal components of the gauge fields are treated in the effective theory as additional scalar fields with masses of order gT : a mass scale we shall refer to as "heavy". For the 2HDM, a diagonalisation of the scalar potential in the three-dimensional theory shows that in the presence of a mixing term, one of the doublets will always be heavy in the vicinity of the critical temperature. We will proceed to construct the theory at "light" scale g 2 T in two steps: First, we integrate out the temporal scalar fields, and second, we integrate out the heavy doublet, using matching both times. This allows us to construct a final effective theory that contains only the light bosonic (zero) modes. Information about the second doublet and the properties of the full theory are contained in the parameters of this final effective theory, which turns out to have the same form as the one derived from the SM in [28]. Existing lattice results can thus be used to make conclusions about the nature of the EWPT in the 2HDM.
We emphasize that since the resulting final effective theory will be super-renormalizable, non-perturbative simulations in the three-dimensional theory can reasonably be performed (see [32]) in contrast to performing direct simulations in the full four-dimensional theory. In addition, combining the perturbative method of DR with non-perturbative threedimensional simulations is very efficient, as the DR procedure is free of IR problems and can be performed accurately using perturbation theory, while the latter is used to study the dynamics of the light bosonic modes, which are the source of IR problems in perturbative studies. Thus, this approach automatically implements the required resummations needed at finite temperature, and therefore mitigates the IR problems, unlike the conven-tional four-dimensional approach with a finite-T effective potential. Finally, we highlight that to the order O(g 4 ) that we work, the parameters of the effective theory-in terms of T and other physical quantities-are independent of the choice of gauge. Thus, with this three-dimensional approach, gauge invariance can be seen to hold explicitly, unlike in many analyses using the perturbative effective potential. For a discussion of gauge invariance in a perturbative study of the EWPT, see [33].
Dimensional reduction has previously been used to study the EWPT in the SM [28,34], MSSM [29,35,36], as well as the 2HDM in Ref. [37,38]. Technical details presented in the appendices of this work have also recently been used in Ref. [39], where dimensional reduction has been applied to the SM augmented with a real triplet, previously studied perturbatively in Ref. [40]. Similar techniques are currently being applied to the realsinglet extension of the SM as well [41], while this model has already been analysed-in limited regions of parameter space-with the three-dimensional approach in [42,43] (for perturbative analyses of the EWPT in the singlet extension, see [44,45]). For a compact and illuminating review on the use of dimensional reduction in studies of the EWPT, see [46].
Despite the fact that dimensional reduction and a non-perturbative analysis of the EWPT have been succesfully applied to the SM and MSSM, these methods are not widely used for BSM models with an extended scalar sector. Our hope is that this detailed work will make the DR method more transparent, as it is a valuable alternative to the widely used, fully perturbative method. This paper is organised as follows. In section 2 we introduce the model in Euclidean space, introduce our conventions and present the form of the effective three-dimensional theories. In section 3 we collect together the main results of this paper: the matching relations between the full theory and the effective theories. In addition, we describe how the four-dimensional MS parameters in the Lagrangian can be related to physical quantities. In section 4 we summarise key aspects of our study and outline future prospects. Phenomenological implications of our results are discussed in a companion paper [47]. Many details of the derivation of our results are relegated to the appendices.

2HDM and three-dimensional effective theories
In this work, we work in a D = d + 1 = 4 − 2ǫ dimensional Euclidean spacetime.

Full four-dimensional theory
The Lagrangian of four-dimensional theory reads L = L gauge +L ghost +L gauge fixing +L fermion +L scalar +L Yukawa +δL +L resummation , (2.1) where the gauge field, ghost, fermion, scalar and Yukawa sector Lagrangians are defined as follows: In the gauge sector we have the SU(2) L , U(1) Y and SU(3) c gauge fields A a µ , B µ , and C α µ appearing inside the field strength tensors G a µν , F µν and H α µν . The associated gauge couplings are g, g ′ , and g s . The only ghost field of relevance is the SU(2) L ghost η a . Left-handed doublet and right-handed singlet lepton fields are denoted ℓ A and e A , with A being the flavour index, while q A refers to left-handed doublet quark fields. u A and d A are right-handed singlet up-and down-type quark fields, respectively. The scalar sectors consists of the doublet fields φ i 1 , φ i 2 and the corresponding charge-conjugated fields where τ 2 is the second Pauli matrix. Finally, in the Yukawa sector we use the approximation that the top quark is the only fermion that couples to the scalar doublets, and furthermore declare that it only couples to one of them, namely φ 2 . The relation Q = I 3 + Y 2 between electric charge Q and isospin I 3 defines the hypercharge of the fields as follows: We work in Landau gauge, as this choice significantly simplifies many diagrammatic calculations. We emphasize, however, that the final matching results we present are gauge independent to the order we are working. This can be verified by performing the calculation in general covariant gauge and explicitly verifying the cancellation of the gauge parameter between terms from field normalisation and correlators. In particular, the effective theory approach allows us to bypass the issues related to gauge dependence that are present in perturbative analyses relying on the use of the finite-T effective potential [33,48].
Consistent construction of the effective theory requires thermal resummation in order to remove problematic contributions originating from two-loop integrals with mixed Matsubara n = 0 and n = 0 modes [28,49]. We implement this in L resummation by adding and subtracting one-loop thermal masses, denoted byΠ, as well as a thermal mixing mass term Π 12 , to zero modes of the scalar fields φ 1 , φ 2 . Schematically where the Euclidean four-momentum is defined as P = (ω n , p), where ω n = 2nπT , and we use the notation m ≡ √ m 2 +Π. Temporal components of the gauge fields are treated similarly; their thermal masses are just the corresponding Debye masses. Terms with minus signs, −Π, are treated as counterterm-like interactions; hence they are called thermal counterterms, despite being UV-finite. Thermal masses are listed explicitly in Appendix C.3. Terms with one-loop resummed masses µ 2 and +m 2 D , +m ′ D 2 contribute to propagators.
This procedure is done both for the doublets φ 1 , φ 2 and the temporal scalar fields A a 0 , B 0 . The temporal gluon field C α 0 does not require resummation at order O(g 4 ). The scalar potential reads: where the parameters µ 2 11,22 , λ 1,2,3,4 are real and µ 2 12 , λ 5,6,7 are, in general, complex. Perturbative expansions of correlation functions, required for dimensional reduction, are organised in terms of the gauge coupling g. We assume that all mass parameters are heavy, i.e. scale as µ 2 ∼ g 2 T 2 , and that all quartic couplings scale as λ ∼ g 2 . Gauge couplings g, g ′ , g s , as well as the top quark Yukawa coupling g Y , are assumed to scale as g.
The Lagrangian can be simplified by imposing a Z 2 symmetry. An exact Z 2 symmetry requires ρ = λ 6 = λ 7 = µ 2 12 = 0, while a soft violation of Z 2 symmetry is achieved with ρ = λ 6 = λ 7 = 0, but µ 2 12 = 0 (see discussions in [50,51]). In Ref. [52] it is described how a treatment of a true hard violation of Z 2 symmetry is inconsistent without kinetic mixing terms with complex coupling ρ. In spite of this, we perform the dimensional reduction following Refs. [37,38] and set ρ = 0 without imposing the full Z 2 symmetry; rather, we keep λ 6 and λ 7 in our calculation. However, when turning to numerical analysis, we restrict ourselves to the case of soft violation of the Z 2 symmetry.
Renormalisation of the four-dimensional theory is discussed in Appendix C.2.

Three-dimensional effective theories
We denote the fields of the effective theories with the same symbols as those of the fourdimensional theory, even though their normalisation is different and will affect the mapping between the full and effective theories. For a generic field, the relation between the fourdimensional and three-dimensional fields reads where Π ψ (P ) is the self-energy of the field, a prime denotes a derivative with respect to P 2 , and δZ ψ is the field renormalisation counterterm. The effective-theory gauge couplings are denoted by g 3 and g ′ 3 . The Lagrangian of the effective theory (again in Landau gauge) then has the schematic form temporal + δL (3) .
We include the SU(2) L and U(1) Y gauge fields in the gauge sector part, where only spatial Lorentz indices are summed over. The spatial SU(3) c gluon fields can be neglected. The counterterm part δL (3) plays an important role in determining relations between the continuum and lattice three-dimensional theories and is needed for the calculation of lattice counterterms [53]. In a continuum three-dimensional theory with dimensional regularisation, one-loop correlators are finite, while two-loop contributions to self-energies are divergent. The three-dimensional theory is super-renormalisable, and the from two-loop mass counterterms one can solve the exact running of the mass parameters in terms of the three-dimensional theory renormalisation scale Λ 3 . The form of L (3) scalar is the same as in the four-dimensional theory, but we denote the couplings with an additional subscript, emphasising that they are couplings of a threedimensional theory. The temporal scalar contribution reads Here the covariant derivative of an isospin triplet is D r A a 0 = ∂ r A a 0 + g 3 ǫ abc A b r A c 0 and for the temporal gluon field we have used the usual derivative instead of the covariant deriva- from the effective theory. Spatial gluons do not couple to the scalar fields, and self-interactions of temporal gluons and their interactions with other temporal scalars would have a very small contribution to quantities of interest, such as scalar mass parameters of the light scale effective theories.
Furthermore, since the mass parameters can be close to zero near the phase transition, contributions of the type 1/µ 2 need to be considered carefully. These appear in two-loop calculation of scalar two-point correlators. In order to cancel potential IR divergences, we apply a procedure analogous to the thermal resummation in the four-dimensional theory by adding and subtracting one-loop corrections from temporal scalar fields to fundamental scalar masses. Terms with plus signs contribute to the masses in scalar propagators, while terms with minus signs are treated as (counterterm-like) interactions, i.e.
where m ≡ √ m 2 +Π. Note that we do not need to include a counterterm interaction for the mixing mass parameter µ 2 12,3 , as the one-loop correction from the temporal scalar fields is of higher order. Explicit expressions for these mass corrections are given in Appendix C.3.
After the heavy temporal scalars have been integrated out, their effect is encoded in the parameters and fields of a new theory, and we havē where the parameters are denoted with a bar asḡ 3 ,ḡ ′ 3 ,μ 2 11,3 , etc. In this theory, the phase transition takes place close to the point where the mass matrix has a zero eigenvalue, and in the diagonal basis the other mass parameter is then generically heavy. By performing a unitary transformation (see Appendix A), one can remove the mixing mass term, and the resulting theory is given by where the scalar potential reads and φ and θ are light and heavy, respectively. Note that in general the diagonalisation procedure generates non-zero couplings λ 6 and λ 7 even in the case of a softly broken Z 2symmetry.
The heavy doublet θ can be integrated out in a similar fashion to the temporal scalars, and we obtain a final effective theory, which has the same form as the effective theory constructed for the SM in Refs. [28,34]: where the couplings are denoted with a hat asĝ 3 ,ĝ ′ 3 , and (2.14) Note that this method of three-step DR is analogous to that of Ref. [35] in the MSSM. Couplings are RG invariant, and the mass parameter runs at two-loop order. Due to super-renormalisability, the running ofμ 2 can be solved exactly from two-loop mass renormalisation, and the corresponding β function receives no additional corrections at higher loop orders. In certain regions of parameter space, it is possible for both doublets to be light in the vicinity of the electroweak phase transition, in which case the final effective three-dimensional theory is given by Eq. (2.10). In a non-perturbative study of that theory, both doublets would need to be treated dynamically on a lattice. We leave this case for a future study, however, and in this work we limit our analysis to the regions of parameter space where the second doublet is heavy and can be integrated out. In this case, we use the DR matching relations that map the four-dimensional theory to the effective threedimensional theory of Eq. (2.13), and recycle the existing non-perturbative results of [7]. Note that in this work, non-perturbative effects from the U(1) gauge field were neglected. For simplicity, we follow this reference; a non-perturbative analysis with U(1) field is presented in [54] and shows no significant difference from the case where only the SU(2) field is considered.
In the final effective theory, one of the four parametersĝ 3 ,ĝ ′ 3 ,μ 2 ,λ can be used to measure all the dimensionful quantities as well as to fix the RG scale for the mass parameter (the couplings are RG invariant). We follow Ref. [7] and chooseĝ 3 . Then, the dynamics of the system is determined by three dimensionless ratios, x ≡λ (

2.15)
Properties of the phase transition, however, depend essentially on only one parameter: As a justifiable approximation, non-perturbative effects of the U(1) gauge field can be neglected by setting z = 0, and in practice, y ≈ 0 on the critical line, close to its leading order value. This means that the character of the transition is described only by magnitude of the parameter x. Results of the Monte Carlo simulations [7] show that for a first-order transition, 0 x 0.11. The transition gets weaker as x increases, and above x ≈ 0.11 only a smooth cross-over remains.
With a DR mapping between the four-dimensional 2HDM and the SM-like effective three-dimensional theory, we can scan the physical parameter space, searching for x < 0.11 and y = 0 to find regions of first order transitions and the corresponding critical temperatures. Note that if x < 0 for some physical input parameters, the three-dimensional theory is not bounded from below and simulations are not possible. This indicates that our DR procedure has broken down, either because of neglected higher order corrections to the matching relations, or neglected dimension six (hereafter 6-dim.) or higher-dimensional operators.

Dimensional reduction
For dimensional reduction, the required ingredients are the matching relations between the full four-dimensional and the effective three-dimensional theories, the one-loop β functions and the relations between MS parameters and physical quantities. In the first step of the dimensional reduction, i.e., when the superheavy scale is integrated out, matching relations are calculated up to O(g 4 ) in our power counting. This accuracy requires one-loop accuracy for couplings and two-loop for mass parameters. However, it suffices to evaluate the Debye masses m D , m ′ D , m ′′ D at one loop only [28]. In the second step of DR, when integrating out the heavy scale, it is convenient and numerically reasonable to perform calculations to same loop order as in the first step of dimensional reduction. One-loop β functions are required to make the matching relations renormalisation scale independent at order O(g 4 ).
For the relations of the MS parameters and physical quantities, we use tree-level relations, despite the fact that for consistent O(g 4 ) accuracy, one should use one-loop corrected relations. This would require performing one-loop vacuum renormalisation of physical quantities, which is a non-trivial task and is left as future work. In the special case of the Inert Doublet Model, similar one-loop vacuum renormalisation calculations can be found in Ref. [55]. The tree-level relations to physical quantities and the matching relations are listed below, while the β functions are listed in the Appendix C.2. In addition, details of the derivation are given in the appendices (see also [56]).
While the analysis of the final three-dimensional theory is non-perturbative and exact within the accuracy of Monte Carlo simulations, the DR procedure is perturbative. Therefore, the validity of perturbation theory at each step of the DR should be estimated. Perturbative errors arise from two sources: Firstly, there are higher order corrections to the parameters of the effective theories. Secondly, 6-dim. and higher-dimensional operators have been neglected in the effective theories. We discuss these higher order operators in Section 3.2.3. Since the first order phase transition requires large values of the couplings [47], it is particularly important to estimate the validity of the DR. For the same reason, we expect the one-loop corrected relations to physical quantites to be of importance.
In contrast to a perturbative treatment with the effective potential, the three-dimensional approach effectively handles the IR behaviour of the finite-T theory: the DR is free from IR problems, while IR physics is correctly captured by the non-perturbative Monte Carlo simulations of the three-dimensional theory. In the DR procedure, by using thermal mass resummed propagators and correspoding thermal counterterms, we are explicitly able to show that at two-loop level, products of zero-mode and non-zero mode contributions in the correlators vanish. Due to this cancellation, one could neglect the effect of the zero modes at two-loop level as only the non-zero modes contribute to the final result. However, keeping the zero modes and explicitly verifying this cancellation serves as a valuable cross-check of our calculations, even though it technically complicates the calculation.

Physical quantities in the full four-dimensional theory
We relate the MS parameters of the 2HDM to physical parameters at tree level, parametrising the complex Higgs fields as . (

3.2)
Here Θ is an angle which simply rotates the field bases (see the discussions in [52] and [57]). In a pure-Higgs theory, the scalar potential and thus the physical Higgs masses are independent of Θ; however, changes to Θ are not a symmetry of the Lagrangian considered in this work, for here the Higgses are coupled to fermions. Despite the lack of symmetry, we may still set Θ = 0, since we are only matching the Higgs masses at tree level, and it is only through loop effects that the parameter Θ affects the masses. The top quark is chosen to couple to φ 2 . Note that our DR is performed using the approximation that the lepton couplings to the scalars can be discarded, and thus no distinction between Type-I and Type-II 2HDM is made.
In this section-and for the final numerical analysis-we explicitly discard the Z 2 hardbreaking couplings λ 6 , λ 7 from our scalar potential of Eq. (2.4). The vacuum expectation values (vevs) are chosen so that they extremize the potential. At tree-level, this is achieved by imposing equations which lead to the following conditions for the mass parameters: where λ 345 ≡ λ 3 + λ 4 + Re λ 5 , and Θ has been set to zero. With this choice of Θ, the vevs v 1 , v 2 are real and constrained experimentally by the relation v 2 where v = 246 GeV is the SM vev. Mixing of the two vevs is parametrised by the angle β, and we use the shorthand notation t β ≡ tan β = v 2 /v 1 . Furthermore, we shall also restrict our analysis to the region of parameter space where λ 5 , and thus µ 2 12 , are real. Re µ 2 12 then defines a natural mass scale of the theory, where we have chosen to simplify the notation by denoting µ 2 ≡ −Re µ 2 12 . Physical states are obtained from the φ ± k , φ k and η k by a diagonalisation and include two CP-even scalars h, H 0 , a CP-odd pseudoscalar A 0 and the charged Higgses H ± . Three of the eight degrees of freedom are absorbed into Goldstone bosons. Mass eigenstates are then related to φ ± k , φ k , η k as Here, α is defined as the mixing angle between the CP-even scalars, and we have introduced the shorthand notation s α , c α , s β , c β ≡ sin α, cos α, sin β, cos β. We shall assume that h is the observed Higgs boson with mass m h = 125 GeV. The quantity c β−α ≡ cos(β − α) is phenomenologically important, as the choice c β−α = 0 corresponds to the alignment limit where h couples to SM particles exactly like the Standard Model Higgs [58]. Physical tree-level masses are found by mass matrix diagonalisation and have been calculated in Refs. [1,58], so we will not list them here. Inverting the eigenvalue relations allows us to write the MS parameters in terms of masses m h , m H 0 , m A 0 , m H ± and mixing parameters t β , c β−α , µ, which are what we input into our parameter space scans. These relations are listed in Appendix B. Precision tests of the 2HDM suggest that m H ± should be close to either m H 0 or m A 0 [59,60], and we choose to set m H ± = m A 0 as this is the region where the perturbative analysis of Ref. [4] mainly found strong first-order phase transitions.
For the gauge couplings and top Yukawa coupling, at tree-level where we have denoted g 2 0 ≡ 4 √ 2G f m 2 W , with G f being the Fermi constant related to the lifetime of the muon.
Tree-level stability [52,58,[61][62][63] and unitarity [64][65][66] requirements set additional constraints on the potential parameters. The relevant equations are listed in Appendix B. Points in the parameter space not satisfying these requirements are considered unphysical and are discarded from our scans.

Matching relations between the full four-dimensional and effective threedimensional theories
In this section, we list the relations between the parameters of the full and the dimensionallyreduced effective theories. The recipe to obtain the matching relations is presented in [28,42] and requires calculation of four-point correlators for couplings and two-point correlators for mass parameters (see Appendix C.1 for computational details). In the companion paper [47], we have listed matching relations in the case of a softly broken Z 2 symmetry and with Im(λ 5 ) = 0, but here we include λ 6 , λ 7 and Im(λ 5 ). We use the following notation: where γ is Euler-Mascheroni constant. Results below are generalised versions of relations earlier presented in [37,38]

Integration over the superheavy scale
Matching relations for the first step of DR are listed in this section. 14) ) (3.39) The SM result for the three-dimensional mass parameter reads This result can also be found from Ref. [28], apart from the two-loop contributions involving g ′ , as it was assumed to scale as g ′ ∼ g 3/2 . In the 2HDM, full results for the scalar mass parameters read: Re(λ 6 λ 7 ) and µ 2 11,3 =µ 2 11 (Λ) + Re(λ 6 λ 7 ) and finally µ 2 12,3 =µ 2 12 (Λ) + T 2 4 λ * 7 (Λ) + λ 6 (Λ) The renormalisation scale in the three-dimensional theory, Λ 3 , as well as other threedimensional parameters, appear above as exact solutions of the RG equations for the mass parameters of the effective theory. The three-dimensional scale can fixed as Λ 3 = g 2 3 [28].

Integration over the heavy scale
In this section, we list the results for parameter matching in the simplified three-dimensional theories where the heavy temporal scalars and the heavy second doublet have been integrated out.
Integrating out the heavy second doublet In order to integrate out the heavy doublet, we diagonalize the scalar Lagrangian by means of a unitary transformation. Relations between couplings of the off-diagonal and diagonalised theories are given in Appendix A.
In order to make a connection to the existing lattice results of Ref. [7], we fix the renormalisation scale of the final effective theory theory as Λ ′′ 3 =ĝ 2 3 .

Effects from 6-dim. operators
In our O(g 4 ) DR procedure, we omitted the effects coming from higher-order operators.
For the validity of DR, it is important to estimate the effect of these neglected operators; indeed, we observe regions in the parameter space where the lattice parameter x becomes negative, signaling the breakdown of DR (see the Fig. 2 in the companion paper [47]). By using the effective potential and the background field method (see Appendix C.1), it is trivial to obtain the three-dimensional coefficients for the most simple 6-dim. operators, of which we analyse (φ † 1 φ 1 ) 3 and (φ † 2 φ 2 ) 3 . The magnitude of these correlators provide a rough estimate of the validity of DR, yet we emphasise that a full evaluation of the other 6-dim. operators for a more comprehensive estimate is a formidable task.
Coefficients of the aforementioned operators in the effective potential read In the SM case in Ref. [28], in Eq. (201) it was shown that the dominant 6-dim. contribution comes from the top quark (compared to the Higgs self-coupling and gauge contributions) and that the relative shift caused by this 6-dim. operator to the Higgs vev value in the effective theory is of the order of one percent. Therefore, we get a rough estimate of the validity of DR by investigating the ratios If these ratios are large, we can expect that even the first step of the dimensional reduction fails.

Discussion
In this work, we have derived three-dimensional high temperature effective theories for the 2HDM. We emphasize that each of these theories is able to reproduce long-distance physics of the full 2HDM and can be studied on the lattice. An advantage of this three-dimensional approach is its accurate treatment of the IR physics at high temperatures, which cannot be accurately described by purely perturbative methods. For this reason, results obtained from the three-dimensional effective theories can be expected to be more realible than those obtained with perturbation theory alone. However, as the dimensional reduction procedure is perturbative, and first order transitions require large couplings, the validity of the perturbation expansion must be carefully analysed in order to reach reliable results. In the context of dimensional reduction, this requires careful analysis of higher-order corrections to the parameters of the effective theory, and estimates for neglected higher-order operators. In the context of a finite-T effective potential, perturbativity of the one-loop contribution must be analysed by including full 2-loop corrections and investigating relative convergence. These highly-important analyses are left as future prospects.
Using existing lattice results, we have performed a non-perturbative analysis on the three-dimensional theory where the second doublet has been integrated out. In the alignment limit of the 2HDM, our results are presented in a companion paper [47]. Beyond this, a more thorough analysis is still underway, including wider scans of the parameter space in order to identify the regions of first order transitions. Moreover, in the near future we will present non-perturbative determinations of thermodynamic quantities of interest, such as the latent heat and surface tension, in addition to the character and strength of the transition that were analysed in the companion paper [47].
For an accurate statement about thermodynamic properties of the EWPT in BSM theories-and consequences for gravitational-wave signatures and EWBG-it is crucial to perform a non-perturbative analysis. While our present work has tackled this challenge, there remains much to be done, as only the character, strength and the critical temperature of the transition have been determined so far in our work. Non-perturbative evaluation of other thermodynamic quantities of interest, such as the nucleation temperature, bubble nucleation rates and sphaleron transition rates, are left for future work. Non-perturbative evaluation of these quantities is crucial for reliable, quantitative prediction of gravitationalwave signatures and for analysing the viability of baryogenesis.

B.1 Stability and unitarity constraints
For the potential defined by 2.4, boundedness from below is achieved if These have been obtained in Refs. [52,58,[61][62][63]. Furthermore, to guarantee that the v = 246 GeV vacuum is a global minimum of the potential, we impose the following constraint (in the case of a softly-broken Z 2 symmetry): Refs. [28,42]. For pure scalar correlators at one-loop level, we apply both of these methods as a cross-check of the correctness of our calculation. At two-loop, we only need scalar self energies, and we have calculated those diagrammatically, rather than evaluating the two-loop effective potential as was done in [28] in the case of the SM. For this calculation, one has to evaluate a decent amount of individual diagrams; we choose not to present intermediate results in a diagram-by-diagram form, and only the final results have been given in Section 3.2. However, in Appendix C.5 we provide a list of the required integrals for this calculation.

C.1 Effective potential
An economic way of calculating the scalar correlators needed for dimensional reduction is to use the effective potential. In order to calculate the effective potential, one decomposes the scalar fields into quantum and classical fields, φ i → φ i + ϕ i . The quantum fields φ i are integrated over in the path-integral formalism in the usual way, while the classical fields ϕ i are not. The functional form of the effective potential can be found by expanding the resulting potential in the background fields. Following Ref. [28], we include one-loop contributions from scalars, gauge bosons and fermions to the effective potential: Note that for one-loop matching, we only need the non-zero mode contributions of these sumintegrals. We can extract the mass matrix from the quadratic parts in the gauge, scalar, and fermion fields. This matrix can be diagonalized and, using different choices of background fields, we extract the contributions to the correlators at zero external momentum. For simplicity we only consider the Z 2 symmetric case here, so that λ 6 = λ 7 = 0 and furthermore µ 12 = 0. In addition, we set Im(λ 5 ) = 0. The effective potential can be expanded in the background fields, (C.5) By determining the above coefficients V , we can find the correlators needed for dimensional reduction. We shift the scalar fields φ i → φ i + ϕ i , and focus first on the bilinear scalar terms, We encounter a complication when trying to distinguish the V 3 , V 4 and V 5 contributions. In order to separate them, we make three different choices for the background fields Case 2 : For each case, we diagonalise the mass matrix, and evaluate the scalar part of the effective potential by using background-field dependent mass squared eigenvalues, and the sumintegral of Eq. (C.71). By expanding to order O(g 2 ) in mass parameters and to O(g 4 ) in couplings, we find expressions of the form We immediately obtain coefficients V 11 , V 22 , V 1 and V 2 from the above expansions, and the remaining V 3 , V 4 and V 5 can be solved from the linear system of coefficients of v 2 1 v 2 2 , v 2 1 ω 2 0 and v 2 1 ω 2 + .
In the gauge sector, the covariant derivative D µ φ † i D µ φ i , where i = 1, 2, couples the scalar and gauge fields. When the scalar fields are shifted by background fields, we get the bilinear terms The bilinear mixing terms between gauge bosons and Goldstone bosons are removed by fixing a gauge, using the usual Faddeev-Popov gauge-fixing procedure. One must also include ghost fields, with ghost masses and new ghost interactions proportional to the gauge fixing parameter ξ. However, in the Landau gauge, ξ = 0, the ghost masses and interactions vanish. As there is no bilinear mixing between gauge bosons and ghost, we can safely go to the Landau gauge in the ghost sector. The part bilinear in the gauge fields is given by the last two terms above. We follow the same procedure as outlined above, making use of the three choices of background fields. In Case 1, we get two massive charged, one massive neutral and one massless gauge bosons with squared mass eigenvalues , (C.14) In Case 2, we similarly find and in finally Case 3, we find four massive gauge bosons, with eigenvalues By using Eq. (C.3) we can evaluate the coefficients of background fields in each case. Finally, we include contributions from the top quark, that only couples to φ 2 , and thus only affects V 22 and V 2 .
By collecting all contributions from the scalars, gauge bosons, and top quark at oneloop, coefficients of the effective potential of Eq. (C.5) take the form V 11 =µ 2 11 + T 2 12 9 4 g 2 + 3 4 g ′2 + 6λ 1 + 2λ 3 + λ 4 , (C.20a) We have included both tree-level and one-loop contributions. From these coefficients, one can identify the required correlators and the corresponding counterterms.

C.3 One-loop thermal masses
Here we collect the one-loop thermal masses that are needed for thermal counterterms in the four-dimensional theory, In the effective theory containing temporal scalar fields A 0 , B 0 and C 0 , the analogous mass corrections readΠ Contributions from temporal gluons are of higher order, and have been omitted. After the temporal scalars have been integrated out, mass correction for the φ field in the diagonalised theory is given bȳ

C.4 Normalisation of fields
In this appendix, we collect the relations between four-and three-dimensional fields:

C.5 Collection of integrals
The Euclidean four-momentum is denoted as P = (ω n , p) for bosons, where ω n = 2nπT , and as P = (ν n , p) for fermions, where ν n = (2n + 1)πT . In dimensional regularization, spatial integration is performed in d ≡ 3 − 2ǫ dimensions. Shorthand notation for the combined Matsubara sum and space integration: bosons: (sum over nonzero modes), fermions: (C.51) All integrals relevant for O(g 4 ) DR are listed below.

C.5.2 Four-dimensional sum-integrals
Two-loop sum-integrals We introduce the following shorthand notation: m ≡ m 2 + m 2 T , (C.73) S(P, m) ≡ 1 P 2 + m 2 + δ P 0 m 2 T , (C.74) Results for the sum-integrals are given up to terms of order O(|m|T ) ∼ O(gT 2 ). Finally, note that the following fermionic sunset sum-integral vanishes: {P },K 1 P 2 K 2 (P + K) 2 = 0. (C.95) These master sum-integrals are required in evaluation of two-loop scalar two-point functions directly in the unbroken phase, without using the effective potential (which was used in [28]). It will turn out that terms with mixed zero mode and non-zero mode contributions, i.e. terms of the form I 3 2 (m)I 4b 1 , are entirely cancelled in resummation, and so the matching relations are obtained solely from the pure non-zero mode parts.