Radiative effects on false vacuum decay in Higgs-Yukawa theory

We derive fermionic Green's functions in the background of the Euclidean solitons describing false vacuum decay in a prototypal Higgs-Yukawa theory. In combination with appropriate counterterms for the masses, couplings and wave-function normalization, these can be used to calculate radiative corrections to the soliton solutions and transition rates that fully account for the inhomogeneous background provided by the nucleated bubble. We apply this approach to the archetypal example of transitions between the quasi-degenerate vacua of a massive scalar field with a quartic self-interaction. The effect of fermion loops is compared with those from additional scalar fields, and the loop effects accounting for the spacetime inhomogeneity of the tunneling configuration are compared with those where gradients are neglected. We find that scalar loops lead to an enhancement of the decay rate, whereas fermion loops lead to a suppression. These effects get relatively amplified by a perturbatively small factor when gradients are accounted for. In addition, we observe that the radiative corrections to the solitonic field profiles are smoother when the gradients are included. The method presented here for computing fermionic radiative corrections should be applicable beyond the archetypal example of vacuum decay. In particular, we work out methods that are suitable for calculations in the thin-wall limit, as well as others that take account of the full spherical symmetry of the solution. For the latter case, we construct the Green's functions based on spin hyperspherical harmonics, which are eigenfunctions of the appropriate angular momentum operators that commute with the Dirac operator in the solitonic background.


I. INTRODUCTION
The possible metastability of the electroweak vacuum is among the most important features of the Standard Model (SM) and may point to its embedding in the framework of more fundamental theories [1][2][3][4][5]. Due to the contributions from SM fermions to the beta function of the Higgs self-coupling, the latter turns negative at values of the Higgs field much larger than in the electroweak vacuum, generating a lower-lying minimum in the effective potential at high scale. The electroweak vacuum may then decay to this global minimum through quantum tunneling. Prior to the discovery of the Higgs boson, the beta function had been computed to next-to-leading order (NLO) accuracy [6][7][8][9][10][11][12][13][14][15]. The discovery of the 125 GeV Higgs boson [16,17] has motivated further assessments of the metastability of the SM electroweak vacuum at nextto-next-to-leading order (NNLO) [18,19]. The Higgs particular, in false vacuum decay, the so-called bounce provides an inhomogeneous background, corresponding to the profile of a nucleated critical bubble, and the beta functions for the coupling constants do not account for the pertaining gradient effects. One-loop radiative corrections due to fluctuations about the bounce, which fully account for the inhomogeneity of the background, have been calculated using the Gel'fand-Yaglom theorem [37,38]. While the latter is a powerful way of evaluating functional determinants, it has shortcomings: If the quantum-corrected bounce cannot be obtained by improving a classical solution via perturbation theory, further elaboration is necessary. This problem applies to situations where the scale of the nucleated bubble depends on radiative effects, i.e. as occurs in approximately scale invariant theories such as the SM [22,23,39], but also when the true vacuum only emerges radiatively in the first place through, e.g., the Coleman-Weinberg mechanism [40][41][42][43]. Furthermore, evaluating the functional determinant using the Gel'fand-Yaglom theorem does not lead to a systematic method of computing higher-order corrections.
These problems can be addressed by constructing Green's functions in the solitonic background, which can then be used to evaluate effective actions [44,45] and to construct self-consistent equations of motion in perturbation theory or resummed variants thereof [46]. Compared to calculations in homogeneous backgrounds, the reduced symmetry makes it harder, however, to advance calculations to high orders. Nevertheless, problems such as the perturbative improvement of bounce solutions and decay rates [46,47], and finding bounces in radiatively generated [43] or classically scale invariant potentials [39] can be addressed this way. So far, Green's functions and loop effects have been calculated for scalar particles. In view of applications to well-motivated scenarios, such as the SM and other situations where false vacuum decay may play a role [48][49][50][51][52][53], gauge and fermion loops should also be included. The fermionic case is of particular relevance because of the pivotal role played by the top quark in electroweak symmetry breaking. Moreover, it corresponds to a simple example of a field with spin in the background of a nucleated spherical bubble and thus is a suitable setting to develop and test the corresponding methods.
In this article, we develop a method of calculating Green's functions for Dirac fermions coupling to spherically symmetric bounce solutions, as well as to planar configurations. The latter limit allows the analysis to be simplified significantly in the case of quasidegenerate vacua. As an application, we consider the archetypal model originally analyzed by Coleman and Callan [35,36] of tunneling between quasi-degenerate vacua in scalar field theory and compute the leading ra-diative corrections to the bounce, as well as the tunneling action. Referring to the tunneling degree of freedom as the Higgs field, these corrections are induced by the selfinteractions of the Higgs boson, as well as its Yukawa couplings to Dirac fermions and dimensionless couplings to extra scalar fields, such that we can compare these different loop effects.
The paper is organized as follows. In Sec. II, we review the calculation by Coleman and Callan [35,36] of the semi-classical tunneling rate and its first quantum corrections in φ 4 theory, upon which we aim to add fermion loop corrections in this work. In particular, we compute the corrections to the classical bounce solution and the decay rate by means of the Green's functions for the fluctuations about the tree-level bounce. In Sec. III, we present a simple way of computing the Green's functions and functional determinants of the Dirac field in the planar-wall approximation that is applicable to the thin-wall limit. A method of deriving fermionic Green's functions for problems where the vacua are not quasidegenerate, such that the thin-wall approximation does not apply, is presented in App. A. In particular, working in four-dimensional Euclidean space, we construct these Green's functions based on the appropriate set of angular momentum eigenstates of spinors on a three-dimensional hypersphere, such that the problem of finding the radial Green's function separates from the angular problem, as is familiar from calculations in three space dimensions with two-dimensional spherical symmetry. The planarwall limit from Sec. III can be obtained as a limiting case of this construction accounting for the full hyperspherical symmetry. Evaluating the loop corrections based on these Green's functions leads to ultraviolet divergences. These need to be canceled by counterterms that must be chosen to satisfy renormalization conditions on the Higgs mass, self-coupling and wave-function normalization, as we work out in Sec. IV. Based on these developments, in Sec. V, we present numerical results, comparing the contributions from the Dirac and scalar fields to the effective action and the equation of motion for the bounce. Our conclusions are stated in Sec. VI.

II. RADIATIVE CORRECTIONS TO FALSE VACUUM DECAY
In this section, we review the pertinent details of the calculation of the decay rate of a metastable vacuum state. The fundamentals for approaching this problem have been developed in the seminal works by Langer, Coleman and Callan [32,33,[35][36][37], up to the inclusion of one-loop corrections to the decay rate in perturbation theory. The calculation of scalar-field loops in the case of vacuum decay in field theory, and based on Green's functions evaluated in the inhomogeneous background of FIG. 1. The classical potential U (Φ) for for the archetypal scalar theory with false vacuum decay, described in Eqs. (1) and (2).
the tunneling soliton, has been introduced in Ref. [46] and carried out to two-loop order in the decay rate in Ref. [47]. In the present work, we extend this methodology to consider the radiative corrections from fermion loops, making direct comparison with the corresponding effects from scalars. We focus, in particular, on the amplitude for the vacuum transition, as well as the deformations of the classical soliton, and present the main elements of the approach. The formulation in terms of the two-particleirreducible effective action is presented and applied in Refs. [42,43]. That approach is useful when the bounce solutions or even the vacua are dominated by radiative effects in the first place, such as in Coleman-Weinberg [40] scenarios of symmetry breaking or in classically scaleinvariant models [39,[41][42][43]. The main technical developments in that case apply to the negative and zero modes of the symmetry-breaking scalar field. Since fermion fluctuations are the main objective, we develop their computation here in a standard perturbation expansion. Based on this, the generalization to effectiveaction techniques can be implemented straightforwardly, and this may be presented elsewhere.

A. Prototypal Higgs-Yukawa model
We consider a Higgs-Yukawa model based on the Euclidean Lagrangian where and κ is the dimensionless Higgs-Yukawa coupling. The Euclidean gamma matrices are obtained from their Minkowskian counterparts through the replacement γ i → iγ k , for k = 1, 2, 3, and γ 0 → γ 4 . For definiteness, we have chosen the same potential as in the archetypal example of tunneling in field theory considered by Coleman and Callan [35,36], cf. Fig. 1. We choose µ 2 , λ, g > 0, such that there are false and true vacua located at ϕ ± ≈ ± v + O(g/ √ λ), where v = 6µ 2 /λ. The constant U 0 = (µv) 2 /2 − gv 3 /3! is chosen such that the false vacuum has vanishing energy density.

B. Leading-order bounce and tunneling rate
The dominant contribution to the tunneling rate is due to the lowest-lying bounce solution ϕ (and the fluctuations about it). The (tree-level) bounce is a solution to the classical equation of motion that satisfies the boundary conditions ϕ| x4→±∞ = ϕ + andφ| x4=0 = 0, where the dot denotes the derivative with respect to x 4 . The denotes the derivative of the classical potential in Eq. (2) with respect to the field ϕ.
Notice that we are interested in a field configuration that starts and ends (in Euclidean time) in the false vacuum -hence the name bounce. For the bounce action to be finite, we also require that ϕ| |x|→∞ = ϕ + . Given the anticipated O(4) invariance of the bounce, it is convenient to work in four-dimensional hyperspherical coordinates, wherein the equation for the bounce takes the form with r 2 = x 2 + x 2 4 . The boundary conditions become ϕ| r→∞ = ϕ + . The solution must be regular at the origin, and we therefore require that dϕ/dr| r = 0 = 0. The solution is a soliton that interpolates between the 'escapepoint' field value ϕ * (which lies close to the true vacuum ϕ − ) at the origin of Euclidean space and the false vacuum ϕ + at infinity, see Fig. 1. Hence, the bounce describes a four-dimensional hyperspherical bubble nucleated within the false vacuum.
The Euclidean amplitude for transitions from the false vacuum (at τ = − T /2) to the false vacuum (at τ = + T /2) is given in terms of the partition function (cf. Eq. (11)) where H is the full Hamiltonian and S[Φ,Ψ, Ψ] is the Euclidean action. The partition function is to be evaluated by expanding around the bounce solution, and N is a normalization factor given by the inverse of the partition function evaluated in the false vacuum. In terms of Z[0, 0, 0], the decay rate is given by [36] In the thin-wall approximation [36] (see also Ref. [54]), applicable when the minima are quasi-degenerate, i.e. when the cubic coupling g is very small, we may neglect the damping term in Eq. (4), as well as the contribution from the cubic self-interaction gϕ 3 . In this approximation, Eq. (4) has the well-known kink solution [55] where γ = µ/ √ 2. Since we will later use ϕ generally to denote the background scalar field variable, we add a superscript "(0)" to have it uniquely denote the classical bounce ϕ (0) . The radius of the critical bubble R is obtained by extremizing the bounce action with respect to R, which gives and Note that this result for the bounce action contains contributions from both the surface tension of the bubble and the latent heat released by the transition to the true vacuum.

C. One-loop corrections to the action
For non-vanishing external sources, the Euclidean partition function is given by It can be used to compute the one-point functions We proceed by expanding Φ(x) = ϕ(x) + 1/2Φ (x), Ψ(x) = 0 + 1/2Ψ (x) and Ψ(x) = 0 + 1/2Ψ (x), 1 such that the action S[Φ,Ψ, Ψ] is given to quadratic order by where we use the shorthand notation S[ϕ] ≡ S[ϕ,ψ = 0, ψ = 0]. This expansion defines the tree-level inverse Green's functions with δ 4 (x − y) being the four-dimensional Dirac delta function. The Klein-Gordon and Dirac operators in the background of ϕ can then be obtained as Substituting Eq. (13) into Eq. (5), we encounter functional determinants of the operators (15) at quadratic order in the fluctuations. The spectrum of the Klein-Gordon operator G −1 (ϕ (0) ; x) requires special treatment [36]. It contains four eigenmodes φ µ = B −1/2 ∂ µ ϕ (0) with zero eigenvalues. These are the Goldstone modes resulting from the spontaneous breakdown of spacetime translational invariance. The functional integral over the four zero eigenmodes is traded for an integral over the collective coordinates of the bounce, giving a factor [56] where V is a three-space volume. In addition, there is one negative eigenvalue because the bounce actually corresponds to a maximum of the Euclidean action with respect to dilatations of the critical bubble. Along the direction of the negative mode, the path integral can be evaluated by analytic continuation [36,37], leading to an overall factor of −i|λ 0 | −1/2 /2 from the square root of the determinant. Putting everything together, we can therefore express the factor arising from the Gaußian integrations as where the superscript (5) indicates that the five lowest eigenvalues discussed above are omitted.
We still need to normalize by the factor N , that is we must divide by the square root of the determinant of the Klein-Gordon operator evaluated at the homogeneous false vacuum ϕ + ≈ v [36]. For this purpose, we explicitly extract the five lowest eigenvalues 4γ 2 from the determinant of G −1 (v; x) that are not canceled when building the quotient of the determinants [57]; namely, The logarithms of these determinants appear as additive corrections to the classical action B[ϕ (0) ] and, within the present approximations, they can be interpreted as the one-loop contribution to the effective action. Together with the corresponding fermion terms, we denote these as where and throughout this paper the determinant of the Dirac operators is understood to be taken over both the coordinate and spinor spaces. Here, and in the remainder of this work, we suppress the argument in B S,D when evaluated at the classical bounce, in accordance with what is implied for B. We note that these expressions still need to be renormalized, which we will carry out in Sec. IV.

D. One-loop corrections to the bounce
Before we assemble the preceding results into an expression for the decay rate, we consider, in addition, the first quantum corrections to the bounce. To this end, we write ϕ = ϕ (0) + δϕ and expand Eq. (12a) to next-toleading order. We then obtain where we define the tadpole functions from scalar and fermion loops as in which tr s indicates the trace over the (suppressed) spinor indices. Alternatively, rather than the expansion of Eq. (12a) based on the partition function, the one-loop tadpole functions also descend from the functional determinants by functional differentiation; specifically, We note that fermion loops at coincident points do not appear in diagrammatic expansions of Green's functions in Yukawa theory but they do appear in the expansion of the one-point function, cf. Eq. (22b). To clarify how the coincident fermion propagator emerges in the tadpole function, we carry out the variation of the one fermionloop effective action in App. B explicitly. At one-loop accuracy, the bounce satisfies 2 i.e. the tadpoles can be interpreted as corrections to the equation of motion.
When substituting the quantum-corrected bounce into the action S[ϕ] and the one-loop terms B where and where, for the first identity, we have used the fact that the classical bounce ϕ (0) is the stationary point of the classical action and, for the second identity, we have used Eq. (21). Expanding we obtain from Eq. (20) Comparing with Eq. (26) and using Eq. (23), we see that the total correction to the expansion in δϕ is [46] B (2) ≡ δS + δB Diagrammatically, the contributions to B (2) correspond to (one-particle reducible) dumbbell graphs, cf. Fig. 3 (b) for the fermion contributions. They therefore constitute a subset of the two-loop corrections. However, this subset can be the dominant two-loop contribution for theories with a large number of degrees of freedom propagating in the loop [46]. For the numerical examples in Sec. V, we consider such a setup with a large number of fermion and scalar fields coupling to the Higgs degree of freedom Φ.

E. Radiatively-corrected decay rate
We can now summarize all the contributions that we include in the approximation of the tunneling rate per unit volume as where B, B (1) S and B (1) D are evaluated at the classical bounce. Note that the exponent in this expression corresponds to an approximation of the effective action, up to the modes that are not positive definite but lead to a vanishing contribution in the planar-wall limit, cf. Eq. (35) and the associated comments.
In summary, given the leading-order approximation ϕ (0) to the bounce, Eq. (7), and B for the action, Eq. (8), we apply the following procedure in order to calculate the radiative corrections to the bounce and to the decay rate: • We first invert Eqs. (14a) and (14b) to find the Green's functions G(ϕ (0) ; x, y) and D(ϕ (0) ; x, y).
• The Green's functions are used in order to calculate B (1) S and B (1) D according to Eq. (20) and the discussion in Sec. IV C.
• The Green's functions are also used to obtain the tadpole functions (22), which, in turn, yield the corrections to the bounce (21).
• When substituted into the tree and one-loop action, the radiative corrections to the bounce yield the quadratic correction B (2) in Eq. (30) by means of Eq. (26), corresponding to dumbbell graphs. We account for these contributions in the calculation of the decay rate because they can be relevant when a large number of fields is running in the loops.
• Finally, the pieces B, B S,D and B (2) can be put together to obtain the decay rate per unit volume (31).

III. GREEN'S FUNCTIONS, LOOP-IMPROVED EFFECTIVE ACTION AND BOUNCE IN THE PLANAR-WALL APPROXIMATION
In this section, we calculate the Green's functions for the scalar and fermion fields in the background of the classical bounce and in the planar-wall limit. A method of calculating the Green's function of the Dirac operator in a spherically symmetric background in terms of all of its spinor components is worked out in App. A. In the approximation of a planar wall, it is calculationally simpler, however, to obtain the functional determinant, as well as the fermion tadpole Π D , from the Green's functions of the corresponding Klein-Gordon-like operators, obtained from "squaring" the Dirac operator, as we discuss below.

A. Planar-wall limit of Green's function and functional determinants
First, we review the calculation of the Green's function for the Klein-Gordon operator, as presented in Ref. [46]. Equation (14a) for the scalar Green's function can be recast to the standard form In a spherically symmetric background, this equation can be solved by separating the angular part by means of a partial-wave decomposition. This reduces the problem to one of finding the hyperradial function G j (ϕ (0) ; r, r ), where j is the quantum number of angular momentum, see, e.g., Ref. [46]. However, to keep things simple, we apply here the planar-wall approximation: When the radius of the bubble wall (9) is very large compared to µ −1 , i.e. when g √ λµ, we can treat the geometry as planar. The bounce is then a function of the coordinate z ⊥ perpendicular to the wall only and has no dependence on the parallel coordinates z on the three-dimensional hypersurface, cf. Fig. 2. In the remainder of this article (save App. A), we employ this planar-wall approximation.
Without loss of generality, we choose z ⊥ ≡ z = x 4 − R and z = x. The Fourier transform of the Green's function with respect to the coordinates x, satisfies the equation and we refer to it as the three-plus-one representation. It can be found numerically or analytically, and the analytic solution is derived in Refs. [39] and [46]. Starting from the representation of the Green's function as a spectral sum, one can show [58][59][60][61] that it can be used to calculate the functional determinant directly, and this approach was used in the planar-wall limit for scalar fluctuations in Ref. [43]. Specifically, the func-tional determinant can be evaluated as where, making use of the isotropy parallel to the bubble wall, we have written the coincident Green's function G(ϕ (0) ; z, z, |k|) ≡ G(ϕ (0) ; z, z, k). The d 3 x integration is parallel to the bubble wall, and Λ is a three-momentum cutoff. Note that, in these integrals, the negative and zero modes correspond to |k| = 0 and therefore lead to a vanishing contribution.
In analogy with the scalar case, the Green's function for the Dirac operator in a spherically symmetric background can be solved by a separation ansatz. The angular solutions in this approach are spin hyperspherical harmonics of definite total angular momentum, as is explained in detail in App. A. Considerable simplifications are, however, possible in the planar-wall approximation.
In the planar-wall approximation, ϕ (0) depends on z (or x 4 ) only, and we must define the local Dirac mass m D (z) ≡ κ ϕ (0) (z). We therefore need to generalize the well-known method of evaluating the determinant of the Dirac operator to the situation where the mass may vary in one direction of spacetime. Noticing that the determinant is invariant under chiral conjugation, we multiply by γ 5 from the left and right of the inhomogeneous Dirac operator. This yields where we have used the anticommutation relations between γ 5 and γ µ and the multiplicativity property of the determinant. We can thus express Now, after employing the representation of gamma matrices where σ i are the Pauli matrices, we obtain This is a determinant in block form, and we can make use of the relation which applies even when A and B do not commute with each other, provided A and B are square matrices of the same dimensions. We therefore arrive at where, in the second step, the determinant over the remaining 2 × 2 block has been performed. When m D is constant, this reduces to the well-known result for the fermion determinant. The result of this computation can be substituted into Eq. (31) for the decay rate.
With the above results, we can now rewrite the fermion radiative correction (20b) to the tunneling action as where B (1) The problem has thus been reduced to one of evaluating the determinants of scalar fluctuation operators. This can be accomplished by first calculating the Green's functions D ± (ϕ (0) ; x), as well as their Fourier transforms with respect to the coordinates parallel to the bubble wall, as in Eq. (33). From these, the quantities B (1) D± can be obtained in analogy with the determinant for the scalar field in Eq. (35), i.e.

B
(1) The fermionic tadpole is now most straightforwardly obtained by functionally differentiating the one-loop determinant as in Eq. (23): We emphasize that δ/δϕ denotes the functional derivative and, after making use of the Jacobi formula, we obtain from Eq. (44) that These expressions can be substituted into Eqs. (21) and (22b) in order to obtain the fermion contribution to the one-loop correction δϕ to the bounce. We also note from Eqs. (22b), (46) and (47) that there follows the following identity for the spinor trace: In the homogeneous limit ϕ (0) = const., D + and D − coincide, and we quickly recover the familiar result for the one-loop fermion determinant -including the factor of 4 from the spinor trace -after summing over ±.
We conclude this part of the discussion by noting that it is also possible to compute D(ϕ (0) ; x, x ) and all of its spinor components directly by solving the equations derived in App. A. While this is calculationally more cumbersome, we have checked the analytical and numerical results for the fermion loop effects reported on in this work using both methods. B. One-loop correction to the bounce in the planar-wall limit The one-loop corrections to the bounce in Eq. (21) satisfy the equation of motion as we can quickly confirm by acting on Eq. (21) with the tree-level Klein-Gordon operator. This can be solved using the Green's function G(ϕ (0) ; z, z , k) at k = 0, for which there is an analytic solution [46]: where u = tanh(γz) (51) and ν = 2 1 + k 2 /(4γ 2 ). Proceeding in this way, we obtain As can be seen from Eq. (50), G(u, u , ν) is singular as ν → 2, i.e. |k| → 0. However, since Π(u ) is an even function, the part of the integrand multiplied by G(u, u , 2) is an odd function in u . The integral remains finite because the singularity of G(u, u , 2) turns out to reside in the even part. We can therefore replace G(u, u , 2) in Eq. (52) with the odd part of the Green's function [46] which can be expressed as [46] G odd (u, u ) = ϑ(u − u ) 1 32 γ in the domain 0 ≤ u, u ≤ 1. For given Π(u), we can thus compute the first correction to the bounce.

IV. RENORMALIZATION
The calculation of radiative corrections to the decay rate and to the bounce solutions requires the renormalization of the scalar operators. We therefore add the following counterterms to the Lagrangian in Eq. (1): The one-loop corrections entail, in particular, the usual ultraviolet divergences from the quartic scalar and Yukawa interactions. In position space, these appear when taking the coincident limit of the Green's functions and, in momentum space, when performing the loop integrals. In the present setup, we use the mixed threeplus-one representation of momentum space parallel to the three-dimensional hypersurface and position space in the perpendicular direction. As a regulator, we introduce a three-momentum cutoff Λ. The counterterms must be uniquely specified by certain renormalization conditions. For the problem of metastable vacua, we find it most useful to impose these conditions on the derivatives of the effective action evaluated at the false vacuum. In the one-loop approximation, the effective action is given by (56) When ignoring derivative operators, the effective action in the false vacuum coincides with the Coleman-Weinberg effective potential. We therefore use the derivatives of the latter to define the renormalization conditions for µ 2 and λ, and these are worked out in Sec. IV A.
The radiative effects also lead to corrections to the wave-function normalization. These are logarithmically divergent for the fermion loops and finite for the scalar loops. In Sec. IV B, we extract these terms analytically by a gradient expansion of the Green's functions. From the latter, we can calculate the leading gradient corrections to the effective action, which can then be identified with the wave-function normalization, where we impose as a renormalization condition the standard unit residue at the single-particle pole. In Sec. IV C, we then summarize these results by combining the counterterms with the one-loop determinants and tadpole functions. In this way, renormalized results for B (1) , B (2) and the tadpole function are obtained, which lead directly to the renormalized decay rate and the correction δϕ to the bounce.
A. Renormalization of coupling constants using the Coleman-Weinberg potential The Coleman-Weinberg effective potential is obtained by evaluating the effective action (56) for configurations ϕ that are constant throughout spacetime. The latter is given by The spacetime and four-momentum variables x and k are understood to have natural units; the overall factor of in the one-loop corrections appears only for bookkeep-ing purposes. 3 It can be directly read from the above expression that We note that we could just as well compute these quantities from the Green's functions in the three-plus-one representation for the limit of a constant field background, using Eqs. (35) and (45) (cf. Eqs. (77) and (81) later). Factoring out the integral over the spacetime fourvolume, one gets the Coleman-Weinberg effective potential where we reiterate that the factor of again appears only for bookkeeping purposes.
In order to apply a regularization scheme that is congruent with our calculations of the decay rate and the corrections of the bounce (both of which are performed in the background that is spatially inhomogeneous in the direction perpendicular to the bubble wall), we make a three-plus-one decomposition of momentum space. Performing the integral over k 0 , we obtain Evaluating the remaining integral up to a cutoff Λ, this yields (dropping the bookkeeping factors) 3 Were we to proceed otherwise, we would need to include additional inverse powers of in the logarithms to ensure the mass terms have the correct dimensions. Note that these factors would then cancel the overall factor of when one arrives at, e.g., Eq. (61).
The renormalized Coleman-Weinberg potential is obtained when adding the counterterms δµ 2 and δλ in accordance with Eq. (55) as where δµ 2 and δλ are specified by the following renormalization conditions: These counterterms are given explicitly in Eqs. (84) below, along with the remaining one for the wave-function correction.

B. Wave-function renormalization through adiabatic expansion of the Green's functions
While the effective potential offers a convenient way to define the counterterms for the coupling constants, it does not lead to conditions on the renormalization of the derivative operator, i.e. the wave-function normalization. Our objective is to express this additional counterterm in an analytic form. This can be achieved by performing a gradient expansion of the Green's functions around a constant background field configuration. One may also refer to this calculation as an adiabatic or Wentzel-Kramers-Brillouin (WKB) expansion.
We first construct the adiabatic expansion of the scalar Green's function, which satisfies with This can be solved by the ansatz where and where we impose The latter enforce the boundary condition that the Green's function vanishes at infinity.
In order to isolate the leading gradient effects, we make the WKB ansatz When substituted into Eq. (67), this leads to where a prime on W (z) (and in the following on m S (z)) denotes a derivative with respect to z. In the gradient expansion, the zeroth-order approximation to W is given by W (0) 2 = M 2 . Substituting this back into Eq. (70), we obtain Demanding continuity and the correct jump in the first derivative to reproduce the delta function in Eq. (64) at the coincident point z = z , we obtain the matching conditions such that where is the Wronskian. Putting these results together, we find that the Green's function is approximated to second order in gradients by In terms of the local squared mass m 2 S (z) = U (ϕ (0) (z)), this reads The first term on the right-hand side can be identified with the Green's function in the three-plus-one representation and a homogeneous background, i.e. for field values that are constant throughout spacetime. We explicitly define this as G hom (ϕ (0) ; z, z , k) = 1 2 where we show the general expression that also holds away from the coincident points in order to highlight the exponentially decaying behavior of the Green's functions for large separations. When discussing parametric examples in Sec. V, we will use G hom in order to compare loop contributions without gradient effects with those that include them. The last two terms in Eq. (76) are the leading gradient corrections. These corrections are ∼ 1/|k| 5 such that the three-momentum integral leads to a finite correction to the derivative operators of the field Φ, i.e. the wave-function normalization. It is most straightforward to determine the correction to this operator from its contribution to the effective action, which we can compute by substituting the result (76) into Eq. (35): where we recall that B S is related to the functional determinants through the definitions (20). While this result is finite, i.e. independent of the momentum cutoff, the corresponding contribution from the Dirac field is not, as we will see below.
We compute the fermionic contributions to the wavefunction normalization analogously.
The propagators D ± (ϕ (0) ; z, z , k) satisfy the same equation as G(ϕ (0) ; z, z , k), i.e. Eq. (64), but with Solving that equation in the WKB approximation, we find Again, in order to isolate the gradient effects in Sec. V, we define the contributions that arise in homogeneous backgrounds without gradients as We note that, while the three-momentum integrals over the derivative terms in Eq. (80) are finite, the trace of the coincident Dirac propagator (48), and therefore the tadpole correction, cf. Eqs. (22b), (46) and (47), contains an extra derivative with respect to z. This generates a logarithmically divergent contribution that can be identified with the divergent part of the wave-function normalization. Computing the corrections to the derivative operators in the effective action using Eq. (45), we find We use this result along with Eq. (78) in order to specify the renormalization condition where the Lagrangian is implied to contain the counterterms as per the definition (55). This relation fixes the counterterm δZ.

C. Renormalized bounce, effective action and decay rate
We can now summarize the one-loop counterterms as The renormalized tadpole correction can then be defined as which should replace Πϕ in the equation of motion (24) for the bounce. Finally, to the one-loop contributions to and the renormalized result for B (2) is obtained when replacing Π with Π ren in Eq. (30).

V. QUANTUM CORRECTIONS TO THE BOUNCE AND DECAY RATE FROM FERMION AND SCALAR LOOPS
We now apply the methods elaborated in the previous sections in order to compute the leading radiative corrections to the bounce and to the decay rate numerically. In presenting the results, we isolate the contributions from the scalar and the fermion loops.
The term B (2) D turns out to contain the dominant twoloop contributions from fermions to the effective action when enhanced by multiple copies of the fermion fields, i.e. Ψ i with i = 1, . . . , N Ψ . To illustrate this, we present a diagrammatic representation of the fermionic corrections to the bounce action in Fig. 3. Diagram (a) is the one-loop term B (1) D of order (κ 4 /λ 2 )N Ψ . Note that the dependence on λ comes from the background bounce ϕ (0) ∼ 1/ √ λ appearing in the Dirac operator. Diagram (b) is the main contribution to B (2) D , which is of order (κ 6 /λ 2 )N 2 Ψ , cf. Eqs. (30) and (52). In addition to diagram (b), there is also a contribution of order (κ 6 /λ 2 )N Ψ represented by diagram (c), which we do not calculate. Diagram (c) is therefore suppressed by a relative factor of 1/N Ψ compared to diagram (b) and may hence be neglected, cf. also Ref. [46].
In contrast to the Higgs field Φ, the vacuum expectation values of the fermion fields do not change through the bubble wall, such that their role may be referred to as spectator. In order to directly compare effects from fermion and scalar loops, we introduce additional scalar species χ i with i = 1, . . . , N χ that also couple as spectators to the Higgs field. In this setup, the Lagrangian without counterterms is For simplicity and in order to have the scalar and fermionic spectators have similar properties, we take m χi = 0. The developments for the effective action and the tadpole corrections of the previous sections generalize straightforwardly to loops of the fields χ i , yielding oneloop corrections B (1) χ to the action and tadpole functions Π χ .
Since it is of interest to compare the loop contributions from the various species, we decompose the tadpole function and the one-loop corrections to the action as where the subscript S indicates the Higgs scalar field Φ, D indicates the Dirac fermion spectators Ψ i and "spec", the scalar spectator fields χ i . We note that all results presented in this section for the tadpole corrections Π, as well as the corrections to the effective action B (1,2) , are understood to be renormalized, i.e. to include contributions from counterterms. In order to identify the impact of the gradient effects, we can use the propagators G hom (ϕ (0) ; x, x) and D hom (ϕ (0) ; x, x) in place of the Green's functions G(ϕ (0) ; x, x) and D(ϕ (0) ; x, x) in Eqs. (35) and (45) to compute the one-loop action terms in the homogeneous background. We can compute the corresponding spectator corrections, leading altogether to the results Accordingly, substituting the homogeneous Green's functions in Eqs. (22) (and analogously for the scalar spectators), we also obtain

A. Tadpoles and corrections to the bounce
The tadpole correction is Π ϕ (0) , where the classical bounce ϕ (0) is an odd function and Π is an even function about the center of the bubble wall. In Fig. 4, we show the total Π along with the contributions from the different species that run in the loops.
We compare these full one-loop results with the tadpole functions without gradient corrections. For each value of u, we calculate these assuming a constant background ϕ = ϕ (0) (u). For the field Φ, the squared mass may become negative for a certain range of u such that the resulting tadpole function acquires an imaginary part that we do not show in the diagrams. The imaginary part arises if we assume a constant field configuration, because there is a continuum of negative modes contributing to the Gaußian integrals in the tachyonic region where −µ 2 + (λ/2)ϕ 2 < 0. It is, however, an artefact since constant configurations do not correspond to extremal points (or saddle points) of the effective action away from the minima of the potential at ± v. In turn, as explained in Ref. [46], imaginary parts, save the one associated with the negative mode, do not appear when calculating the loop diagrams for fluctuations around the full bounce solutions that also account for gradients, i.e. when expressing the loops in terms of the Green's function solutions in the bounce background.
The parameters chosen for Fig. 4 are a set which we also use as a benchmark point for the remaining numerical results presented in this section. The bumps around u = ± 0.4 in the graph of Π S without gradient contributions are due to the transition to the tachyonic region in the classical potential. This causes a divergence in the derivative of Π S . Note that in the graph of the full one-loop result for Π S , which accounts for gradient corrections, these bumps are absent because the field gradients counteract the tachyonic instability. For the fermionic and scalar spectator fields, there are no tachyonic regions across the bubble wall. For the fermionic and scalar spectator fields, we can see from Fig. 4 that the gradient effects suppress the tadpole corrections when compared with the corrections without gradients. For the Higgs field, the tadpole correction appears to be enhanced compared to the real part of the loop correction without gradients. One should note, however, that this is not directly comparable with corrections from the spectator fields because of the tachyonic modes and the imaginary part in the Coleman-Weinberg potential that has no direct physical interpretation.
The parameters in Eq. (92) have been chosen such that there is a substantial amount of accidental cancellation between the fermion and scalar loop contributions. In Fig. 4, we see this in the plot of the total tadpole correction Π, from which we also conclude that, in such a situation, the relative impact of the gradient corrections can be of order one. For all tadpole corrections Π, we The tadpole Π is the total result, ΠS is the contribution from scalar fields Φ, Πspec is from scalar spectator fields χi and ΠD is from Dirac spectator fields Ψi. Graphs with solid lines show the tadpole corrections Π including gradient effects; dotted graphs show the corrections Π hom where these effects are ignored. The correction δϕ to the bounce as a function of γ(r− R) with (δϕ, solid) and without (δϕ hom , dotted) gradients in separate panels for the contributions δϕS from the Higgs field, δϕspec from the scalar spectators, δϕD from the Dirac spectators and for the total correction δϕ. note that these are largest around the point u = 0, corresponding to ϕ (0) = 0. From Eq. (21), we also observe that it is the combination ϕ (0) Π(ϕ (0) ; y) that acts as a source for the radiative correction. We can therefore expect that the relative impact on the radiative correction δϕ to the bounce is suppressed. In Figs. 5 and 6, we compare the radiative corrections to the bounce with and without gradient effects. These corrections are calculated using Eq. (47), where we substitute Π S,D,spec and Π hom S,D,spec in order to obtain the particular contributions δϕ S,D,spec and δϕ hom S,D,spec . These are linear responses to the particular one-loop tadpoles, such that δϕ = δϕ S + δϕ D + δϕ spec and accordingly for the corrections with the superscript "hom" that exclude gradient effects. We observe that, for the Dirac spectator fields Ψ i (subscript "D") and the scalar spectators χ i (subscript "spec"), inclusion of the gradient effects smooths the field profile, that is the turning points in the functions δϕ D,spec are softened relative to those in δϕ hom D,spec . For the correction δϕ S from Higgs loops, there appears to be an opposite effect. However, we recall that we have dropped the imaginary parts from Π hom S such that this contribution is not directly comparable with the one from the spectator fields. We also note the larger total correction from the Dirac spectators when compared with their scalar counterparts. The apparent relative factor of minus four can, of course, be accounted for by the number of Dirac spinor degrees of freedom and the opposite sign of the fermion loop. [The squared mass of the Dirac spectator fields across the wall is κ 2 ϕ(z) 2 , for the scalar spectators it is α ϕ(z) 2 /2. For the values of α and κ chosen in Eq. (92), these are therefore coincident, such that the relative factor of minus four is actually exact when ignoring gradients.] The relative impact of the gradient corrections on δϕ turns out to be about twice as large for the Dirac spectators when compared with the scalar spectators.

B. Corrections to the action
We now investigate the gradient effects on the action and the tunneling rate. For this purpose, we define the contributions to the action per unit area of bubble surfacē and accordingly for the individual contributions from the Higgs, Dirac and scalar spectator fields, as well as for the loop correctionsB (1,2)hom derived from the Green's functions in the homogeneous background. For the benchmark parameters (92), the results for the various corrections to the action are compared in Table I. With or without gradient effects, we see that the corrections B (1) for the scalar loops are negative (leading    (1) andB (2) (i) without gradients (i.e. based on the Green's functions in homogeneous backgrounds) and (ii) with gradients (i.e. based on the Green's function in the background of the tree-level soliton). These quantities are computed for the benchmark point (92). We draw attention to the fact that the values ofB (2) differ in sign for cases (i) and (ii), leading to a relative increase in the tunneling rate when the gradients are included. For completeness, the value ofB -which does not differ between cases (i) and (ii) -is 2.828.
to an enhanced decay rate) and for fermion loops, these are positive (leading to a suppressed decay rate). Since our renormalization conditions fix λ at the false vacuum, fermion loops lead to a running toward smaller |λ| and scalar loops toward larger |λ| around the center of the bounce, where ϕ = 0. From the tree-level action (10), we see that larger (smaller) λ lead to a faster (slower) decay rate, which qualitatively explains this numerical observation. As for the relative contributions from field gradients to B (1) shown in Table I, we note that for Dirac spectators, these are larger by a factor of about four in comparison with the contribution from scalar spectators. In Figs. 7, 8 and 9, we show the results of varying the dimensionless couplings λ, κ and α separately about the benchmark point given in Eq. (92). Only the Dirac spectator correctionB (1) D depends on κ and only the scalar spectator correctionB (1) spec depends on α, which is why these are the functions that we present in Figs. 8 and 9. When comparing the Dirac and scalar spectator contributions in Figs. 8 and 9, one should bear in mind that the squared masses depend on the couplings as κ 2 ϕ(z) 2 and α ϕ(z) 2 /2, respectively. This explains the stronger curvature in the plot for Dirac spectators. Since a change in the self-coupling λ of the Higgs field changes the shape of the bounce, this is of relevance for all three types of one-loop corrections,B (1) spec . From Fig. 7, we see, however, thatB (1) S is insensitive to the value of λ. This is because, in our perturbation expansion, one loop corresponds to one power in λ, such thatB (1) S is one order higher in λ compared toB, which itself scales as 1/λ, cf. Eq. (10). In principle, there are further logarithmic dependencies on λ [cf., e.g., the effective potential (61)] but these turn out to cancel when performing the integral over dz, as is shown analytically in Refs. [46,47,57].
In accordance with the above remarks on the effect of fermion and scalar loops on the decay rate, we can understand the dependence ofB D increases and therefore the decay rate decreases. This is due to a higher barrier from fermion fluctuations in the effective potential for larger Yukawa couplings. Similarly, as α increases,B (1) spec decreases and therefore the decay rate increases. This implies that there is a lower barrier in the effective potential for larger couplings α. To illustrate this point, the dependence of the barrier in the Coleman-Weinberg po- tential (61) on κ and α is shown in Fig 10. We emphasize, however, that, while the Coleman-Weinberg effective potential can be used in order to interpret the leading effects from radiative corrections, it does not include the subleading gradient effects on the particles running in the loops. We note that the bounce evolves rapidly about z = 0, where the gradient effects on the tadpoles are largest, such that their relative impact on the one-loop action, which is an integral quantity over z, is comparably small. In Fig. 11, we explicitly isolate the gradient effects by plotting δB (1) =B (1) −B (1)hom for the one-loop contributions from the particular species.
In Fig. 12, we compare the dependence of B (2) on the various coupling constants with and without gradient effects. The relative difference between the cases with and without gradient effects is of order one. When recalling the dependence of B (2) on the tadpole functions given in Eq. (30), we see that this sizable difference is because of the large relative impact of gradient effects on Π due to the cancellation from fermion and scalar loop effects, cf. Fig. 4. We reiterate, however, that this cancellation is coincidental due to the parameter choice in Eq. (92).

VI. CONCLUSIONS
We have extended the Green's function method proposed in Ref. [46] to fermion fields in order to deal with the radiative effects on false vacuum decay in Higgs-Yukawa theory. The bounce background provides an inhomogeneous mass term m D = κϕ (0) for the fermions via their Yukawa coupling. We have therefore constructed the Green's functions for the Dirac operator in a Euclidean spacetime-dependent background with O(4) spherical symmetry. In the thin-wall approximation, however, we have found that it is possible to reduce the functional determinant of the Dirac operator to functional determinants of two scalar operators of second order in derivatives, which simplifies explicit calculations. We have then used this procedure to calculate the radiative corrections to the bubble profile and the decay rate  numerically in the thin-wall approximation.
The results show that, while the quantum corrections from scalar loops broaden the bubble wall and lead to a faster decay rate, fermion loops, on the other hand, have the opposite effect. These conclusions could have been reached also when ignoring the gradient effects. However, essential for this work is that we numerically calculate these quantum corrections, fully capturing the gradient effects beyond the Coleman-Weinberg effective potential. In particular, the profile of the radiatively improved bounce solutions are found to be smoother when the gradient effects are included. In the planar-wall approximation, the corrections to the bounce profile from the gradient effects are small. This is because these corrections are proportional to (Π − Π hom )ϕ (0) . While Π − Π hom peaks at the bubble wall u = 0, the classical bounce ϕ (0) vanishes there, leading to a suppression of the effect on the decay rates [46]. However, these conclusions reached for the archetypal example of tunneling between quasi-degenerate vacua in scalar theory with quartic interactions may change when considering phenomenological models leading to different bubble profiles. One interesting application of this work would therefore be a study of the gradient effects on top-quark loops in electroweak phase transitions and the decay of the electroweak vacuum. These tasks require a numerical calculation and renormalization of the fermionic Green's functions in a spherical background, based on the methods developed in the present work, and further, in the case of a thermal phase transition, the inclusion of finite-temperature effects. The differences in the coupling of the scalars and fermion fields to the inhomogeneous background may also be of interest in the study of solitons in supersymmetric theories.
These form a basis {N i , N i } of the Lie algebra o(4): where ijk is the Levi-Civita tensor. We can introduce another basis {A i , A i }: The Eigenstates of A 2 and A 2 have eigenvalues j(j + 1) and j (j + 1), respectively. Eigenstates of the total orbital angular momentum operator (A7) have eigenvalues p(p+2)+q 2 = 2[j(j +1)+j (j +1)] and transform as {j/2, j /2} representations with partition numbers ( j+j 2 , j−j 2 ). Under spatial reflection: The spatial reflection symmetry of O(4) therefore forces j = j for the total angular momentum eigenstates. Thus, the total orbital angular momentum eigenstates transform as {j/2, j/2} representations of O(4) with partition numbers (j, 0), and 1 2 M · M has eigenvalues 4j(j + 1).
In the case of the spin group Spin(4) SU L (2) ⊗ SU R (2), we have the generators with the six independent operators forming a basis of the Lie algebra so(4): As before, we can introduce another basis {S L i , S R j }: The sets of operators S L ≡ {S L i } and S R ≡ {S R i } form two commuting su(2) algebras, reflecting the fact that the so(4) algebra is equal to the direct sum of two su(2) algebras, i.e. so(4) = su L (2) ⊕ su R (2), and, correspondingly, that the SO(4) group is locally isomorphic to the direct product of two SU (2) groups, i.e. SO(4) SU L (2) ⊗ SU R (2). The operators S L act on the subspace generated by the projector P L = 1 2 (I 4 −γ 5 ), and the operators S R act on the subspace generated by the projector P R = 1 2 (I 4 + γ 5 ). We call spinors that transform in these two subspaces left-and right-handed, respectively. Lefthanded spinors transform as {1/2, 0} representations and are labeled by the partition numbers (1/2, +1/2); righthanded spinors transform as {0, 1/2} representations and are labeled by the partition numbers (1/2, −1/2).
In terms of the above generators, we can write the Dirac operator in the form wherex µ ≡ x µ /|x| ≡ x µ /r. Hence, provided we can find the Green's functionD of the operator in square brackets, the Green's function of the complete Dirac operator will be given by D = γ ·xD. Our aim then is to find the eigenstates of M · Σ, i.e. the simultaneous eigenstates of 1 2 M 2 , 1 2 Σ 2 and the total angular momentum K 2 ≡ Now we consider the eigenstates of M · Σ. We need to consider the coupling between the representation of O(4) and that of Spin(4). We first couple the states transforming under O(3) with those under SU L (2), leading to the quantum numbers J and m J , and then proceed similarly for O (3) and SU R (2), leading to the quantum numbers J and m J . Finally, we couple the two resulting states, replacing the quantum numbers m J , m J by L, M . There is a unitary transformation relating these two couplings [62]: and, correspondingly, the partition numbers take values P ≡ J + J = j + s, |j − s| , (A19a) Q ≡ J − J = ± s, ∓ s .
The states have the following eigenvalues of the total angular momentum: In the coordinate representation, the orbital angular momentum states are given by the four-dimensional hyperspherical harmonics Y j m (e r ) = e r |(j, 0), , m , where the four-dimensional unit vector is e r = (cos χ, sin χ cos θ, sin χ sin θ sin ϕ, sin χ sin θ cos ϕ) .
(A23) In terms of these coordinates, the hyperspherical harmonics in four dimensions are given by [66] Y j m (e r ) = 2 Γ( + 1) 2(j + 1)(j − )! π(j + + 1)! We now turn out attention to the problem of finding the fermionic Green's function for a radially varying mass m(r). The Euclidean-space Dirac equation takes the form We can proceed by making the following ansatz for the solution: D(x, x ) = λ a λ (r, r ) + b λ (r, r ) γ · x D λ (e r , e r ) , Given the completeness of the eigenstates in Eq. (A17), this leads to the equation γ · x x · ∂ − J r 2 a λ (r, r ) + m(r)b λ (r, r ) D λ (e r , e r ) + m(r)a λ (r, r ) + ∂ ∂r + J + 3 r r b λ (r, r ) D λ (e r , e r ) = δ(r − r ) r 3D λ (e r , e r ) .
Here, we have also used Eq. (A13) and the relation The radial equation only depends on the orbital angular momentum quantum number j and the partition number P . The two linearly independent components of Eq. (A30) allow us to eliminate b λ such that we obtain a second-order equation. Substituting for the allowed values for J from Eq. (A21) leads to (A32) Note that, for a constant mass, this coincides with the spherical Klein-Gordon equation, as one would expect. We are ultimately interested in the coincident limit of the fermionic Green's function. It is therefore useful to evaluate the dyadic product of the eigenstates in Eq. (A17) and to trace over all but P and j. This can be simplified dramatically, since we may take χ = 0 without loss of generality. Doing so, = 0 in all non-vanishing contributions, and we find the result Q,L,M,± | e r |(j, 0), (s, ±s); (P, Q), L, M | 2 = 1 4π 2 (j + 1) K I 4 , where K =      j + 2 , P = j + s , j > 0 , j , P = j − s , j > 0 ,