Adiabatic regularization with a Yukawa interaction

TL;DR

This paper extends adiabatic regularization to include Yukawa interactions between fermions and scalar fields in an expanding universe, providing explicit renormalized quantities and confirming consistency with the heat kernel method.

Contribution

It introduces a method to regularize fermionic fields with Yukawa interactions in curved spacetime, including explicit expressions and counterterms.

Findings

Derived explicit expressions for renormalized stress-energy tensor and bilinear in FLRW spacetime.

Identified covariant counterterms for Yukawa interactions in curved spacetime.

Confirmed the approach's consistency with the heat kernel method in the massless limit.

Abstract

We extend the adiabatic regularization method for an expanding universe to include the Yukawa interaction between quantized Dirac fermions and a homogeneous background scalar field. We give explicit expressions for the renormalized expectation values of the stress-energy tensor $\langle T_{\mu\nu} \rangle$ and the bilinear $\langle \bar\psi\psi\rangle$ in a spatially flat FLRW spacetime. These are basic ingredients in the semiclassical field equations of fermionic matter in curved spacetime interacting with a background scalar field. The ultraviolet subtracting terms of the adiabatic regularization can be naturally interpreted as coming from appropriate counterterms of the background fields. We fix the required covariant counterterms. To test our approach we determine the contribution of the Yukawa interaction to the conformal anomaly in the massless limit and show its consistency with the heat kernel method using the effective action.

Equations

$$S_m=\int d^{4}x\sqrt{-g}\left\{\frac{i}{2}[\bar{\psi }{\underline{\gamma }}^{\mu }{\nabla }_{\mu }\psi -({\nabla }_{\mu }\bar{\psi }){\underline{\gamma }}^{\mu }\psi )]-m\bar{\psi }\psi -g_Y\Phi \bar{\psi }\psi \right\},$$

$$(i{\underline{\gamma }}^{\mu }{\nabla }_{\mu }-m-g_Y\Phi )\psi =0,$$

$$T_{\mu \nu }^m:=\frac{2}{\sqrt{-g}}\frac{\delta S_m}{\delta g^{\mu \nu }}=\frac{V_{\nu a}}{\det V}\frac{\delta S_m}{\delta V_a^{\mu }}=\frac{i}{2}\left[\bar{\psi }{\underline{\gamma }}_{(\mu }{\nabla }_{\nu )}\psi -({\nabla }_{(\mu }\bar{\psi }){\underline{\gamma }}_{\nu )}\psi \right].$$

$${\nabla }_{\mu }{T}_m^{\mu \nu }=g_Y\bar{\psi }\psi {\nabla }^{\nu }\Phi .$$

$${\nabla }_{\mu }{T}_m^{\mu \nu }+\frac{1}{\sqrt{-g}}\frac{\delta S_m}{\delta \Phi }{\nabla }^{\nu }\Phi =0,$$

$${\nabla }_{\mu }⟨T_m^{\mu \nu }⟩=g_Y⟨\bar{\psi }\psi ⟩{\nabla }^{\nu }\Phi .$$

$$S=S_g+S_m=\frac{1}{16\pi G}\int d^{4}x\sqrt{-g}R+\int d^{4}x\sqrt{-g}\left\{\frac{1}{2}{g}^{\mu \nu }{\nabla }_{\mu }\Phi {\nabla }_{\nu }\Phi -V(\Phi )\right\}+S_m,$$

$$G^{\mu \nu }+8\pi G({\nabla }^{\mu }\Phi {\nabla }^{\nu }\Phi -\frac{1}{2}{g}^{\mu \nu }{\nabla }^{\rho }\Phi {\nabla }_{\rho }\Phi +g^{\mu \nu }V(\Phi ))=-8\pi G{T}_m^{\mu \nu },$$

$$\square \Phi +\frac{\partial V}{\partial \Phi }=-g_Y\bar{\psi }\psi .$$

$$G^{\mu \nu }+8\pi G({\nabla }^{\mu }\Phi {\nabla }^{\nu }\Phi -\frac{1}{2}{g}^{\mu \nu }{\nabla }^{\rho }\Phi {\nabla }_{\rho }\Phi +g^{\mu \nu }V(\Phi ))$$

$$\square \Phi +\frac{\partial V}{\partial \Phi }$$

$${\nabla }_{\mu }({\nabla }^{\mu }\Phi {\nabla }^{\nu }\Phi -\frac{1}{2}{g}^{\mu \nu }{\nabla }^{\rho }\Phi {\nabla }_{\rho }\Phi +g^{\mu \nu }V(\Phi ))=(\square \Phi +\frac{\partial V}{\partial \Phi }){\nabla }^{\nu }\Phi ,$$

$${\nabla }_{\mu }⟨T_m^{\mu \nu }{⟩}_{ren}=g_Y⟨\bar{\psi }\psi {⟩}_{ren}{\nabla }^{\nu }\Phi =-(\square \Phi +\frac{\partial V}{\partial \Phi }){\nabla }^{\nu }\Phi .$$

$$\left({\partial }_0+\frac{3}{2}\frac{\dot{a}}{a}+\frac{1}{a}{\gamma }^0\vec{\gamma }\vec{\nabla }+i(m+s(t)){\gamma }^0\right)\psi =0,$$

$$\left({\partial }_t+\frac{3\dot{a}}{2a}+i{\gamma }^0\vec{\gamma }\frac{\vec{k}}{a}+i{\gamma }^0(m+s(t))\right){\psi }_{\vec{k}}=0.$$

$$\displaystyle\psi_{\vec{k},\lambda}(t)=\frac{1}{a^{3/2}(t)}\left({\begin{array}[]{c}h^{I}_{{k}}(t)\xi_{\lambda}(\vec{k})\\ h^{II}_{{k}}(t)\frac{\vec{\sigma}\vec{k}}{k}\xi_{\lambda}(\vec{k})\\ \end{array}}\right)\ ,$$

$$\xi_{+1}(\vec{k})=\frac{1}{\sqrt{2k(k+k_{3})}}\left({\begin{array}[]{c}k+k_{3}\\ k_{1}+ik_{2}\\ \end{array}}\right)\ ,\hskip 28.45274pt\xi_{-1}(\vec{k})=\frac{1}{\sqrt{2k(k+k_{3})}}\left({\begin{array}[]{c}-k_{1}+ik_{2}\\ k+k_{3}\\ \end{array}}\right)\ ,$$

$$h_k^{II}=\frac{ia}{k}\left(\frac{\partial h_k^I}{\partial t}+i(m+s(t))h_k^I\right),h_k^I=\frac{ia}{k}\left(\frac{\partial h_k^{II}}{\partial t}-i(m+s(t))h_k^{II}\right).$$

$$\displaystyle u_{\vec{k},\lambda}(t)=\frac{e^{i\vec{k}\vec{x}}}{\sqrt{(2\pi)^{3}a^{3}(t)}}\left({\begin{array}[]{c}h^{I}_{{k}}(t)\xi_{\lambda}(\vec{k})\\ h^{II}_{{k}}(t)\frac{\vec{\sigma}\vec{k}}{k}\xi_{\lambda}(\vec{k})\\ \end{array}}\right)\ .$$

$$\displaystyle v_{\vec{k},\lambda}(t)=Cu_{\vec{k},\lambda}(t)=\frac{e^{-i\vec{k}\vec{x}}}{\sqrt{(2\pi)^{3}a^{3}(t)}}\left({\begin{array}[]{c}h^{II*}_{{k}}(t)\xi_{-\lambda}(\vec{k})\\ h^{I*}_{{k}}(t)\frac{\vec{\sigma}\vec{k}}{k}\xi_{-\lambda}(\vec{k})\\ \end{array}}\right)\ .$$

$$|h_k^I{|}^2+|h_k^{II}{|}^2=1.$$

$$\psi (x)=\int d^3\vec{k}\sum _{\lambda }\left[B_{\vec{k}\lambda }{u}_{\vec{k}\lambda }(x)+D_{\vec{k}\lambda }^{†}{v}_{\vec{k}\lambda }(x)\right].$$

$$h_k^I(t)=\sqrt{\frac{\omega (t)+m}{2\omega (t)}}{e}^{-i{\int }^t\omega (t^{\prime })d{t}^{\prime }},h_k^{II}(t)=\sqrt{\frac{\omega (t)-m}{2\omega (t)}}{e}^{-i{\int }^t\omega (t^{\prime })d{t}^{\prime }},$$

$$h_k^I(t)=\sqrt{\frac{\omega (t)+m}{2\omega (t)}}{e}^{-i{\int }^t\Omega (t^{\prime })d{t}^{\prime }}F(t),h_k^{II}(t)=\sqrt{\frac{\omega (t)-m}{2\omega (t)}}{e}^{-i{\int }^t\Omega (t^{\prime })d{t}^{\prime }}G(t),$$

$$\Omega$$

$$F$$

$$G$$

$$(\omega -m)G=\Omega F+i\dot{F}+\frac{iF}{2}\frac{d\omega }{dt}\left(\frac{1}{\omega +m}-\frac{1}{\omega }\right)-(m+s)F,$$

$$(\omega +m)F=\Omega G+i\dot{G}+\frac{iG}{2}\frac{d\omega }{dt}\left(\frac{1}{\omega -m}-\frac{1}{\omega }\right)+(m+s)G,$$

$$(\omega +m)F{F}^{\ast }+(\omega -m)G{G}^{\ast }=2\omega .$$

$$f^{(1)}$$

$$f^{(2)}$$

$$(\omega -m)G^{(1)}=(\omega -m)F^{(1)}+{\omega }^{(1)}-s+\frac{i}{2}\frac{d\omega }{dt}\left(\frac{1}{\omega +m}-\frac{1}{\omega }\right),$$

Full text

Keywords: quantum field theory in curved spacetime, adiabatic regularization, Yukawa interaction, semiclassical gravity, cosmology, inflation, preheating.

I Introduction

A major problem in the theory of quantized fields in curved spacetimes parker-toms ; Waldbook ; fulling ; birrell-davies is the computation of the expectation values of the stress-energy tensor components. These calculations are rather convoluted, as they involve products of fields at coincident spacetime points, which are ultraviolet (UV) divergent even for free fields. In cosmological scenarios, this is connected to the fundamental phenomenon of particle creation by time-dependent backgrounds parker66 ; parker68 . A nonadiabatic expansion of the Universe generically induces particle creation of both bosonic and fermionic species, which leads to new UV divergences in their quadratic expectation values, not present in the absence of expansion.

A very efficient renormalization method, specifically constructed to deal with the UV divergences of a free field in an expanding universe, is adiabatic regularization. Originally, this technique was introduced to tame the divergences of the mean particle number of a scalar field in a Friedmann-Lemaitre-Robertson-Walker (FLRW) universe parker66 , and was later extended to get rid of the divergences of the stress-energy tensor parker-fulling ; Bunch80 ; Anderson-Parker . The key ingredient of the adiabatic scheme is the asymptotic expansion of the field modes, in which increasingly higher-order terms in the expansion involve increasingly higher-order time derivatives of the metric (the scale factor). Due to dimensional reasons, this is equivalent to an UV asymptotic expansion in momenta. This way, one can expand adiabatically the integrand of the unrenormalized bilinear, identify the UV-divergent terms, and subtract them directly to obtain a finite, covariant expre ssion. The renormalized expectation value is hence expressed as a finite integral in momentum space, depending exclusively on the mode functions defining the quantum state. The particular form of the adiabatic expansion depends on the spin of the quantized field. For free scalar fields, the well-known WKB expansion provides an adequate solution (see, for instance, parker-toms ; fulling ; birrell-davies ). For spin-1/2 fields, however, the adiabatic expansion takes a different form LNT ; RNT (see also Ghosh:2016 ). The adiabatic method has been proven to be equivalent to the DeWitt-Schwinger point-splitting scheme christensen76 ; Christensen78 for both scalar fields Birrell ; Anderson-Parker , and spin-$1/2$ fields rio1 . The method is specially suitable for numerical calculations, Hu-Parker ; Anderson and also for analytic approximations Anderson-Eaker ; Molina-Paris-Anderson-Ramsey . The adiabatic regularization has also been use d to scrutinize the two-point function defining the variance and power spectrum in inflationary cosmology and related issues inflation-r ; Woodard ; Markkanen . We would like to note that, in cosmological perturbation theory, other fields apart from the inflaton are regarded themselves as first order. Then, in order to stay at linear order in cosmological fluctuations at the level of the equations of motion, it is necessary to study these fields as propagating in the perfectly homogeneous FLRW background. This is what is customarily done in fact for the scalar curvature and tensor metric perturbations. Interactions between these metric perturbations and the additional fields are of higher order in cosmological perturbation theory. Therefore, when adiabatic regularization (or any other renormalization method) is applied to inflationary cosmology, the additional fields should be considered, at leading order, as quantum fields propagating in the homogeneous FLRW spacetime.

Particle creation also takes place if the quantized field, either a boson or a fermion, is coupled to a classical background scalar field evolving nonadiabatically in time. In this case, the interaction term acts in the boson/fermion equation of motion as a time-dependent effective mass, which excites the field and increases its mean particle number. The most paradigmatic example of this is probably preheating after inflation KLS . In the same way as before, new UV divergences associated to this particle creation appear in the expectation values of the different bilinears, not present in the absence of interaction, which must be appropriately removed to obtain finite quantities. Although renormalization of expectation values for interacting fields is generally much more complicated, adiabatic regularization can be generalized to include interactions to classical scalar background fields. In this case, the adiabatic expansion of the field modes used to identify the UV- divergent terms depends on both the scale factor and the background field, as well as their respective time derivatives. If the quantized field is a scalar with a Yukawa-type coupling, the adiabatic expansion is still of the WKB form Anderson-Molina-Paris-Evanich-Cook ; Molina-Paris-Anderson-Ramsey . However, a generalization of the adiabatic scheme for other interacting species in an expanding universe is absent in the literature. Here, we will try to partially fill in this gap.

In this work, we extend the adiabatic regularization method to Dirac fields living in a FLRW universe and interacting, via the standard Yukawa coupling, with an external scalar field. In this approach, the Dirac field is quantized, while both the metric and the background field are regarded as classical. This kind of system appears for example in fermionic preheating, in which the inflaton acts as a background scalar field oscillating around the minimum of its potential, and decays nonperturbatively into fermions due to its Yukawa interactions fermionic-preheating . Another example is the decay of the Standard Model (SM) Higgs after inflation, in which the Higgs condensate oscillates around the minimum of its potential, and transfers part of its energy into all the massive fermions of the Standard Model, coupled to the Higgs with the usual SM Yukawa couplings figueroa ; enqvist-meriniemi-nurmi (another part being transferred to the SM gauge bosons enqvist-meriniemi-nurmi ; SMHiggsLattice ). In certain models, the Higgs decay may also lead to the reheating of the Universe higgs-preheating . In this work, we will not focus on a particular scenario, but consider arbitrary time-dependent scale factors and background fields. The main objective is to provide well-motivated and rigorous expressions for the renormalized expectation values of the fermion stress-energy tensor $\langle T_{\mu\nu}\rangle$ and the bilinear $\langle\bar{\psi}\psi\rangle$. In the semiclassical equations of motion, these are the quantities that incorporate the backreaction of the created matter onto the background fields. To check the validity of the adiabatic method, we will also compute the contribution of the Yukawa interaction to the conformal anomaly in the massless limit, and check its consistency with the heat-kernel method using the effective action.

The paper is organized as follows. In Sec. II we give a general overview of the problem, introducing all necessary notation and equations of motion of the system. In Sec. III we develop the adiabatic expansion of the fermion field modes, subject to a Yukawa interaction with a scalar background field. In Sec. IV we derive general expressions for the renormalized expectation values of the stress-energy tensor $\langle T_{\mu\nu}\rangle$ and the bilinear $\langle\bar{\psi}\psi\rangle$. In Sec. V we further analyze the adiabatic regularization program by determining the covariant counterterms associated to the UV divergences. In Sec. VI we apply the method to calculate the conformal anomaly, and include a discussion concerning the ambiguity of the coefficients of the anomaly on the renormalization scheme. We also show that our results are compatible with the heat-kernel method. Finally, in Sec. VII we summarize our results and conclude. The paper is accompanied by three appendixes. In Appendix A we compute the conformal anomaly of a quantized scalar field coupled to a scalar background field. In Appendix B we apply adiabatic regularization to a simple, analytically solvable example. Finally, in Appendix C we gather the terms of the fermionic adiabatic expansion up to fourth order.

In this work, we take the FLRW metric as $ds^{2}=dt^{2}-a^{2}(t)d\vec{x}^{2}$, and we use the Dirac-Pauli representation for the Dirac gamma matrices, \gamma^{0}=\scriptsize\left({\begin{array}[]{cc}I&0\\ 0&-I\\ \end{array}}\right), \vec{\gamma}=\scriptsize\left({\begin{array}[]{cc}0&\vec{\sigma}\\ -\vec{\sigma}&0\\ \end{array}}\right), with $\vec{\sigma}$ the usual Pauli matrices. We also assume natural units $\hbar=1=c$.

II Semiclassical equations for a quantized Dirac matter field with Yukawa coupling

We consider the theory defined by the action functional $S=S[g_{\mu\nu},\Phi,\psi,\nabla\psi]$, where $\psi$ represents a Dirac field, $\Phi$ is a scalar field, and $g_{\mu\nu}$ stands for the spacetime metric. We decompose the action as $S=S_{g}+S_{m}$, where $S_{m}$ is the matter sector

[TABLE]

and $S_{g}$ is the gravity-scalar sector, which will be presented in the next subsection. Here, $\underline{\gamma}^{\mu}(x)$ are the spacetime-dependent Dirac matrices satisfying the anticommutation relations $\{\underline{\gamma}^{\mu},\underline{\gamma}^{\nu}\}=2g^{\mu\nu}$, related to the usual Minkowski ones by the vierbein field $V_{\mu}^{a}(x)$ defined through $g_{\mu\nu}(x)=V_{\mu}^{a}(x)V_{\nu}^{b}(x)\eta_{ab}$. On the other hand, $\nabla_{\mu}\equiv\partial_{\mu}-\Gamma_{\mu}$ is the covariant derivative associated to the spin connection $\Gamma_{\mu}$, $m$ is the mass of the Dirac field, and $g_{Y}$ is the dimensionless coupling constant of the Yukawa interaction. In (1), both the metric $g_{\mu\nu}(x)$ and the scalar field $\Phi(x)$ are regarded as classical external fields. The Dirac spinor $\psi(x)$ will be our quantized field living in a curved spacetime and possessing a Yukawa coupling to the classical field $\Phi$. The Dirac equation is

[TABLE]

and the stress-energy tensor is given by birrell-davies

[TABLE]

The presence of the Yukawa interaction with the external field $\Phi$ modifies the standard conservation equation. We have, instead,

[TABLE]

These equations can be easily seen as the consequence of the invariance of the action functional $S$ under spacetime diffeomorphisms $\delta x^{\mu}=\epsilon^{\mu}(x)$: $\delta\Phi=\epsilon^{\mu}\nabla_{\mu}\Phi$, $\delta g_{\mu\nu}=2\nabla_{(\mu}\epsilon_{\nu)}$. One gets

[TABLE]

which reproduces (4). We will assume that the quantum theory fully respects this symmetry. Therefore, we demand

[TABLE]

II.1 Adding the gravity-scalar sector

The complete theory, including the gravity-scalar sector in the action, can be described by

[TABLE]

where $S_{m}$ is the action for the matter sector given in (1). We will reconsider the form of the action in Sec. V, in view of the counterterms required to cancel the UV divergences of the quantized Dirac field. However, let us work for the moment with the action (7). The Einstein equations are then

[TABLE]

and the equation for the scalar field is

[TABLE]

The semiclassical equations are obtained from (8) and (9) by replacing $T_{m}^{\mu\nu}$ and $\bar{\psi}\psi$ by the corresponding (renormalized) vacuum expectation values $\langle T_{m}^{\mu\nu}\rangle_{ren}$ and $\langle\bar{\psi}\psi\rangle_{ren}$,

[TABLE]

These equations are consistent with the Bianchi identities $\nabla_{\mu}G^{\mu\nu}=0$ , since

[TABLE]

and, from (6) and (11), we have

[TABLE]

When the spacetime is an expanding universe [$ds^{2}=dt^{2}-a^{2}(t)d\vec{x}^{2}$], and $\Phi$ is a homogeneous scalar field $\Phi=\Phi(t)$ (e.g. an inflaton), Eqs. (10) and (11) describe the backreaction on the metric-inflaton system due to matter particle production and vacuum polarization codified in the renormalized vacuum expectation values $\langle T_{m}^{\mu\nu}\rangle_{ren}$ and $\langle\bar{\psi}\psi\rangle_{ren}$. It is then important to elaborate an efficient method to compute these quantities in this cosmological setting.

III Adiabatic expansion for a Dirac field with Yukawa coupling

In a spatially flat FLRW spacetime, the time-dependent gamma matrices are related with the Minkowskian ones by $\underline{\gamma}^{0}(t)=\gamma^{0}$ and $\underline{\gamma}^{i}(t)=\gamma^{i}/a(t)$, and the components of the spin connections are $\Gamma_{0}=0$ and $\Gamma_{i}=(\dot{a}/2)\gamma_{0}\gamma_{i}$. The Dirac equation with the Yukawa interaction $i\underline{\gamma}^{\mu}\nabla_{\mu}\psi-m\psi=g_{Y}\Phi\psi$, taking $\Phi$ as a homogenous scalar field $\Phi=\Phi(t)$, is then

[TABLE]

where we have defined $s(t)\equiv g_{Y}\Phi(t)$. If we expand the field $\psi$ as $\psi=\int\frac{d^{3}\vec{k}}{(2\pi)^{3/2}}\psi_{\vec{k}}(t)e^{i\vec{k}\vec{x}}$, and we substitute it into (14), we obtain the following differential equation for $\psi_{\vec{k}}$:

[TABLE]

In order to solve this equation, it is convenient to write the Dirac field in terms of two two-component spinors of the generic form

[TABLE]

where $\xi_{\lambda}$ with $\lambda={\pm}1$ are two constant orthonormal two-spinors ($\xi_{\lambda}^{\dagger}\xi_{\lambda^{\prime}}=\delta_{\lambda,\lambda^{\prime}}$), eigenvectors of the helicity operator $\frac{\vec{\sigma}\vec{k}}{2k}\xi_{\lambda}=\frac{\lambda}{2}\xi_{\lambda}$. The explicit forms of $\xi_{+1}$ and $\xi_{-1}$ are

[TABLE]

where $\vec{k}=(k_{1},k_{2},k_{3})$ and $|\vec{k}|=k$. The time-dependent functions $h_{k}^{I}$ and $h_{k}^{II}$ satisfy the first-order coupled equations

[TABLE]

Given a particular solution $\{h_{k}^{I}(t)$, $h_{k}^{II}(t)\}$ to Eqs. (20), one can construct the modes

[TABLE]

Equation (23) will be a solution of positive-frequency type in the adiabatic regime. A solution of negative-frequency type can be obtained by applying a charge conjugate transformation $C\psi=-i\gamma^{2}\psi^{*}$ (we follow here the convention in Peskin-Schroeder )

[TABLE]

The Dirac inner product is defined as $(\psi_{1},\psi_{2})=\int d^{3}xa^{3}\psi_{1}^{\dagger}\psi_{2}$. The normalization condition for the above four-spinors, $(u_{\vec{k}\lambda},v_{\vec{k}\,^{\prime}\lambda^{\prime}})=0$, $(u_{\vec{k}\lambda},u_{\vec{k}\,^{\prime}\lambda^{\prime}})=(v_{\vec{k}\lambda},v_{\vec{k}\,^{\prime}\lambda^{\prime}})=\delta_{\lambda\lambda^{\prime}}\delta^{(3)}(\vec{k}-\vec{k}\,^{\prime})$, reduces to

[TABLE]

Since the Dirac scalar product is preserved by the cosmological evolution, the normalization condition (27) holds at any time. This ensures also the standard anticommutation relations for the creation and annihilation operators [$\{B_{\vec{k},\lambda},B_{\vec{k^{\prime}},\lambda^{\prime}}^{\dagger}\}=\delta^{3}(\vec{k}-\vec{k^{\prime}})\delta_{\lambda\lambda^{\prime}}$, $\{B_{\vec{k},\lambda},B_{\vec{k^{\prime}},\lambda^{\prime}}\}=0$, and similarly for the $D_{\vec{k},\lambda}$, $D_{\vec{k^{\prime}},\lambda^{\prime}}^{\dagger}$ operators], defined by the Fourier expansion of the Dirac field operator

[TABLE]

III.1 Adiabatic expansion

We now compute the adiabatic expansion of a Dirac field living in a FLRW spacetime, and possessing a Yukawa interaction term with a classical background field. We know that, in the adiabatic limit, and in the absence of interaction, the natural solution of the field modes $h_{k}^{I}$ and $h_{k}^{II}$ is

[TABLE]

where $\omega=\sqrt{\frac{k^{2}}{a^{2}}+m^{2}}$ is the frequency of the field mode. This will constitute the zeroth-order term of the adiabatic expansion. Mimicking the ansatz introduced in LNT , we write the $h_{k}^{I}$ and $h_{k}^{II}$ functions as

[TABLE]

where $\Omega(t)$, $F(t)$, and $G(t)$ are time-dependent functions, which we expand adiabatically as

[TABLE]

Here, $F^{(n)}$, $G^{(n)}$, and $\omega^{(n)}$ are terms of $n$th adiabatic order (we explain exactly what we mean by that below). By substituting (30) into the equations of motion (20) and the normalization condition (27), we obtain the following system of three equations,

[TABLE]

To obtain the expressions for $\Omega^{(n)}$, $F^{(n)}$, and $G^{(n)}$, we introduce the adiabatic expansions (31) into (32), and solve order by order. As usual, we consider $\dot{a}$ of adiabatic order one, $\ddot{a}$ of adiabatic order two, and so on. On the other hand, we consider the interaction term $s(t)$ of adiabatic order one, so that the zeroth-order term in (30) recovers the free field solution in the adiabatic limit, defined in (29). Similarly, time derivatives of the interaction increase the adiabatic order, so that $\dot{s}$ is of order two, $\ddot{s}$ of order three, and so on. With this, a generic expression $f^{(n)}$ of adiabatic order $n$ (e.g. $f^{(n)}=F^{(n)},G^{(n)},\Omega^{(n)}$) will be written as a sum of all possible products of $n$th adiabatic order formed by $s$, $a$, and their time derivatives. For example, functions of adiabatic orders one and two will be written respectively as

[TABLE]

with $\alpha_{n}\equiv\alpha_{n}(m,k,a)$ and $\beta_{n}\equiv\beta_{n}(m,k,a)$. The assignment of $s$ as adiabatic order one is consistent with the scaling dimension of the scalar field, as it possesses the same dimensions as $\dot{a}$.

III.1.1 First adiabatic order

By keeping only terms of first adiabatic order in (32), the system of three equations gives

[TABLE]

We now treat independently the real and imaginary parts by writing $F^{(1)}=f_{x}^{(1)}+if_{y}^{(1)}$ and $G^{(1)}=g_{x}^{(1)}+ig_{y}^{(1)}$. We obtain for the real part

[TABLE]

which has as solutions

[TABLE]

On the other hand, the imaginary part of the system gives

[TABLE]

These two equations are not independent. The obtained solution for $g_{y}^{(1)}$ and $f_{y}^{(1)}$ is

[TABLE]

where $A$ is an arbitrary first-order adiabatic function. We will choose the simplest solution

[TABLE]

obeying the condition $F^{(1)}(m,s)=G^{(1)}(-m,-s)$. Therefore, the adiabatic expansion will also preserve the symmetries of Eqs. (20) with respect to the change $(m,s)\rightarrow(-m,-s)$. We have checked that physical expectation values are independent to any potential ambiguity in this kind of choice.

III.1.2 Second adiabatic order

In the same way, the second-order terms of (32) give

[TABLE]

where the first-order terms have already been deduced above. Taking the real part of these equations, we obtain

[TABLE]

which has as solutions

[TABLE]

and $g^{(2)}_{x}(m,s)=f^{(2)}_{x}(-m,-s)$. On the other hand, taking the imaginary part of the equations, we have

[TABLE]

As before, this system contains an arbitrariness in its solution,

[TABLE]

where now $B$ is a linear combination of second-order adiabatic terms. By imposing again the condition $F^{(2)}(m,s)=G^{(2)}(-m,-s)$, one finds

[TABLE]

and $g_{y}^{(2)}(m,s)=f_{y}^{(2)}(-m,-s)$.

III.1.3 Third and fourth adiabatic order

The same procedure can be repeated for all orders. The real part of the expansion is totally determined by the system of equations (32), while every imaginary part contains an arbitrariness that can be solved by fixing the condition $F^{(n)}(m,s)=G^{(n)}(-m,-s)$. The third- and fourth-order terms of the expansion are explicitly written in Appendix C.

IV Renormalization of the stress-energy tensor $\langle T_{\mu\nu}\rangle$ and the bilinear $\langle\bar{\psi}\psi\rangle$

The classical stress-energy tensor in a FLRW spacetime has two independent components. For a Dirac field, they are (no sum on $i$),

[TABLE]

We define the vacuum state $|0\rangle$ as $B_{\vec{k},\lambda}|0\rangle\equiv D_{\vec{k},\lambda}|0\rangle\equiv 0$, and denote any expectation value on this vacuum as e.g. $\langle T_{\mu\nu}\rangle\equiv\langle 0|T_{\mu\nu}|0\rangle$. In the quantum theory, the vacuum expectation values of the stress-energy tensor take the form (see for example RNT )

[TABLE]

and

[TABLE]

The above formal expressions contain quartic, quadratic, and logarithmic UV divergences, which turn out to be independent of the particular quantum state. These divergences are similar to those described in Utiyama-DeWitt . To characterize them, one plugs in (47)-(48) the adiabatic expansion of $h_{k}^{I}$ and $h_{k}^{II}$, given in Eq. (30). We shall see that, in the presence of a Yukawa interaction, all adiabatic orders up to the fourth one generate UV divergences. This is different to what happens in the case of a free field, where the divergences only appear at zeroth and second adiabatic orders RNT . In general, adiabatic renormalization proceeds by subtracting those adiabatic terms from the integrand of the expectation values, producing a formal finite quantity. There are two important considerations regarding these subtractions. First, they must refer to all contributions of a given adiabatic term of fixed (adiabatic) order, othe rwise general covariance is not maintained. And second, one subtracts only the minimum number of terms required to get a finite result parker-toms .

We now proceed to calculate the renormalized expressions for the energy density and pressure.

IV.1 Renormalized energy density

We start by performing the adiabatic expansion of the energy density in momentum space (48)

[TABLE]

where $\rho_{k}^{(n)}$ is of $n$th adiabatic order. The adiabatic terms producing UV divergences (after integration in momenta) are

[TABLE]

where we have used the notation $a^{(3)}\equiv d^{3}a/dt^{3},a^{(4)}\equiv d^{4}a/dt^{4}$, etc.

We note that if we turn off the Yukawa coupling, we recover the results obtained in RNT . The Yukawa interaction produces new contributions and, in particular, we have now non-zero terms at first and third adiabatic orders. The physical meaning of them will be given later on. Note here that in the UV limit, $\rho_{k}^{(0)}\sim k$, $(\rho_{k}^{(1)}+\rho_{k}^{(2)})\sim k^{-1}$, and $(\rho_{k}^{(3)}+\rho_{k}^{(4)})\sim k^{-3}$. This indicates that subtracting the zeroth-order term will cancel the natural quartic divergence of the stress-energy tensor, subtracting up to second order will cancel also the quadratic divergence, and subtracting up to fourth order will cancel the logarithmic divergence. Therefore, defining the adiabatic subtraction terms as

[TABLE]

the renormalized 00 component of the stress-energy tensor is

[TABLE]

This integral is, by construction, finite.

IV.2 Renormalized pressure

The method proceeds in the same way for the pressure. The renormalized ii component of the stress-energy tensor is given by

[TABLE]

where $p_{k}^{(0-4)}\equiv p_{k}^{(0)}+p_{k}^{(1)}+p_{k}^{(2)}+p_{k}^{(3)}+p_{k}^{(4)}$, and

[TABLE]

The corresponding adiabatic terms for the pressure are

[TABLE]

As before, we see that in the UV limit, $p_{k}^{(0)}\sim k$, $(p_{k}^{(1)}+p_{k}^{(2)})\sim k^{-1}$, and $(p_{k}^{(3)}+p_{k}^{(4)})\sim k^{-3}$. Subtracting the zeroth-order term eliminates the quartic divergence, subtracting up to second order removes the quadratic divergence, and subtracting up to fourth order removes the logarithmic divergence. If the Yukawa interaction is removed, we recover again the results in RNT . The interaction produces also nonzero contributions to the first and third adiabatic orders.

IV.3 Renormalization of $\langle\bar{\psi}\psi\rangle$

We are also interested in computing the renormalized expectation value $\langle\bar{\psi}\psi\rangle_{ren}$. The formal (unrenormalized) expression for this quantity is

[TABLE]

We define the corresponding terms in the adiabatic expansion as $\langle\bar{\psi}\psi\rangle_{k}=\langle\bar{\psi}\psi\rangle_{k}^{(0)}+\langle\bar{\psi}\psi\rangle_{k}^{(1)}+\langle\bar{\psi}\psi\rangle_{k}^{(2)}+\langle\bar{\psi}\psi\rangle_{k}^{(3)}+...$ . Due to the Yukawa interaction, ultraviolet divergences arrive till the third adiabatic order. In general, we have

[TABLE]

From here, we obtain

[TABLE]

The adiabatic prescription leads then to

[TABLE]

In this case, we observe that in the UV limit, $(\langle\bar{\psi}\psi\rangle_{k}^{(0)}+\langle\bar{\psi}\psi\rangle_{k}^{(1)})\sim k^{-1}$, and $(\langle\bar{\psi}\psi\rangle_{k}^{(2)}+\langle\bar{\psi}\psi\rangle_{k}^{(3)})\sim k^{-3})$. Subtracting up to first order eliminates the quadratic divergence, and up to third order removes the logarithmic one.

Our results can be generically implemented together with numerical methods to compute the renormalized expectation values $\langle T_{\mu\nu}\rangle_{ren}$ and $\langle\bar{\psi}\psi\rangle_{ren}$. On the other hand, we would like to briefly comment that a higher-order adiabatic expansion also serves to generate asymptotic analytical expressions for the renormalized stress-energy tensor in some special situations. This happens in spacetime regions where the relevant modes always evolve adiabatically. For instance, if we approximate the form of the exact modes $\{h_{k}^{I},h_{k}^{II}\}$ by their higher-order adiabatic expansion, we can find in a very straightforward way an analytic approximation for the renormalized quantities in the adiabatic regime, as in the example given in Appendix B. Outside the adiabatic regime one should use numerical methods to find the exact modes and plug them in the generic renormalized expressions obtained above.

V Ultraviolet divergences and renormalization counterterms

The ultraviolet divergent terms of the adiabatic subtractions can be univocally related to particular counterterms in a Lagrangian density including the background gravity-scalar sector. By writing

[TABLE]

the equations of motion for the scalar field are

[TABLE]

From (70), we can write the identity

[TABLE]

where $\langle\bar{\psi}\psi\rangle_{ren}$ is finite and the remaining integrals at the right-hand side of (73) are the adiabatic subtraction terms. As we shall see, the ultraviolet divergences of the adiabatic subtraction terms can be removed by counterterms of the form: $\delta{Z}\Box\Phi$, $\delta{\lambda_{1}}$, $\delta{\lambda_{2}}\Phi$, $\delta{\lambda_{3}}\Phi^{2}$, $\delta{\lambda_{4}}\Phi^{3}$, $\delta{\xi_{1}}R$, and $\delta{\xi_{2}}R\Phi$. To deal with the UV-divergent subtraction terms we use dimensional regularization Bunch80 . We can check that the (covariantly) regulated divergences take the same form as the above covariant counterterms. For $\langle\bar{\psi}\psi\rangle^{(0)}$ we have ($n$ denotes the spacetime dimension)

[TABLE]

where we will retain only the poles at $n=4$. This divergence can be absorbed by $\delta\lambda_{1}$. Additionally, we also have

[TABLE]

This divergence of adiabatic order one can be absorbed by $\delta\lambda_{2}$. The divergences of adiabatic order two

[TABLE]

can also be eliminated by $\delta\xi_{1}$ and $\delta\lambda_{3}$. Finally, the three divergences of adiabatic order three

[TABLE]

are absorbed by $\delta Z$, $\delta\xi_{2}$ and $\delta\lambda_{4}$.

On the other hand, the tensorial equations are

[TABLE]

and we find similar cancellations. However, two extra divergences appear. Focusing, for simplicity, at zeroth adiabatic order, we have

[TABLE]

At first adiabatic order we encounter the following divergences

[TABLE]

At second adiabatic order we find these divergences

[TABLE]

At third adiabatic order we get the following divergences ($H\equiv\dot{a}/a$)

[TABLE]

Finally, at fourth adiabatic order the divergences are

[TABLE]

All the above divergent expressions arising from the Yukawa interaction can be written covariantly as

[TABLE]

and can be consistently removed (including also the divergences for $\langle\bar{\psi}\psi\rangle$) by the renormalization parameters

[TABLE]

[TABLE]

[TABLE]

We remark that the set of needed counterterms is all possible counterterms having couplings with non-negative mass dimension, up to Newton’s coupling constant. This is also in agreement with the results in perturbative Quantum Field Theory in flat spacetime. The renormalizability of the Yukawa interaction $g_{Y}\varphi\bar{\psi}\psi$ of a quantized massive scalar field $\varphi$ with a massive quantized Dirac field $\psi$ requires us to add terms of the form $\frac{\lambda Z_{\lambda}}{4!}\varphi^{4}$, $\frac{\kappa Z_{\kappa}}{3!}\varphi^{3}$, and also a term linear in $\varphi$ Srednicki . The presence of a curved background would require us to add the terms $\xi_{1}R\varphi$ and $\xi_{2}R\varphi^{2}$. We note that a term of the form $\xi_{2}R\varphi^{2}$ is required by renormalization for a purely quantized scalar field $\varphi$ if a self-interaction term of the form $\frac{\lambda}{4!}\varphi^{4}$ appears in the bare Lagrangian density bunch-parker ; bunch81 . Here we have found that the Yukawa interaction demands the presence of the renormalized terms $\xi_{1}R\varphi$ and $\xi_{2}R\varphi^{2}$ (as well as the terms $\lambda_{i}\varphi^{i}$), even if they are not present in the bare Lagrangian density. Similar counterterms have been identified in the approach in Ref. Baacke-Patzold .

Therefore, the tentative semiclassical equations presented in Sec. II should be reconsidered to include the above-required counterterms. In terms of the renormalized parameters we have

[TABLE]

and

[TABLE]

where the potential $V(\Phi)$ should contain the terms

[TABLE]

Obviously, additional terms, not required by renormalization, can be added to the potential if one adopts an effective field theory viewpoint. Some of the renormalized parameters ($\Lambda$, $\xi_{1}$, $\lambda_{1}$, $\cdots$) could take, by fine-tuning, zero values. We do not consider these issues in this work.

VI Conformal anomaly

In this section we will analyze the massless limit of the theory and work out the conformal anomaly. In the massless limit the classical action of the theory enjoys invariance under the conformal transformations

[TABLE]

with

[TABLE]

Variation of the action yields the identity

[TABLE]

which, in our case, turns out to be $g^{\mu\nu}T_{\mu\nu}-g_{Y}\Phi\bar{\psi}\psi=0$. At the quantum level the theory will lose its conformal invariance as a consequence of renormalization [which respects general covariance and hence (6)] and generates an anomaly

[TABLE]

$C_{f}$ is independent of the quantum state and depends only on local quantities of the external fields.

To calculate the conformal anomaly in the adiabatic regularization method, we have to start with a massive field and take the massless limit at the end of the calculation. Therefore,

[TABLE]

Since the divergences of the stress-energy tensor have terms of fourth adiabatic order, the adiabatic subtractions for $\langle\bar{\psi}\psi\rangle$ should also include them. The fourth-order subtraction term, which produces a nonzero finite contribution when $m\to 0$, is codified in $\langle\bar{\psi}\psi\rangle^{(4)}$. The term $m\langle\bar{\psi}\psi\rangle_{ren}$ vanishes when $m\to 0$. The remaining piece produces the anomaly [recall (64)-(65)]

[TABLE]

Applying the adiabatic expansion computed in Sec. III and doing the integrals we obtain

[TABLE]

Since $C_{f}$ is a scalar, we must be able to rewrite the above result as a linear combination of covariant scalar terms made out of the metric, the Riemann tensor, covariant derivatives, and the external scalar field $\Phi$. Our result is

[TABLE]

The same result is obtained by using the results of Sec. IV. $C_{f}$ can be reexpressed as [recall (56)-(57)]

[TABLE]

Performing the integrals we get exactly (101) and hence (102).

In Appendix A we have computed the conformal anomaly for a massless scalar field $\phi$ with conformal coupling to the scalar curvature $\xi=1/6$, and with a Yukawa-type interaction of the form $g_{Y}^{2}\Phi^{2}\phi^{2}$. Adiabatic regularization predicts the following conformal anomaly

[TABLE]

In the absence of Yukawa interaction ($h=0$, $g_{Y}=0$) we reproduce the well-known trace anomaly for both scalar and spin-$1/2$ fields (restricted to our FLRW spacetime) birrell-davies . We recall that the trace anomaly is generically given for a conformal free field of spin $0,1/2$, or $1$ in terms of three coefficients

[TABLE]

where $C_{\mu\nu\rho\sigma}$ is the Weyl tensor and $G=R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma}-4R_{\mu\nu}R^{\mu\nu}+R^{2}$ is proportional to the Euler density. The coefficients $a$ and $b$ are independent of the renormalization scheme and are given by Duff ; HHR

[TABLE]

where $N_{s}$ is the number of real scalar fields, $N_{f}$ is the number of Dirac fields, and $N_{v}$ is the number of vector fields. Our results with $g_{Y}=0$ fit the values in (VI). [We note that in the FLRW spacetime of adiabatic regularization the Weyl tensor vanishes identically]. In contrast, the coefficient $c$ depends in general on the particular renormalization scheme Wald78 . A local counterterm proportional to $R^{2}$ in the action can modify the coefficient $c$. For instance, for vector fields the point splitting and the dimensional regularization method predict different values for $c$.

When the Yukawa interaction is added, the general form of the conformal anomaly is

[TABLE]

Now, the coefficients $f$ and $g$ are independent of the renormalization scheme but $d$ and $e$ are not. The finite Lagrangian counterterms required by the renormalizability of the Yukawa interaction obtained in previous sections,

[TABLE]

might alter the values of the coefficients $d$ and $e$, but not the coefficients $f$ and $g$. Note that, due to classical conformal invariance, one should consider only those counterterms having dimensionless coupling parameters. Therefore, our results for the $f$ and $g$ coefficients are

[TABLE]

Finally, to show explicitly that the above coefficients are independent of the renormalization scheme, we will compute them using the heat-kernel method given in parker-toms , by means of the one-loop effective action.

VI.1 Consistency with the heat-kernel results

The conformal anomaly for a field $\phi^{j}(x)$ obeying the second-order wave equation

[TABLE]

is given by

[TABLE]

where $E_{2}(x)$ is the second Seeley-DeWitt coefficient. The minus sign is for bosons and the plus sign is for fermions. These coefficients are local, scalar functions of $Q(x)$ and the curvature tensor. $E_{2}$ is given by

[TABLE]

where $W_{\mu\nu}=[\nabla_{\mu},\nabla_{\nu}]$. For a single massless scalar field with $\xi=1/6$ and an interaction of the form $h^{2}\phi^{2}\Phi^{2}$ we have

[TABLE]

and $W_{\mu\nu}=0$. For a spatially flat FLRW universe we get

[TABLE]

in full agreement with the result (104) obtained using adiabatic regularization.

For a single massless Dirac field with a Yukawa interaction we have [$(i\gamma^{\mu}\nabla_{\mu}-g_{Y}\Phi)\psi=0$]

[TABLE]

and

[TABLE]

Using the properties of the trace of products of gamma matrices, we get

[TABLE]

The above result reproduces the coefficients $f$ and $g$ obtained from adiabatic regularization. We note that there is a mismatch in the coefficients $d$ and $e$. These are, however, the coefficients that might depend on the renormalization scheme.

VII Summary

When a quantum field is coupled to a classical, nonadiabatic time-dependent background, it gets excited, and undergoes a regime of particle creation. In this case, new UV-divergent terms appear in the expectation values of its quadratic products, which must be appropriately removed to obtain a physical, finite quantity. In cosmological scenarios, adiabatic regularization provides an appropriate solution to this challenge: by means of an adiabatic expansion of the field modes, one can identify the covariant UV-divergent terms of the corresponding bilinear, and subtract them directly from the unrenormalized quantity. The background may be the expansion of the Universe itself, as in the case of inflation, or a classical homogeneous scalar field. The adiabatic scheme can be applied in both situations, for both bosonic and fermionic species.

In this work, we have developed the adiabatic regularization method for spin-1/2 fields in an expanding universe, coupled to a classical background scalar field with a Yukawa interaction term. The results of this work are a natural generalization of the studies initiated in LNT ; RNT , and broaden significantly the range of applicability of the adiabatic method. We have computed the adiabatic expansion of the spin-1/2 field modes up to fourth adiabatic order, and used it to obtain expressions for the renormalized expectation values of the stress-energy tensor $\langle T_{\mu\nu}\rangle_{ren}$ and the bilinear $\langle\bar{\psi}\psi\rangle_{ren}$. These quantities are fundamental ingredients in the study of the semiclassical equations of fermionic matter interacting with a background field, as they codify the backreaction effects from the created matter on the metric/background fields. Therefore, it is essential to develop an efficient renormalization scheme to correctly quantify the effects of this backreaction. All expressions obtained are generic, depending only on the background scalar field and scale factor time-dependent functions. This constitutes probably the major advantage of the adiabatic renormalization scheme. We leave the method prepared to perform numerical computations in future investigations.

Finally, we have tested the overall theoretical construction of the adiabatic scheme by justifying the method in terms of renormalization of coupling constants, as well as by computing the conformal anomaly. Our calculation of the conformal anomaly with the Yukawa interaction has been proved to be fully consistent with the generic results obtained via the one-loop effective action. Therefore, by considering such a system, we have also improved our general understanding of quantum field theory in curved spacetimes.

Acknowledgments

F.T. thanks Daniel G. Figueroa and Juan Garcia-Bellido for useful discussions. This work is supported by the Grants. No. FIS2014-57387-C3-1-P , No. FPA2015-68048-C3-3-P, No. MPNS of COST Action No. CA15117, and the Severo Ochoa Programs SEV-2014-0398 and SEV-2012-0249. A.d.R. is supported by the FPU Ph.D. fellowship FPU13/04948, A. F. is supported by the Severo Ochoa Ph.D. fellowship SEV-2014-0398-16-1, and F.T. is supported by the Severo Ochoa Ph.D. fellowship SVP-2013-067697.

Appendix A Scalar field with a Yukawa-type coupling

In this appendix we compute the conformal anomaly of a quantized real scalar field $\phi$, coupled to another background scalar $\Phi$ with a Yukawa-type interaction. This result will be used in Sec. VI. The interaction term can be chosen of the form $g\Phi\phi^{2}$ or $h^{2}\Phi^{2}\phi^{2}$. Although the adiabatic regularization can be equally applied in both cases, we will focus on the latter case, since the coupling constant $h^{2}$ is dimensionless and the classical theory inherits the conformal invariance. Therefore, the action functional of the scalar matter field is given by

[TABLE]

As before, the scalar field lives in a spatially flat FLRW metric $ds^{2}=dt^{2}-a^{2}(t)d\vec{x}^{2}$, and we assume that the external field is homogeneous $\Phi=\Phi(t)$. In this case, the equation of motion is

[TABLE]

where we have introduced the notation $s(t)\equiv h\Phi(t)$, similar to the one used for the spin-$1/2$ field in the main text. The quantized field is expanded in Fourier modes as

[TABLE]

where $f_{\vec{k}}(x)=e^{i\vec{k}\vec{x}}h_{k}(t)$, and $A_{\vec{k}}^{\dagger}$ and $A_{\vec{k}}$ are the usual creation and annihilation operators. Substituting (120) into (119) we find

[TABLE]

where $\sigma(t)=(6\xi-\frac{3}{4})(\frac{\dot{a}^{2}}{a^{2}})+(6\xi-\frac{3}{2})(\frac{\ddot{a}}{a})$, and $\omega_{k}(t)=\sqrt{\frac{k^{2}}{a(t)^{2}}+m^{2}}$. The adiabatic expansion for the scalar field modes is based on the usual WKB ansatz

[TABLE]

which satisfies automatically the Wronskian condition $h_{k}\dot{h}_{k}^{*}-h_{k}^{*}\dot{h}_{k}=2i$. One can substitute the ansatz into Eq. (121), and solve order by order to obtain the different terms of the expansion. The function $W_{k}(t)$ obeys the differential equation

[TABLE]

Note that here, $s(t)\equiv h\Phi$ is assumed of adiabatic order one as in the fermionic case. One obtains systematically $\omega^{(odd)}=0$ for all terms of odd order in the expansion. At second adiabatic order one gets

[TABLE]

and at fourth adiabatic order, the result is

[TABLE]

Expressions for the subtraction terms in conformal time have been obtained in Molina-Paris-Anderson-Ramsey . Here we will briefly sketch the renormalization counterterms associated to the UV divergences of the stress-energy tensor and the variance $\langle\phi^{2}\rangle$. We follow a strategy similar to the one used in Sec. V. The Lagrangian density with the required renormalization counterterms is

[TABLE]

One can check that the above counterterms are enough to absorb all the UV divergences that emerge in the quantization of the scalar field. We note that, due to the symmetry $\Phi\to-\Phi$ of the matter Lagrangian, counterterms of the form $R\Phi$, $\Phi$, $\Phi^{3}$ are absent. However, a higher-derivative term of the form $R^{2}$ is now necessary, which did not appear for the Dirac field in a FLRW spacetime.

We assume the conformal coupling to the curvature $\xi=1/6$. For a massive field we have

[TABLE]

Classical conformal invariance is obtained when $m^{2}=0$. In adiabatic regularization the conformal anomaly is computed by taking the massless limit

[TABLE]

where $(W_{k}^{-1}(t))^{(4)}=\omega^{-3}(\omega^{(2)})^{2}-\omega^{-2}\omega^{(4)}$ is the fourth-order term in the adiabatic expansion of $W_{k}^{-1}$. Note that here, $\langle\phi^{2}\rangle^{(4)}$ is evaluated including fourth-order adiabatic subtractions. This is different to the physical vacuum expectation value $\langle\phi^{2}\rangle_{ren}$, which has to be evaluated with subtractions only up to second order. This is why only the purely fourth-order adiabatic piece contributes to the anomaly. The explicit expression of $(W_{k}^{-1})^{(4)}$ for arbitrary $\xi$ is

[TABLE]

The integral in comoving momenta is finite and independent of the mass. Assuming now $\xi=1/6$, the result is

[TABLE]

We can rewrite the expression in terms of covariant scalar terms as

[TABLE]

which is the result given in Eq. (104).

Appendix B A simple example

In this appendix we consider a simple mathematical example to illustrate how the adiabatic method works. We compute the bilinear $\langle\bar{\psi}\psi\rangle_{\rm ren}$ of a Dirac field, coupled to a background scalar field evolving in Minkowski spacetime [$a(t)=1$] as

[TABLE]

For convenience, we have absorbed the Yukawa coupling $g_{Y}$ in the dimensionless constant $\mu$. To avoid the mathematical instability at $t\rightarrow 0$, we will only consider times in the range $-\infty<t<0$. This model has three convenient aspects which simplify significantly the analysis. First of all, the field mode equations (20) have an analytical solution in terms of the well-known Whitakker functions, so we do not have to solve the equation numerically. Second, at time $t\rightarrow-\infty$ we have $s,\dot{s}\cdots\rightarrow 0$, so that the system is adiabatic initially, and there is no ambiguity when imposing initial conditions to the field modes. And third, as we shall see, the system behaves in such a way that, as long as we are well before the instability, $\langle\bar{\psi}\psi\rangle_{\rm ren}$ can be approximated by the fourth order in its adiabatic expansion, giving a final renormalized bilinear that can be easily integrated.

It is useful to define a new dimensionless time $z\equiv mt$ and momenta $\kappa\equiv k/m$. The field equations (20) for $h_{k}^{I}$ and $h_{k}^{II}$ in terms of these variables become

[TABLE]

and from these, we obtain the second-order uncoupled equations

[TABLE]

Let us also define a dimensionless frequency $\omega_{\kappa}\equiv\sqrt{\kappa^{2}+1}$, so that $\omega=\sqrt{k^{2}+m^{2}}=m\omega_{\kappa}$. The general solution for $h^{I}_{k}(t)$ is a linear combination of the first and second kind Whittaker functions $M_{\alpha,\lambda_{1}}\left(2i\omega_{\kappa}t\right)$ and $W_{\alpha,\lambda_{1}}\left(2i\omega_{\kappa}t\right)$, where $\alpha\equiv\frac{-i\mu}{\sqrt{\kappa^{2}+1}}$ and $\lambda_{1}\equiv-\frac{1}{2}i(2\mu-i)$. The solution for $h^{II}_{k}(t)$ is similar, with the change $\lambda_{1}\to\lambda_{2}\equiv-\frac{1}{2}i(2\mu+i)$, so we have

[TABLE]

Note that $h_{k}^{I}$ and $h_{k}^{II}$ must obey the constraint (27), so there is only 1 degree of freedom in the fermion solution, which is determined when imposing the initial conditions. To fix the constants in the linear combinations, we impose the adiabatic behavior (29) at $z\to-\infty$, getting

[TABLE]

The final solution is then

[TABLE]

The renormalized expectation value $\langle\bar{\psi}\psi\rangle_{ren}$ is given, from (70), by

[TABLE]

where we have that the adiabatic contributions of order $n$ go as $\langle\bar{\psi}\psi\rangle^{(n)}\propto c_{n}(\mu,\kappa)z^{-n}$, with $c_{n}$ time-independent functions of $\mu$ and $\kappa$. The above integral is finite, as one can easily check from the asymptotic expansion of the Whittaker function $W_{\alpha,\lambda}(x)$.

We can compute analytically the leading term at $z\to-\infty$ by performing the adiabatic expansion of $\langle\bar{\psi}\psi\rangle$ up to fourth order, and subtracting from it the zeroth, first, second, and third orders. Therefore, the leading behavior at very early times is

[TABLE]

Computing the integral, we finally get

[TABLE]

Substituting (131) in this expression, and setting $a=1$, we finally obtain

[TABLE]

where we have written the solution in terms of $z$. In Fig. 1 we show $\frac{1}{m^{3}}|\langle\bar{\psi}\psi\rangle|_{ren}$ as a function of time, comparing the exact result (137) with the approximation (140). At very early times $z\to-\infty$ we have, as expected, $\langle\bar{\psi}\psi\rangle_{ren}\sim 0$. We observe that the approximation holds quite well, except when the instability is approached.

Appendix C Adiabatic expansion

In this appendix, we provide the terms of the adiabatic expansion of the spin-1/2 field modes up to fourth order. Although the first- and second-order terms have already been written in Section III, we copy them here for convenience. As introduced in Eqs. (30) and (31), the adiabatic expansion takes the form

[TABLE]

The terms $G^{(n)}$ can be obtained from $F^{(n)}$ with the relation $G^{(n)}(m,s)=F^{(n)}(-m,-s)$, so we do not explicitly write them here. We denote by $f_{x}^{(n)}$ and $f_{y}^{(n)}$ to the real and imaginary parts of $F^{(n)}$ respectively, so that $F^{(n)}=f_{x}^{(n)}+if_{y}^{(n)}$.

The first-order terms are

[TABLE]

The second-order terms are

[TABLE]

The third-order terms are

[TABLE]

Finally, the fourth-order terms are

[TABLE]