Influence of backreaction of electric fields and Schwinger effect on inflationary magnetogenesis

We study the generation of electromagnetic fields during inflation when the conformal invariance of Maxwell's action is broken by the kinetic coupling $f^{2}(\phi)F_{\mu\nu}F^{\mu\nu}$ of the electromagnetic field to the inflaton field $\phi$. We consider the case where the coupling function $f(\phi)$ decreases in time during inflation and, as a result, the electric component of the energy density dominates over the magnetic one. The system of equations which governs the joint evolution of the scale factor, inflaton field, and electric energy density is derived. The backreaction occurs when the electric energy density becomes as large as the product of the slow-roll parameter $\epsilon$ and inflaton energy density, $\rho_{E}\sim \epsilon \rho_{\rm inf}$. It affects the inflaton field evolution and leads to the scale-invariant electric power spectrum and the magnetic one which is blue with the spectral index $n_{B}=2$ for any decreasing coupling function. This gives an upper limit on the present-day value of observed magnetic fields below $10^{-22}\,{\rm G}$. It is worth emphasizing that since the effective electric charge of particles $e_{\rm eff}=e/f$ is suppressed by the coupling function, the Schwinger effect becomes important only at the late stages of inflation when the inflaton field is close to the minimum of its potential. The Schwinger effect abruptly decreases the value of the electric field, helping to finish the inflation stage and enter the stage of preheating. It effectively produces the charged particles, implementing the Schwinger reheating scenario even before the fast oscillations of the inflaton. The numerical analysis is carried out in the Starobinsky model of inflation for the powerlike $f\propto a^{\alpha}$ and Ratra-type $f=\exp(\beta\phi/M_{p})$ coupling functions.

Observed intergalactic magnetic fields can have either astrophysical or primordial origin, and both magnetogenesis scenarios are currently under discussion. Although astrophysical mechanisms based on a Biermann battery [19][20][21] have been proposed to generate the "seed" magnetic fields and different types of dynamo can enhance them [22][23][24][25], it is problematic to embed the magnetic fields with a large correlation length into the cosmic voids. Therefore, the primordial origin of the large scale magnetic fields seems to be more realistic. In particular, the phase transitions in the early Universe may lead to the magnetic fields of the necessary strength [26][27][28][29][30][31]. However, their coherence length is determined by the horizon size at the moment of phase transition and is much less than Mpc today. Then the most natural is the inflationary magnetogenesis, proposed in Refs. [32,33], which can easily attain very large coherence length.
Since Maxwell's action is conformally invariant, the fluctuations of the electromagnetic field do not undergo enhancement in the conformally flat inflationary background [34]. In order to generate electromagnetic fields we need to break the conformal invariance. This can be done by introducing the interaction with scalar or pseudoscalar inflaton fields or curvature scalar (see, e.g., the pioneer works [32,33,35,36]). In our study we consider the kinetic coupling of the electromagnetic field to the scalar inflaton via the term f 2 (φ)F µν F µν , which was first introduced by Ratra [33] and has been revisited many times for different types of coupling functions [37][38][39][40][41][42][43][44][45].
This model modifies the standard electromagnetic Lagrangian, multiplying it by the function of the inflaton field. As it was mentioned in [40], if one rescales the electromagnetic potential in order to absorb f 2 , electric charges of particles effectively will depend on f −1 . Therefore, to avoid the strong coupling problem during inflation one needs to require f ≥ 1. Since the inflaton field and the scale factor change monotonously during inflation it is natural to assume that the coupling function is a decreasing function during inflation which attains large values in the beginning.
For decreasing coupling functions, it is well known that the electric energy density dominates the magnetic one [39,40,44] and if we try to generate the magnetic field strong enough to be in accord with the observations, the electric field appears to be even stronger and its energy density exceeds that of the inflaton field. This is known as the backreaction problem. In previous studies, the authors tried to avoid this problem, because it does not allow us to solve the background equations and Maxwell equation separately [39,40,42,43]. However, it is interesting what happens when the backreaction becomes important and whether the amplification really stops in this regime. These are open questions in the literature.
The cosmological Schwinger effect relates the particle production by an electric field and the exponential expansion of the Universe and contains interesting features which are absent in its flat-space counterpart, namely (i) the infrared hyperconductivity in the bosonic case when the conductivity becomes very large in the limit of small mass of charged particles and (ii) the negative conductivity in the weak-field regime eE ≪ H 2 which can, in principle, lead to the enhancement of the electric field.
The induced current of created particles obtained by direct averaging of the corresponding current operator contains ultraviolet divergences, which can be regularized using adiabatic subtraction [48,52] or the point-splitting method [53]. Although these techniques remove the divergent parts, the finite part is not uniquely defined and depends on the applied subtraction scheme. Adiabatic subtraction to the second order in the adiabaticity parameter and pointsplitting regularization define the minimal subtraction scheme which is commonly used. However, in the massive limit, when the particle's mass is much greater than the Hubble parameter, the current and conductivity contain exponentially unsuppressed terms, which is not consistent with the standard Bogolyubov calculations. In order to eliminate this discrepancy the authors of Ref. [57] proposed a new maximal subtraction scheme, which normalizes the behavior of the current in the massive case. However, the main features in the small mass regime remain the same and the infrared hyperconductivity in the scalar case is observed [48,53].
In addition, the charged particles show negative conductivity in the weak-field regime, which can lead to an instability and an avalanchelike enhancement of electric field up to a certain critical value [52,58]. In the very recent article [59] Stahl considers a possibility to enhance the quantum fluctuations of the electromagnetic field even without interaction with an inflaton. However, the negative conductivity is rather speculative [57,58] and may be an artifact of used subtraction schemes.
An attempt to combine the generation of the electromagnetic field due to kinetic coupling with the inflaton, the Schwinger effect and the backreaction into a consistent picture was made in Ref. [62]. It was shown that the expressions for the Schwinger current in the time-dependent electric background in the strong-field regime have the same functional dependence as in the case of a constant electric field. However, the backreaction on the background evolution was calculated perturbatively only in the first order, which is valid only at the early stages of inflation. The impact on the generated magnetic field was not discussed as well.
This paper is organized as follows. We derive a self-consistent system of equations which describes the joint evolution of the scale factor, inflaton field, and electric field energy density in Sec. II, where we take into account the backreaction of generated fields on the background evolution and the Schwinger effect which is important at the late stages of inflation. In Sec. III, we study the regime when the backreaction becomes important and analyze the main features of the electromagnetic field power spectra generated in this regime. The results of numerical calculations of the power spectra of generated fields and the present value of the magnetic field in the Starobinsky inflation model for two types of coupling functions are represented in Sec. IV. The summary of the obtained results is given in Sec. V.

II. SELF-CONSISTENT SYSTEM OF EQUATIONS
We consider a spatially flat Friedmann-Lamaître-Robertson-Walker Universe with metric tensor and use the natural system of units where = c = 1, M p = (8πG) −1/2 = 2.4 · 10 18 GeV is a reduced Planck mass, and e = √ 4πα ≈ 0.3 is the absolute value of the electron's charge. The action which describes inflaton field φ, electromagnetic field A µ , and charged field χ (either bosonic or fermionic) reads where V (φ) is the inflaton effective potential; f (φ) is the kinetic coupling function; and L charged (A, χ) is a gaugeinvariant Lagrangian of the charged field χ interacting with the electromagnetic field A µ . The corresponding Euler-Lagrange equations have the form where the 4-current is defined as usual The right-hand side of Eq. (3) describes the backreaction of created electric fields on the evolution of the inflaton field. Evolution of the Universe is driven by the total energy density of all fields. It can be calculated as the 00-component of the stress-energy tensor. The latter is defined as usual as where L 0 is the Lagrangian density of the inflaton and electromagnetic fields which is represented by the first three terms in the square brackets in Eq. (2). In the simplest case the inflaton field is spatially homogeneous. If the coupling function f (φ) always decreases in time, then it is well known [39,40] that the electric component of the created electromagnetic field dominates the magnetic one and leads to the backreaction problem. Therefore, we take into account the presence of electric field F 0i = a(t)E i (t) and neglect the magnetic component F ij = a 2 ε ijk B k ≈ 0. Then the energy density reads where ρ χ is the energy density of the charged particles produced by the Schwinger effect. The Friedmann, Klein-Gordon, and Maxwell equations take the following form: where overdots denote the derivatives with respect to cosmic time and the prime denotes the derivative with respect to the inflaton field. It has to be mentioned that we consider the classical electric field, which is generated from quantum fluctuations due to interaction with the inflaton. It is useful to analyze the mode composition of this field. When the modes are inside the horizon they oscillate in time and have to be treated as quantum fluctuations. However, when they cross the horizon the mode functions start to behave monotonically and can be chosen real. According to Ref. [63], this corresponds to the quantum-to-classical transition and these modes can be treated classically. Therefore, the contribution to the electric field is made only by the modes, which are outside the horizon. Since their wavelength is larger than the horizon, i.e. the largest observable scale, the corresponding electric field can be considered homogeneous.

A. Equation for electric field energy density
It is convenient to rewrite the Maxwell equation (10) in terms of electric field energy density ρ E = f 2 E 2 /2: The "classical" part of the electric energy density at a certain moment of time is determined by the modes, which crossed the horizon from the beginning of inflation till the moment under consideration: where k i is the momentum of the mode which crosses the horizon at the beginning of inflation. However, Eq. (11) does not take into account the fact that the number of relevant modes with wavelength larger than the horizon changes in time. In order to deal with this, we should introduce an additional term which describes the modes crossing the horizon at a given time t and starting to contribute to the total energy density of the electric field: The power spectra of electric and magnetic fields are expressed through the mode function of the electromagnetic field A(k, t) in the standard way [39]: This mode function satisfies the following equation (in the conformal time η): and the initial conditions for the modes inside the horizon have the form of the Bunch-Davies vacuum When a certain mode crosses the horizon it changes its dependence on time from oscillatory to monotonous. We can assume that at the moment of horizon crossing its behavior is still approximately described by Eq. (17). Then and the "boundary" term (13) takes the form Finally, the equation which governs the behavior of electric field energy density is given bẏ It is possible to derive Eq. (20) in an alternative way by inspecting the time dependence of the electric power spectrum. In the deeply subhorizon regime k ≫ aH the mode function takes the value of the Bunch-Davies initial conditions (17). Far outside the horizon k ≪ aH, the solution of Eq. (16) is Matching solutions (17) and (21) at the moment of horizon crossing t k when k = aH, we find where all the quantities with subscript k must be taken at the moment of horizon crossing t k . Then, using Eqs. (14)-(15), we can write the power spectra of generated fields as follows: Integrating Eq. (23) over modes outside the horizon at a given moment of time, we obtain the energy density of the electric field: where S(k) is a function which only depends on momentum k and is time independent. Then, obviously, the time derivative of the energy density equalsρ Taking into account that t k (k = k H ) ≡ t, we obtain Eq. (20) without the Schwinger term, which was also absent in Eq. (16) and could be added in the final equation phenomenologically. However, it is impossible to derive a similar equation for the magnetic energy density, because the terms proportional to C 1 and C 2 in Eq. (24), in general, have different behaviors in time and their relative contribution could change during the time evolution. In other words, the time dependence of the magnetic power spectrum cannot be extracted in the form of a universal function like in the case of electric field energy density (23).
At the moment of horizon crossing there is only the term in Eq. (24), proportional to C 1 because the integral vanishes. If af 2 is an increasing function, then the second term is a decaying mode and it can be neglected. In this case it is possible to obtain the equation for the magnetic energy density. However it is not interesting, because the coupling function decreases slower than a −1/2 and the generated electromagnetic fields are too small to explain the observational data [39]. We do not take into account the possibility of a growing coupling function because it causes the strong coupling problem in the beginning of inflation.
To derive the boundary term in Eq. (20) we used the decomposition of the electromagnetic field operator over the set of Fourier modes. However, in the presence of the Schwinger current, strictly speaking, this is not advantageous because this current has a nonlinear dependence on the electric field. Nevertheless, we will show in the following sections that the Schwinger effect is important only at the late stages of inflation when the electric field energy density is large and the boundary term is negligible in comparison with other terms in Eq. (20). On the contrary, at the beginning of inflation when the electric energy density is close to zero, the boundary term is very important because it generates the initial value of the electric field which will be enhanced due to the interaction with the inflaton. At that time the Schwinger effect is negligibly small and can be excluded from the consideration. This justifies our derivation.

B. Schwinger effect
To close the system of equations we need two other equations which determine the Schwinger current j and the energy density of particles ρ χ created via the Schwinger process.
We will use the expressions for the Schwinger currents derived in the minimal subtraction scheme in Refs. [48,52]. The general expressions are rather cumbersome. Therefore, it is more convenient to use their asymptotics governed by the parameters where m is a charged particle's mass, and the effective charge e eff = e/f naturally appears in the second expression.
Since the electric energy density always satisfies ρ E ρ tot = 3H 2 M 2 p when the backreaction is taken into account, the parameter L is damped by the effective charge during almost the entire period of inflation and can become large only at the end, when the inflaton field quickly rolls down to the minimum of its potential and f → 1. The typical value of the Hubble parameter, which can be fixed from the observations of the amplitude of primordial perturbations [64], is H ∼ 10 −5 M p ∼ 10 13 GeV. On the other hand, the lightest charged spin-1/2 particle, the electron, has a mass m ∼ 10 −3 GeV and the lightest charged scalar particle, the pion, has a mass m ∼ 0.1 GeV. Therefore, we are interested in the small mass limit M ≪ 1.
In the weak-field regime, L ≪ 1, the fermionic Schwinger current has the form [52] where γ E ≈ −0.577 is the Euler's constant. It is important to note that the current is always negative in the low mass limit m < H and, in principle, can lead to the enhancement of the electric field [59].
In the bosonic case the current is positive; however, it has a very large value in the low mass regime as it diverges like 1/M 2 at M → 0 [48,53]: In the literature this phenomenon is often called the infrared hyperconductivity. It has to be mentioned that this expression is valid in the extremely low field regime, when L ≪ M ≪ 1. When L becomes larger than M , the current quickly diminishes to the values comparable with the fermionic case (28) and in the region L ∼ 1 the negative conductivity can be observed [53]. Expressing the dissipation term on the right-hand side of Eq. (20) in terms of ρ E , we obtain It has to be compared with the main terms in Eq. (20) which are of order Hρ E . For scalar particles, the weak-field expression (30) becomes important when The corresponding L-parameter must be smaller than M , which gives This implies that at the moment when expression (30) may be important the electric energy density ρ E ∼ H 4 ≪ ρ tot . This can take place only at the early stages of inflation when the coupling function is very large. Therefore, the two inequalities (32) and (33) can be satisfied simultaneously either for the extremely low masses or in the case where the coupling function changes not very strongly during inflation. Both cases are unfavorable from the point of view of magnetogenesis as the generated magnetic fields would be very weak: in the first case the magnetogenesis terminates too early, and in the second case the enhancement can be insufficiently strong to explain the observational data. In any case, before going to the strong-field regime L ≫ 1 one has to check that at least one of the inequalities is violated. This can be done only numerically for each choice of the coupling function. For spin-1/2 particles, as well as for the bosons in the region M ≪ L 1, the weak-field expression (31) could play a role only if e/(πf ) 1. However, this is never achieved during inflation, because e/π < 1 and f (φ) > 1. Although the negative conductivity seems to be advantageous for the generation of an electric field, in the kinetic coupling model it has no effect. Therefore, for fermions, one has to consider only the strong-field regime.
The strong-field expressions for the scalar and spinor Schwinger currents are very similar and have the form where g b = 1 and g f = 2 are the numbers of spin degrees of freedom. Expressing the Schwinger term in Eq. (20) in terms of energy density, we obtain This term becomes comparable with Hρ E only at the late stages of inflation when In this case, the parameter L is really large and expression (35) is applicable. Thus, the strong-field regime could become efficient only in the late stages of inflation, when the inflaton field [and the coupling function f (φ)] changes very quickly. This results in an abrupt decrease of the electric energy density due to a high conductivity of produced plasma. Therefore, the main role of the Schwinger effect is to stop the generation of electromagnetic fields. After the effect turns on the magnetic field evolves according to the magnetic flux conservation law, i.e. B ∝ a −2 . Therefore, the magnetic field power spectrum has to be determined before the Schwinger effect turns on.
It is important to mention that all the expressions for the Schwinger current were obtained in the case of a constant electric field energy density. However, such a behavior cannot be realized in an expanding Universe without some artificial currents which in addition could violate the second law of thermodynamics [61]. Moreover, under the real circumstances, the electric field is time dependent and often grows very quickly. Nevertheless, the induced current due to Schwinger pair production was evaluated in Ref. [62] in the time-dependent electric background generated due to the kinetic coupling with the inflaton. The author considers the powerlike coupling function, f ∝ a α , and pays particular attention to the case α = 2. The strong-field expression (4.21) in Ref. [62] for the induced current in the case α = 2 can be rewritten in our notations as follows: which differs from our Eq. (34) only by a factor of 3/7 which is of the order unity. Therefore, our expression has the correct functional dependence and can be used for the numerical analysis of the Schwinger effect. Finally, we also need an equation governing the evolution of the created particles. It is natural to require that the energy dissipated by electromagnetic fields be transfered into created particles. Then, the energy density of produced particles can be described phenomenologically by the following equation: where the source term is given by Eq. (35) [or by Eq. (30) if the weak-field regime for the scalar particles appears to be effective]. Since we consider the particles' masses smaller than the Hubble parameter, we treat the created particles as ultrarelativistic, which leads to the factor 4H in the above equation.

III. POWER SPECTRA IN THE BACKREACTION REGIME
When the electric field energy density becomes large enough, the backreaction can change the regime of evolution of the inflaton field. Let us estimate the corresponding density when the backreaction starts to influence the background evolution. For this purpose, we consider Eq. (9), in which the right-hand side can be rewritten as follows: We suppose that for a short time near the moment of the backreaction occurrence the coupling function behaves like a power of the scale factor f ∝ a −γ , γ > 0 (we will see that just after entering the new regime of evolution γ ≃ 2 independently of the explicit expression of the coupling function). Then we have Since inflation proceeds most of the time in the slow-roll regime, we can neglect the term with the second derivative in Eq. (39). Then we obtain the following relation: It is obvious that the new regime appears when the second term in the brackets becomes of order unity, because in this case the time derivativeφ has to decrease in order to preserve the product (here we take into account that V ′ φ as a function of φ but notφ feels the changes later, when the field significantly deviates from its original trajectory). In other words, the backreaction becomes important when Therefore, the backreaction occurs even earlier than the electric energy density becomes comparable with the inflaton energy density. As a result, the Friedmann equation (8) does not feel directly the energy density of the generated field but only through the inflaton. We would also like to mention that if the coupling function depended only on the curvature or scale factor but not on the inflaton field directly, there would not be the right-hand side in Eq. (39), and the only place where the electric energy density could cause the backreaction would be the Friedmann equation (8). In this case, the backreaction would occur when the electric energy density became comparable with that of the inflaton, i.e. ρ E ∼ ρ inf . Firstly, the scale factor evolution would change in accordance with the Friedmann equation and after that the scalar field would feel it.
We would like to emphasize that Eq. (42) was derived by using the two approximations: (i) the slow-roll regime (we neglectedφ) and (ii) that the electric energy density only begins to influence the evolution [the additional term in the parentheses in Eq. (41) is of the order unity]. Therefore, Eq. (42) can be used only to estimate the value of the electric energy density when the backreaction becomes important. During further evolution (especially at the end of inflation) the value of ρ E can be determined only numerically and is shown in Figs. 1(a) and 5(a) for the two particular choices of the coupling function. Now, let us look at Eq. (20) without the Schwinger term. When the backreaction regime is established, the electric energy density becomes almost constant and its time derivative is small in comparison with the other terms, ρ E ≪ Hρ E . The boundary term on the right-hand side is important only at the early stages when the energy density is close to zero. When Consequently, ρ E can be neglected. Then the two remaining terms must cancel out their leading contributions that give the condition which immediately implies independently of the initial time dependence of the coupling function.
This fact allows us to estimate the behavior of electric and magnetic power spectra for modes which cross the horizon after the backreaction regime becomes relevant. For such modes, Substituting these expressions into Eqs. (23)- (24), we obtain Since the backreaction regime once established lasts until the end of inflation (or until the Schwinger effect turns on), the last mode which crosses the horizon just before the end has no time to be enhanced. Therefore, we have the following expressions for the power spectra (up to a factor of order unity): The present-day value of the observed magnetic field is determined by all relevant modes which can survive the further evolution of the Universe. Assuming the flux conservation, we have where k diff is the momentum which now corresponds to the cosmic diffusion scale, i.e., the smallest size of magnetic configuration which can survive the cosmic diffusion in the late stages of the Universe's evolution. It could be estimated If the backreaction occurs well before k diff crosses the horizon, then the magnetic power spectrum will have a blue tilt and the main contribution to the magnetic field comes from the upper integration bound The value of a0 ae can be found from the fact that the pivot scale k * crosses the horizon N * e-folds before the end of inflation: Then the present-day strength of the magnetic field is given by the following expression: This equation sets the upper bound on the present value of the generated magnetic field in models where the backreaction regime is established well before the shortest relevant mode crosses the horizon. However, this estimate is incorrect in the case where the backreaction occurs after k diff crosses the horizon. This situation will be considered numerically in the next section.

IV. NUMERICAL RESULTS
In this section we analyze the process of electromagnetic field generation numerically. In order to be specific, we consider the Starobinsky model of inflation, which is favored by the latest results of the Planck Collaboration [64]. The inflaton potential in this model has the form where the parameter µ can be fixed by the requirement that the amplitude of primordial scalar perturbations at the moment when the pivot mode k * crosses the horizon agrees with the CMB observations [64] and equals From Eqs. (8) and (9) in the slow-roll regime, taking Substituting it into Eq. (53), we find the value of the parameter µ, where φ * is the value of the inflaton field when the pivot scale k * crosses the horizon. Now, for fixed k * and P R , the only parameter of the model is φ * . It can be found numerically from the requirement that the horizon crossing of the mode k * happens at N * e-folds before the end of inflation.
In order to solve the system of Eqs. (8), (9) and (20) numerically, we impose the following initial conditions: where φ 0 has to be larger than φ * because we wish to consider all physically relevant modes which, of course, include the pivot scale k * .
The system of Eqs. (8), (9), and (20) determines how electric energy density evolves in time and does not allow one to extract the time evolution of the magnetic field which is the main purpose of our analysis. However, in the case ρ B ≪ ρ E , it is possible to calculate the magnetic power spectrum perturbatively. Using the solution of Eqs. (8), (9) and (20) as a background, we numerically integrate Eq. (16) with the initial conditions (17) for all modes which cross the horizon during the inflation stage (until the moment when the Schwinger effect turns on). Then, using the mode functions, we calculate the magnetic field power spectrum (15) and compare it with the background electric energy density. If ρ B ≪ ρ E , our approach is self-consistent. After that the present-day value of the magnetic field can be calculated using Eq. (48).
We determine the time evolution of the magnetic energy density numerically by using Eq. (16), where the Schwinger effect is taken into account only indirectly through the evolution of the inflaton field. The latter is influenced by the electric field density whose dynamics is directly affected by the Schwinger term. A more self-consistent approach would require the inclusion of the Schwinger source term directly on the right-hand side of Eq. (16). However, including such a term is a nontrivial problem, because the Schwinger current depends on the total electric energy density and, therefore, contains contributions from all modes which undergo enhancement. Due to this fact all modes are coupled and there is no possibility to solve Eq. (16) for each separate mode.
There are some general physical arguments which can help to determine the behavior of the magnetic energy density. The Schwinger effect produces charged particles which form plasma. Due to the high conductivity of the latter the electric field significantly decreases and the magnetic field becomes "frozen-in" and starts to evolve in accordance with the flux conservation law, namely B ∝ a −2 . We use this fact to calculate the present-day value of the magnetic field (48).
We consider the two types of coupling functions which are the most popular in the literature, the function which scales like f ∝ a α and the Ratra coupling function, f = exp(βφ/M p ). We would like to mention that for both coupling functions and relevant values of the parameters α and β the infrared hyperconductivity in the bosonic Schwinger current plays no role for the masses m > m π ∼ 0.1 GeV; therefore, in our analysis we use the strong-field expressions for the Schwinger current (35).
A. Coupling function f ∝ a α Here we choose the coupling function which behaves as f ∝ a α during the slow-roll inflation. By using Eqs. (8) and (9), we can express it in the slow-roll regime in terms of the inflaton field as follows: Since the coupling function primarily depends on the inflaton field φ rather than on the scale factor and it is difficult to think of the mechanism which ensures scaling ∝ a α , we will consider expression (57) as a definition that gives the scaling f ∝ a α in the slow-roll regime. It is instructive to verify Eq. (20) analytically for the above coupling function in de Sitter space-time with constant H. Martin and Yokoyama showed [39] that the electric field contribution dominates over the magnetic field contribution in this case and its power spectrum is given by Then the total energy density of the superhorizon modes equals Substituting (59) into the left-hand side of Eq. (20), we obtaiṅ In the absence of charged fields, the right-hand side for f ∼ a α reads We see that they have similar parametric behavior and differ only by a factor of order unity. Therefore, the approximate expression (19) is acceptable. Further, we will perform a numerical analysis of the system. We set α = −2.5 (so that the backreaction problem really occurs) and show in Fig. 1 the time dependencies of the electric field energy density and the inflaton field. At first, we consider only the situation when the backreaction from electric field is present in Eqs. (8) and (9) and the Schwinger effect is absent (see the red dashed lines in Fig. 1). In other words, this corresponds to the situation when there are no charged fields in the Universe.
It is obvious that at low energy density the evolution is the same as for the case without backreaction. Since the coupling function for α < −2 decreases faster than a 2 all the time, the third term on the left-hand side of Eq. (20), which is negative, is always stronger than the second one, which is positive. Hence their interplay leads to the enhancement of the electric field energy density like ∝ a 2|α|−4 ; see Eq. (59). Each separate mode starts to grow just after it crossed the horizon and the earlier it crossed, the larger amplitude it will have. This leads to a red-tilted power spectrum in the regime without backreaction, reported in Ref. [39] and which can be seen in Fig. 3. However, a different behavior arises when the electric field becomes strong enough to affect the background evolution equations. In this case it slows down the rolling of the scalar field. Also, as it was shown in Sec. III, the coupling function evolves approximately like ∝ a −2 and the two terms in Eq. (20) almost compensate each other. In this regime, the electric field energy density becomes almost constant until the end of inflation. The backreaction prolongs also inflation. Finally, we take into account the Schwinger effect. In numerical calculations we use the strong-field expression for the Schwinger source term given by Eq. (35). Since it is proportional to the cube of effective charge e/f , the Schwinger mechanism is inefficient during the early stages when f is very large. However, it becomes relevant at the late stages and then strongly diminishes the electric field density [see the blue solid lines in Figs. 1(a) and 2].
The energy density of the charged particles created due to the Schwinger process also becomes non-negligible only in the late stages. We plot its time dependence in Fig. 2 by the green dashed-dotted line. After the Schwinger effect turns on it quickly transfers almost all the energy density from the electric field to the charged particles, starting the process of reheating even before the fast oscillations of the inflaton field (the so-called "Schwinger reheating" mentioned in Ref. [55]). It is important to note that the expression for the Schwinger current (35) was derived in de Sitter space-time; therefore, its application during the stage of preheating is not well justified.
Finally, we carry out the perturbative calculation of the magnetic energy density. Using the background solutions for the inflaton and scale factor we numerically solve the mode equation (16) for all modes which crossed the horizon during inflation. Calculating the magnetic energy density at a given moment of time we include only the modes which are outside the horizon at that moment. The corresponding dependence is shown in Fig. 2 by the purple dots. These results confirm the applicability of the perturbative approach because the magnetic energy density is indeed much smaller than the electric one. The time evolution of the magnetic energy density during the preheating stage in the absence of the Schwinger effect is shown by purple dots in Fig. 2. It should be compared with the ρ B ∝ a −4 behavior plotted by the black dotted line which takes place when the Schwinger effect is present.
Then, we analyze the power spectra of generated fields and their dependence on the parameter α. Figure 3 shows the electric [panel (a)] and magnetic [panel (b)] power spectra for three different values of the power α = −2.2, −2.5, and −3 at the time when the Schwinger effect becomes efficient and inflation ends. We see that for modes which cross the horizon before the backreaction becomes relevant the power spectrum behaves like in the case where the backreaction is absent, i.e. dρ E /d ln k ∝ k 4+2α and dρ B /d ln k ∝ k 6+2α [39]. However, for shorter modes, which cross the horizon after the setting of the backreaction regime, the power spectra have similar scaling, independent of α, very close to that predicted in Eq. (47).
Having determined the magnetic power spectrum we can calculate the value of the magnetic field at the present epoch using Eq. (48). We have to include all modes which crossed the horizon after the beginning of inflation and which are longer than the cosmic diffusion scale that is estimated today as r diff /a 0 ∼ 1 A.U. [2]. The dependence of the present-day value of the magnetic field on the parameter α is shown in Fig. 4 for the two values of the number of e-folds at which the pivot scale crosses the horizon, N * = 50 (red dashed line) and N * = 60 (blue solid line). Since we know the value of the pivot scale today and we can calculate its value at the moment of horizon crossing in terms of the Hubble parameter at that moment, we can find the total number of e-folds N tot between these two moments of time. Fixing the number of e-folds from the pivot scale horizon crossing to the end of inflation N * we automatically fix  the number of e-folds from the end of inflation to today; see Eq. (50). Taking into account that the postinflationary evolution of the magnetic field is approximately determined by the flux conservation, B ∼ a −2 , and that its spectral index is n B < 4, we conclude that the smaller N * we take, the smaller the value of the magnetic field we obtain, which is obvious also from Fig. 4. Since the magnetic power spectrum is blue [see Fig. 3(b)], the main contribution to the present-day value of the magnetic field B 0 is due to the modes with k k diff . For α > −2 the backreaction never occurs and the magnetic power spectrum behaves like dρ B /d ln k ∝ k 6+2α ; see e.g. Ref. [39]. Its spectral index, n B = 6 + 2α > 2, increases with the increase of α and the amplitude of the spectrum has to decrease in order to match the unperturbed power spectrum for the modes with k > k max . As a result, the present-day value of the magnetic field is dependent on α and falls when α increases. However, for α < −2 the backreaction starts to play an important role and modifies the spectrum; see Fig. 3. In particular, for α −2.2 it changes the spectrum for the modes, which are longer than the cosmic diffusion scale, and influences the present-day value of the magnetic field. For lower values of α it becomes independent of α and its value is in accordance with the upper value (51). The smaller-than-one-order-of-magnitude discrepancy between the analytically derived limit (51) and the exact numerical value can be explained by the time dependence of the Hubble parameter close to the end of inflation and by the slight deviation of the spectral index from n B = 2. For −2.2 < α < −2, the backreaction becomes important only when the modes with k > k diff exit the horizon. As a result, for all the modes which contribute to B 0 , the spectrum is not described by Eq. (47) and has the spectral index n B = 2α + 6 < 2; i.e., their spectrum is less steep than for α < −2.2. At the same time, the amplitude of the spectrum is larger [see the green dashed-dotted line in Fig. 3(b) which is less steep and is located higher than the red dashed and blue solid lines]. This leads to larger values of B 0 for −2.2 < α < −2 compared to those for α < −2.2. This can be seen as a bump in Fig. 4.
As it was mentioned above, another feature of the generated magnetic field is its blue spectrum. This leads to very small coherence scales, comparable to 1 A.U. at the present time. Since n B ≈ 2, the values of the large-scale magnetic field on the Mpc scale are then 11 orders of magnitude smaller. In any case, the value of the magnetic field is too small to be in accordance with observational data [9][10][11][12][13].
It is important to mention the problem of the models with a red-tilted power spectrum, which, in particular, we have for α < −2. For such models it is very important when exactly the inflation starts, because all modes which cross the horizon during inflation undergo enhancement and it is larger for modes which cross the horizon earlier. For example, the power spectra in Fig. 3 and the present-day values of the magnetic field in Fig. 4 were calculated with the assumption that the inflation lasts only N e = 60 e-folds. However, if we assume that it started much earlier, this will require us to include more modes from the infrared region and will lead to earlier setting of the backreaction regime. As a result, for all modes of physical relevance (which crossed the horizon 50-60 e-folds before the end of inflation) the evolution will occur in the backreaction regime and their power spectra will be given by Eq. (47). In other words, if the inflation lasts much more than 60 e-folds, the power spectra of the generated electromagnetic fields for arbitrary α < −2 will be equivalent to the case α = −2, because the electric power spectrum will be scale-invariant and the magnetic one will have the spectral index n B ≈ 2. In this situation the principal upper bound on the present-day magnetic field (51) is valid.
Our results are in accordance with Ref. [41], where the backreaction of the electric fields was taken into account for the powerlike coupling function. We confirm that for all α < −2 in the backreaction regime the power spectra are equivalent to the "attractor" case α = −2 and this strongly suppresses the generation of magnetic fields.

B. Ratra coupling function f = exp(βφ/Mp)
Another example is given by the coupling function of the Ratra type Magnetogenesis in the Starobinsky model with the kinetic coupling (62) was studied numerically in Ref. [44]. However, not all relevant modes were taken into account and the backreaction was not correctly estimated. Here we determine the electric and magnetic power spectra in this model taking into account the backreaction and Schwinger effect.
In the slow-roll regime and in the absence of backreaction, one can derive the following approximate expressions for the time dependence of the scale factor and the inflaton field [44]: where t e = 3e √ 2/3(φ0/Mp) /(2µ) is the moment of time close to the end of inflation when the slow-roll approximation becomes inapplicable. Equation (63) makes it possible to define the moment of time t k when the mode with comoving momentum k crosses the horizon This expression is more accurate for long-wavelength modes, for which t k ≪ t e = 2 µ N e , N e = 3 4 e √ 2 3 φ 0 Mp . Then, using Eq.(64), we can calculate f k : In addition, Therefore, using Eq.(23), the electric power spectrum reads (setting f = 1 at the end of inflation) where s = 2 − √ 6 β Ne is the so-called "anomalous slope," introduced in Ref. [44]. For N e = 60, we have s = 2 − 0.04β, which is very close to the result found numerically in Ref. [44].
To estimate the magnetic power spectrum we take into account that near the end of inflation af 2 decreases very quickly and the the term proportional to C 2 in Eq. (24) dominates. The corresponding integral can be estimated as Then, which also is in agreement with Ref. [44]. However, the power spectra (68) and (70) are only applicable for relatively low momenta k ≪ a e H and for shorter modes more accurate expressions for f k , t k must be used. Moreover, when the strong backreaction regime takes place the power spectra start to behave like that defined in Eq. (47). In our numerical analysis, we use β = 15 so that the backreaction problem indeed takes place. Figure 5 shows the time dependence of the energy density of the electric field [panel (a)] and the inflaton field [panel (b)]. The green dashed lines correspond to the situation when the back reaction is taken into account and the Schwinger effect is switched off. For comparison we also plot the corresponding dependencies for the case when the electromagnetic field is a spectator field (i.e., it feels the background but does not backreact).
At the early stages of inflation the inflaton field is in the slow-roll regime and the coupling function also changes slowly; namely, it behaves like a power function [see Eq. (64)]. Under such circumstances, the Hubble term in Eq. (20) damps the evolution of the electric field energy density. It would tend to zero like a −4 if there were no boundary term on the right-hand side, which describes the new modes which cross the horizon and give their contribution to the electric energy density. Its presence leads to a dynamical equilibrium with electric energy density ρ E ∼ H 4 /(16π 2 ) being almost constant at early times.
At a certain moment of time, the negative term 2(ḟ /f )ρ E becomes larger in the absolute value than the positive Hubble term 4Hρ E and the electric energy density starts growing. As it was shown in Sec. III, when it reaches ρ E,cr ∼ ǫρ inf , the backreaction becomes very important and a new regime of evolution is established. The requirement that f ∝ a −2 in this regime significantly slows down the inflaton field and it deviates from its original trajectory. Figure 5 shows that the backreaction changes drastically the end of inflation, making the inflaton slide towards the minimum at φ = 0 without oscillations. This means the absence of the preheating stage in this case.
Taking into account the Schwinger effect resolves this problem (see the red dotted lines in Fig. 5). When the Schwinger current becomes significant it quickly diminishes the electric field and its backreaction. Then the inflaton field quickly rolls down to the potential minimum and oscillates. This helps us to restore the preheating stage, during which different particles can be created due to the parametric resonance processes [65]. In addition, the Schwinger effect by itself fills the Universe with created charged particles, implementing the Schwinger reheating scenario. This can be seen from Fig. 6, where the energy density of created particles is shown by the green dashed-dotted line.
At the same time, we should note that the expression for the Schwinger current (35) was derived in de Sitter space-time. Therefore, its application during the stage of preheating is not well justified. The magnetic field energy density can be calculated perturbatively. Figure 6 shows its time dependence in comparison with the electric energy density. The results confirm the applicability of the perturbative approach as the magnetic energy density is indeed much smaller than the electric one. The figure shows clearly that after the backreaction becomes relevant the magnetic field begins to decrease in time. Purple dots show the time evolution of the magnetic field energy density in the absence of the Schwinger effect. However, if the latter is present, then the high conductivity of produced plasma leads to the adiabatic evolution of the magnetic field for which ρ B ∝ a −4 . This is shown in Fig. 6 by the black dotted line.
Using the solutions of Eqs. (8), (9), and (20) as a background we numerically solve Eq. (16) for the electromagnetic mode function for all modes which cross the horizon during inflation and calculate the power spectra of generated electromagnetic fields. They are shown in Fig. 7 for three different values of the parameter β. As it was mentioned above, in the early stages of inflation the electric energy density is small and it does not cause the backreaction. For all the modes which cross the horizon at this stage the spectra have the same scaling with k as unperturbed power spectra described by Eqs. (68) and (70); however, their amplitudes are damped in comparison to the free case due to the backreaction which occurs later. When the coupling function starts to change faster, the power spectra undergo significant enhancement. This is reflected in the peaks in power spectra in Fig. 7 and this also corresponds to the fast increase of the electric and magnetic integral energy densities, which is shown in Fig. 6. Finally, when the backreaction regime occurs the power spectra behave as predicted in Eq. (47). For larger values of β the backreaction occurs earlier and the range of modes, which behave in such a manner, is wider.
Finally, from the known power spectra we can calculate the present-day value of the magnetic field using Eq. (48). Figure 8 shows its dependence on β for two values of the number of e-folds from the pivot scale horizon crossing to the end of inflation, N * = 50 (red dashed line) and N * = 60 (blue solid line). Note that not all modes make their contributions to the present value of the magnetic field since not all of them can survive the cosmic diffusion at the later stages of the evolution of the Universe.
For low values of the parameter β 10 the backreaction does not occur at all and the power spectra behave like in the unperturbed case considered in Ref. [44]. The magnetic field shows a steep growing dependence on β. However, for larger values of β the backreaction has an impact on the spectra and the growth slows down. The magnetic field has maximal value in the case when the peak in the magnetic power spectrum occurs for modes in the vicinity of the cosmic diffusion scale k diff . For N * = 60 this happens for β = 18 − 22. For larger values of β the backreaction regime sets much before the diffusion scale crosses the horizon and the modes which contribute to the magnetic field have power spectrum ∝ k 2 . Therefore, the dependence of B 0 tends to a constant given by Eq. (51). This can be seen as plateaus for large β in Fig. 8. Figure 8 shows that the generated magnetic field can be as large as 10 −18 G. However, its blue spectrum with n B = 4 + s for long wavelength modes and n B ≈ 2 for the shortest ones implies that the coherence length is very small and is comparable with the cosmic diffusion scale, i.e. 1 A.U. For such magnetic fields the lower bound required by the observations [9,[11][12][13] is B min ∼ 10 −18 (1 A.U./1 Mpc) −1/2 G ∼ 5 · 10 −13 G. Thus, the theoretical prediction gives much smaller value of the magnetic field than that experimentally observed.

V. CONCLUSIONS
In this work, we have studied how the backreaction of electric fields and the Schwinger effect influence the inflationary magnetogenesis in the model with the standard kinetic coupling of the electromagnetic field to the inflaton f 2 (φ)F µν F µν and decreasing coupling function f (φ). Such a model leads to the generation of strong electric fields which dominate over the magnetic ones and backreact on the inflationary dynamics. We derived the self-consistent system of equations which describes the joint evolution of the scale factor, inflaton field, and electric field during inflation. We took into account also the cosmological Schwinger effect which turned out to be important at the end of inflation, helping to finish the inflation stage and providing the mechanism of reheating which is complementary to the usual scenario with fast inflaton oscillations. To the best of our knowledge, such a self-consistent description has not been studied in the literature before.
In our analysis of the electromagnetic field evolution, it was assumed that the modes inside the horizon are described by the quantum electromagnetic field state in the Bunch-Davies vacuum. The modes far outside the horizon can be described classically, contributing to the electric and magnetic energy densities which can be associated with the observed large-scale fields. As the Universe expands, more and more modes cross the horizon and contribute. We took this fact into account through an additional boundary term in the equation for the electric energy density. The role of this term is rather subtle. Since the equation for the superhorizon modes is classical and homogeneous, it would imply only the trivial solution in the absence of this term which describes the quantum-to-classical transition.
The backreaction becomes important when the energy density of the electric field becomes as large as the product ǫρ inf of the slow-roll parameter ǫ and the inflaton energy density ρ inf . After this, the electric energy density remains almost constant until the end of inflation. The presence of such backreacting electric fields slows down the rolling of the inflaton and prolongs inflation. For the modes that cross the horizon after the backreaction becomes important, the electric power spectrum is scale-invariant and the magnetic one has a blue tilt with the spectral index n B = 2. This behavior does not depend on the explicit form of the coupling function.
Although the Schwinger effect naively should be very important in the presence of electric fields, we found that this is not the case. The point is that the Schwinger effect is tuned by the value of the effective charge e eff = e/f . At the beginning of inflation when the electric field is weak and the coupling function is very large, f ≫ 1, the Schwinger effect is not efficient and has no impact on the evolution of the inflaton and the electromagnetic fields. It becomes important at the late stage when the inflaton field is close to the potential minimum and the coupling function f → 1. In this case, the Schwinger current leads to a quick dissipation of the electric field due to the high conductivity of the produced plasma. The magnetic field becomes "frozen-in" and evolves in accordance with the flux conservation law. The energy density of the produced particles becomes comparable with the inflaton energy density. Therefore, the Schwinger effect starts to reheat the Universe even before the fast oscillations of the inflaton. Thus, we confirm that the Schwinger reheating previously discussed in Ref. [55] can help in filling the Universe with particles.
However, this does not necessarily mean the completion of reheating and does not exclude the possibility of the standard reheating mechanisms, such as parametric resonance processes (narrow, broad, and stochastic) [65] or perturbative decay of the inflaton, because even after the Schwinger effect turned on the inflaton field still oscillates in its potential minimum. To find out what mechanism is the most effective one has to extend the time evolution to the preheating stage and to include the interaction with other fields. However, the applicability of our theory is not well justified during the preheating stage. First of all, since the electric energy density has been strongly reduced due to the Schwinger process, the large field approximation for the Schwinger current given by Eq. (35) is no longer correct. The low field expressions (30)-(31) depend on the renormalization scheme and have, according to the studies in the literature, such features as the negative conductivity or the infrared hyperconductivity which are rather questionable from the physical point of view. Finally, expressions (30)- (31) were obtained in the de Sitter space-time for the constant electric field and cannot be used during the preheating stage. Therefore, an appropriate study of the Schwinger reheating requires the use of the correct expressions for the Schwinger current.
In our analysis, we used the expressions for the Schwinger current obtained in the case of a constant and homogeneous electric field in de Sitter space-time [48,52]. In Ref. [62] it was shown that in the quasiclassical approximation the current has a behavior similar to that in the case of a constant electric field, which implies that our analysis is correct at least qualitatively. Nevertheless, the most accurate approach describing the self-consistent evolution of the electric field and particles produced due to the Schwinger effect is possible to realize only in the framework of the kinetic theory. However, this lies far away from the scope of our article and should be addressed elsewhere.
The magnetic field energy density is subdominant compared to the electric one. Therefore, one can neglect it when solving the background equations. Using the inflaton and scale factor time dependencies, we integrate the equation for the electromagnetic field mode function and compute the magnetic power spectrum. Such a perturbative approach can be justified because the calculated magnetic energy density is much less than the background electric energy density. Having determined the magnetic power spectrum we estimated the value of the magnetic field at the present epoch. We took into account that after the end of inflation the magnetic field evolves according to the flux conservation law, i.e. B ∝ a −2 , and integrated the power spectrum over all modes, which can survive the cosmic diffusion during the late stages on Universe evolution. The last mode which gives the contribution has the wavelength ∼ 1 A.U. at the present time [2].
Were the Schwinger effect absent, the evolution of the magnetic field energy density could be determined perturbatively during the preheating stage as it was done during inflation. Later, during the reheating stage the Universe becomes filled with all sorts of particles and particularly charged ones. Therefore, due to high conductivity the electric field disappears and the magnetic one starts to evolve like B ∝ a −2 . The role of the Schwinger effect, therefore, is to turn on this regime earlier, near the end of inflation. As a result, the present-day value of the magnetic field may be reduced at most by a factor of e −2∆Nr , where ∆N r is the number of e-folds between the end of inflation and the beginning of reheating, i.e., the duration of preheating.
We studied numerically the generation of electromagnetic fields in the Starobinsky inflationary model for two types of coupling functions: (i) the powerlike f ∝ a α and (ii) the Ratra-type f = exp(βφ/M p ). Magnetogenesis in the kinetic coupling model with a powerlike coupling function was considered previously in Ref. [40], where it was shown that the backreaction takes place for α < −2.2 and should be included self-consistently in analysis. We showed that the backreaction becomes important when ρ E ∼ ǫρ inf ≪ ρ inf and the unperturbed calculation becomes inapplicable even earlier, for α < −2.1. Moreover, we extended the analysis for all negative values of α and showed that the present-day value of the magnetic field cannot exceed 10 −22 G, which is still too small to explain the observational data. We numerically confirm the results of Ref. [41] showing that for all α < −2 in the backreaction regime the electromagnetic power spectra are similar to the case α = −2 (the so-called attractor) and this is the reason for strong suppression of the magnetic field.
For the second type of coupling function, the rapid enhancement of the electromagnetic fields occurs only at the late stage of inflation when the inflaton field starts to deviate from the slow-roll regime. Although the spectral index for the modes which contribute to the observed value of the magnetic field is determined by the free dynamics without backreaction and is equal to n B = 6 − √ 6β/N e , the amplitude of the spectrum is damped due to the backreaction. In the most favorable situation, when the physically relevant modes cross the horizon at the time when the coupling function decreases faster than a −2 , the strength of the present-day magnetic field can reach ∼ 10 −18 G. However, the spectrum of this field has a blue tilt that leads to a small coherence length comparable with the cosmic diffusion scale, i.e. ∼ 1 A.U. For such magnetic fields the lower bound required by the observations [9,[11][12][13] is B min ∼ 10 −18 (1 A.U./1 Mpc) −1/2 G ∼ 5 · 10 −13 G. Therefore, the theoretical prediction lies well below the observational limit.
It is known that the electric fields produced in magnetogenesis models induce the curvature perturbations during inflation and they can give tight constraints on the models from the CMB observations. These constraints on the kinetic coupling models have been studied in the literature (e.g. [66,67]). Since the electric spectrum has a large power on the CMB scale in Fig. 3, this case may be excluded by this argument. However, the calculation of the curvature perturbation power spectrum in the strong backreaction regime deserves a separate investigation and lies outside the scope of our article.
In our work, we considered only the case of kinetic coupling between the inflaton and the electromagnetic field and showed that the backreaction from the electric fields strongly suppresses magnetogenesis. This makes the kinetic coupling model unfavorable in contrast to the case of axion coupling of the pseudoscalar inflaton ϕ via the term L ∝ ϕF µνF µν , whereF µν is the dual electromagnetic tensor [68]. This case seems to be promising because it generates the helical magnetic fields [35,[68][69][70], which are more stable against dissipation and their coherence length can be increased due to the inverse cascade process [71][72][73][74][75][76][77]. Nevertheless, this scenario also suffers from the backreaction problem [78] and it would be interesting to implement our equations in this case. We plan to address this issue in future studies.