The Universe after inflation: the wide resonance case

2

January 1997

PHYSICS

LETTERS 6

Physics Letters B 390 ( 1997)80-86

The Universe after inflation: the wide resonance case S.Yu. Khlebnikov a, 1.1. Tkachev b~c a Department of Physics, Purdue University, West Lafayette, IN 47907, USA b Department of Physics, The Ohio State University, Columbus, OH 43210, USA ’ Institute for Nuclear Research of the Academy of Sciences of Russia, Moscow II 7312, Russia

Received 4 September 1996; revised manuscript received 15 October 1996 Editor: H. Georgi

Abstract We numerically study the decay of massive and massless inflatons # into massive excitations X, via a d2Xz coupling, in the expanding Universe. We find that a wide enough resonance can withstand expansion of the Universe, and we determine the range of parameters in which that happens. For a massive inflaton, the effective production of particles with mass ten times that of the inflaton requires very large values of the resonance parameter q, q 2 108. For these large q, the maximal size of produced fluctuations is significantly suppressed by back reaction, but at least in the Hartree approximation they are still not negligible. For a massless inflaton with a h44/4 potential, the expansion of the Universe completely prevents resonant production of particles with masses larger than &5(O) for values of q smaller than 106. PACS:

98.8O.Cq;05.7O.Ph

The quantitative theory of reheating after inflation is an important problem in inflationary cosmology [ 1 J because the way the reheating proceeds determines whether certain symmetries are restored after reheating and which scenarios of baryogenesis are feasible. It was recently realized that stimulated (resonant) production of bosons can be most significant during reheating [ 2-51. However, despite considerable progress made [ 2,6,7], still little was known quantitatively about the important case when the coupling of the inflaton field (4) to bosonic fields (X) is relatively large (the wide resonance case). It was not clear precisely how wide the resonance must be (for a given mass of X) for it to be effective in the expanding Universe. If the coupling required by this condition turns out to be too large, non-linear effects will set in when X fluctuations have not yet developed large occupation

0370-2693/97/$17.00 PII SO370-2693(96)0

numbers; in that case, the resonance will be prevented altogether. Thus even the existence of wide resonance in the expanding Universe was not firmly established, let alone its quantitative details. The purpose of the present work was to study the wide resonance case by numerical means. The main quantity one wants to know for a given inflationary model is the maximal size of fluctuations of fields produced in inflaton decay. The size of these fluctuations determines which phase transitions took place after inflation [ 8,9], whether enough heavy leptoquarks could be produced to generate the baryon asymmetry [ 7-101, and how strong the cosmological breakdown of supersymmetry [ 11,121 was. The reason why fluctuations produced at reheating could be sufficiently large for the above effects to exist is that, as now recognized, reheating in many cases

Copyright 0 1997 Published by Elsevier Science B.V. All rights reserved. 1419-O

S.Yu.Khlebnikov.1.1. Tkachev/ PhysicsLettersB 390 (1997) 80-86

starts with a stage of parametric resonance (called preheating in Ref. [ 21)) during which the bosonic fluctuations grow extremely (exponentially) fast. In the expanding Universe, parametric resonance is effective only when it is sufficiently wide or when the system is (classically) conformally invariant. Otherwise, the redshift of momenta prevents the resonance from developing [ 131. Indeed, in a conformally non-invariant case, each mode of fluctuating field redshifts out of a resonance band in a finite time At; when the resonance band is narrow, At is too short for the modes to be significantly amplified. The conformally invariant case is realized when a massless inflaton decays into a massless fluctuation field. Production of excitations of a field more massive than the inflaton can only take place in the wide resonance regime. In that case, there can occur an interesting and somewhat counterintuitive phenomenon, when quanta of a field that strongly couples to the inflaton can be efficiently produced even if their mass is much larger than the frequency of the inflaton oscillations. The effect can be thought of as a “many--, few” process, in which many quanta of the inflaton field assemble together to produce a few heavy quanta. If the fields that can be efficiently produced by this mechanism include leptoquarks of a grand unified theory (GUT), the subsequent non-equilibrium decay of these leptoquarks could give rise to a sizable baryon asymmetry [7-lo]. In the present work, we wanted, first of all, to find the regions of parameters, for which a wide resonance develops despite the redshift of fields and momenta in the expanding Universe. We considered two models: in the first one, the inflaton potential was m*q5*/2, in the second, it was A4414 In both cases, we considered decay of the inflaton into a heavy scalar field X of mass Mx, due to coupling $g*+*X*. The width of the resonance in these models is determined by parameter 9: in Model 1 (a massive inflaton) q E g*#* (0) /4m*, where d(O) is the value of the inflaton field when it starts to,oscillate; ’ in Model 2 (the massless inflaton) , q F g74A. In Model 2, we have found no resonance for any Mx > d&b(O) with q up to q = 106. In Model 1 with ’ For a static universe, in the linear regime with respect to fluctuations, this agrees with the definition of the parameter q of the canonical Mathieu equation [ 141.

81

Mx/m = 10, which may be appropriate for scalar leptoquarks, we have found that the redshift is overcome only for very large q, q 2 IO*. Thus, a parametric resonance can withstand expansion of the Universe, though effects of the expansion restrict strongly the range of parameters in which that happens. Our second objective in this work was to determine the maximal size achieved by X fluctuations in Model 1. In general, both the expansion of the Universe and the back reaction play a role in determining the maximal size of fluctuations. A fully non-linear calculation along the lines of Ref. [ 51, which would allow us to include consistently all back reaction effects, was technically unfeasible for q >> 103. So, we took the back reaction into account approximately, using a simple Hartree procedure [ 2,4]. The Hartree approximation may overestimate the maximal size of fluctuations. The interaction term in the equation of motion for C$becomes important no later than (X2) becomes of the order of t$*(0) /q (the brackets denote averaging over space). One might expect that in the fully non-linear calculation, the resonant growth of (X2) will in any event not extend beyond this value. In the Hartree approximation, however, when this value is reached, (X2) merely renormalizes the mass of C#J[ 21, and the resonance growth continues. That further growth may be an artifact of the Hartree approximation in this model and disappear in the fully non-linear picture. We have done a number of fully non-linear calculations, with different values of q, for the conformally invariant case (see below). The set of data points that we have for (X2) at the end of resonance in this problem fits a curve even steeper than 1/q, in disagreement with the Hartree approximation. We are not willing to make any inferences for the present case, however, because the time evolution of fluctuations in expanding Universe is rather model dependent. We neglected any self-coupling of the field X. In the Hartree approximation, an estimate shows that a quartic self-coupling xX4/4 can be neglected during the entire resonance stage when 15 g”. If the maximal size of fluctuations is in fact smaller than obtained in the Hartree approximation, this self-coupling can be neglected at the resonance stage even for larger A. The resulting (X2),,,= as a function of q (the brackets denote averaging over space) is shown in Fig. 1, It is suppressed for large q but is still not negligible.

82

S.Yu. Khlebnikov, I.I. Tkachev /Physics Letters B 390 (1997) 80-84

1 Mi a2(7)

2 ---_

10-p 10-s &

lo-*

q3 10-s fi 10-0 ?5 10-7 10-B 10-o lo-‘0

i 102

Fig. I. Maximal value of the variance of fluctuations in Model 1 as a function of the resonance parameter. Dotted lines correspond to no back reaction and are fitted by a power law. Solid curves are obtained in the Hartree approximation. The dashed straight line is a power law fit to the Hartree data.

If we measure the strength of the (massive) fluctuations with regard to symmetry restoration [ 8,9] using the effective “temperature”: T’&/l2 s (X2),8x, then for q = log, the Hartree approximation gives Ten of order 10-3Mpt. This is of order of the GUT scale, so the process could be relevant for symmetry restoration even in grand unified models. One should keep in mind, though, that the size of fluctuations can be further altered (apart from the redshift) by the subsequent chaotic evolution [ 51, absent in the Hartree approximation. That chaotic evolution is rather fast because it is due to stimulated, i.e. classically enhanced, processes. For example, in applications to GUI baryogenesis, the chaotic evolution will in large part take place before the B-non-conserving decays of X particles. Therefore, some quantitative understanding of it is required before reliable estimates of the baryon asymmetry can be made; this will be a subject of a forthcoming paper. Below we describe our models and calculations in more detail. Model I: Massive inflaton. In conformal variables, and after a convenient resealing, the action for the matter fields is

,2

a2( 0)

x2

1 2w2x2

(1)

I

where r = ma(O) 7, q is the conformal time; (p = 4a(7)/&O)a(O), x = Xa(r)/ 5 = 40)x; $(O)a(O); a(T) is the scale factor of the Universe. We take r = 0 at the moment of the first extremum of 4( 7). Solving the equation of motion for the zeromomentum mode of the inflaton coupled to gravity, together with the Einstein equation for the scale factor, we have found C+(O)= 0.28Mpi. After r = 0, a matter dominated stage quickly sets in, for which one could use f&=(fi$)~+‘)~=(0.29r+1)~.

(2)

We continued, instead, to determine the evolution of the scale factor self-consistently, including the influence of produced fluctuations in the Einstein equations. As our self-consistent solution confirms, for some values of the parameters the matter dominated stage crosses over into a radiation dominated stage before resonance ends. The crossover happens when the interaction potential ~$c$~X~ overtakes the mass term of 4. In the Hartree approximation, the zero-momentum mode is the only non-vanishing mode of the field 9; we keep for it the notation cp.Its equation of motion is ..

@++_ a*(T)P + 49(x2bP =0, a2(o)

(3)

where the brackets denote an average over space. This should be solved with the initial conditions ~(0) = 1, @(0) = 0. The equation of motion for the Fourier components of the field x is jik + W;(T)xk

= 0,

(4)

where u;(7)

M* a*(r) = k2 _ z + Am2 a2(o) +4w2w

(5)

We will sometimes use the notation m2 X

z

xM2 m2

’

(6)

The Fourier components are normalized in such a way that the variance of the fluctuation field is

S.Yu. Khlebnikou, 1.1. Tkacheu / Physics Letiers B 390 (1997) 80-86

1 (X2) = (2*)3 s d3k(xk12.

lo-’ 10-e

Because fluctuations are produced mostly with relatively small values of (comoving) momenta k, the inevitable presence of a momentum cutoff in our calculation, somewhere at larger momenta, will automatically result in a properly renormalized (x2). If we wanted to embed our calculation into a fully quantumfield-theoretical picture, we would identify the parameters we use (m, 2, etc.) with the running parameters of that picture, taken at a normalization point of the order of our momentum cutoff; moderate changes in the normalization point will only cause corrections of relative order $/47r = qm2/4?r42 (0). The initial conditions for xk are randomly distributed with the probability density pkhkl

-- 24220)

mexp

&k(o)

[

83

IXk(O)

I2

10-a _ lo-’ 8 10-1 &

10-O lo-‘0 lo-” lo-la 10-13

Fig. 2. Fluctuation variance in Model 1, as a function of time, for two values of the mass of the fluctuating field; back reaction of fluctuations on the inflaton field was included in the Ha&ee approximation.

1 (8) .

According to Eq. ( 1) , m2/42 (0) is the parameter that regulates the strength of quantum fluctuations compared to classical ones in this model. The initial “velocities” are locked to their corresponding “coordinates”

Universe (the adiabatic limit), q,E( 7) would indeed be the correct resonance parameter at time r. To understand the effects of the redshift and the back reaction, it is convenient to see first what happens when back reaction is neglected. This is indeed the correct limit for sufficiently small q (a more precise criterion depends on m,) . Then, rn,ff = m and qff(7) e

=qe.

(11)

&CO) kk(O)

= t -iwk(o)

+

h(o)lxk,

(9)

where h = b/a. For a derivation of these initial conditions, see Ref. [ 51. The semiclassical description is reliable as long as g2/4?r = qm2/4& (0) << 1. The results of numerical calculations for the realistic case m = 10-6A4pi are presented in Figs. l-3. Three factors can potentially suppress parametric resonance in the expanding Universe, they are: the time dependent mass term 0; M$ in Eq. (5), which redshifts the modes away from the resonance; the time dependent mass term of the field (p in Eq. (3), which makes the oscillations of the background field aperiodic; and the back reaction term 0; q in Eq. (3). A very useful quantity to introduce is the effective q, q&f(7)

=

$-$g, eff

(10)

where 4 is the (decreasing) amplitude of the physical field r$, and mi+r is that field’s effective mass including the Hartree correction. In a slowly expanding

For parametric resonance to be effective at time r, qe&T) should not be much smaller than unity: qe&T) < 1 corresponds to a narrow resonance, which in this model cannot compete with the expansion of the Universe. Because a resonance will need some time 70 to start, the condition for it to begin at all, qeR(70) 2 1, already requires large values of q = qeR(0), By our definition, a resonance starts when (X2) rises above quantum fluctuations. E.g. in Fig. 2, ro x 6; some k-modes can start to grow earlier, however. Consider first the case of massless X, Mx = 0. We have found that even then the resonance starts only for q 2 103, see Fig. 1. With q 2 103, the parametric resonance starts at some point but terminates when qeR drops below unity. Hence, in the expanding Universe, the resonance stops and the maximal size of produced fluctuations is limited in this model2 even without back reaction. We have calculated the maximum value * As opposed to the conformally invariant cuse 12.15.51.

84

SYu. Khlebnikov.1.1.Tkachev/ PhysicsLettersI3390 (1997) 80-86

of (X2) that would be achieved without back reaction for a wide range of values of q. It is plotted in Fig. 1 by the dotted curve and for massless X has the label mx = 0. This curve can be approximated by (X2),,, 0; 9 ‘*, which is plotted in Fig. 1 by a solid straight line. Next, consider the case Mx # 0. The rapidly growing mass term of x is a resonance suppressing factor. In the case of the A44 model, we will see that it works against the resonance extremely effectively. In the present model, however, the frequency of oscillations of the field cpis growing at the same rate as the effective mass of x, compare Eqs. (3) and (5). In a sense, massive x modes can remain tuned to the resonance, in some analogy with the conformal case. Still, for a non-zero Mx, the resonance terminates before qeff drops below unity. This is why the curves with m, # 0 areshifted in Fig. 1 to the right with respect to the mx = 0 case. In the absence of expansion of the Universe and back reaction, the linearized equation for fluctuations, Eq. (4), would be a Mathieu equation. In that case, a necessary condition for the X particle production to fall into the strongest resonance band would be, for mx 2 1: mx 5 q’i4, or equivalently q 2 ( M,y/m)4. As one can see from Fig. 1, this estimate is four orders of magnitude off in the expanding Universe: e.g. with q N lo4 the resonance barely develops even for m x = 1. To take expansion of the Universe into account, we can replace q in the above estimate with qeE at some time before the end of the resonance stage, qe.(71). As evident from Eq. ( 1l), finding qeR(?l) is equivalent to finding $(~i). The quantities 71 and $( 71) in general depend on q. Although some analytical estimates of 4 (71) can be made [ 81, its reliable determination requires our numerical results. We have found that for the range of q we consider in Fig. 1, one can use generically f$‘( 7, ) N 10m4~*(0). Using qee instead of q repairs the estimate: to efficiently produce particles of mass Mx, one needs qed 71) 2 4, or q 2 104(Mx/m)4,in agreement with (M&n) Fig. 1. For Mx/m = 10, we get q 2 log, which for m= lOI3 GeV translates into 2 2 4. x 10m3.Considering the complexity of the phenomena involved, the estimate gz 2 10T4 made for this case in Ref. [ 101 is in a reasonable agreement with our result. The term “wide” parametric resonance may be to some extent misleading. The resonance peak

in the power spectrum is itself narrow, Ak << q.$a(T)/a(O). However, it can be anywhere in the instability band 0 < k,, 5 qifa(T) /a(O), and indeed moves through this band a few times as qeR changes. So, due to time dependence of qeff in expanding Universe, all k from this range become resonant in turn,3 and it is for that reason that one may consider q$a( 7)/a(O) as an effective width of the resonance. As a result of the motion of the resonance peak, all k up to k of order q’i4 become filled up in the power spectrum, although to a varying extent. The motion of the resonance peak is reflected in Fig. 2. For mx = 0,there is no resonance during time 10 5 7 5 12: the resonance peak has moved out of the strongest instability band. At 12 ,$ 7 5 13, there is the last resonant period, which ends when qeK becomes somewhat smaller then unity. In Fig. 1, notice also the sensitivity of the maximal size of fluctuations with respect to small changes in q. This sensitivity is expected, and one of the reasons for it is the dependence of the position of the resonance peak on q. Changing q by an amount of order fi puts the initial peak at a quite different place in the instability band. That changes the entire picture of time evolution and leads to a different value of (X2),,,ax.Small relative changes in q can be simulated by small numerical uncertainties, such as that of the coefficient in the relation &J(O) = 0.28Mp1, so for example Fig. 2 may not correspond to exactly the value q = 104. It anyway makes sense to average out the rapid changes in Fig. 1 and concentrate on the general pattern, As a general trend, we observe a broad maximum of (X2),,,axas a function of q, at fixed Mx. It occurs in the region of q where the redshift and the back reaction play comparable roles in terminating the resonance. Now consider the limit when the maximum value attained by fluctuations is determined by back reaction. This is the case for sufficiently large q. For example, for massless X, this is the case for q 2 104, see Fig. 1. The general trend in (X2),;u, in this region of q is given by

(12)

3It is perhapsnot obvious a priori that this motion of the resonance peak does not kill the resonancealtogether.

85

S.Yu. Khlebnikov.1.1. Tkochev/ PhysicsLettersB 390 (199?J 80-86

7

Fig. 3. Time evolution of the homogeneous inflaton field and the fluctuation variance in Model 1 in the Hartree approximation, for given values of the parameters.

which is plotted in Fig. 1 by the dashed straight line. In this region, we can also estimate (X2),= in the Hartree approximation analytically, and the estimate is in rather good agreement with Eq. ( 12). It proceeds as follows. We first notice that because the decay of C#J into X is a stimulated process, it mostly takes place in a short period close to the end of the resonance, when X fluctuations are already large. Then, we should distinguish between the amplitude of the inflaton oscillations just before the decay, & ~1) , which is determined mostly by the redshift, and its amplitude right after the decay, c$(72). For a wide resonance, &( 72) can be smaller than C$(71) by orders of magnitude. Such a drop is clearly seen in Fig. 3. Because 72 - ~1 is such a short time, the expansion of the Universe during it can be neglected, and from energy conservation we get g262(72)(x2)

N m2+2(71).

The value of r$(71) weakly depends on q. As before, we use the generic estimate c$(~1) N 10m2#( 0) for all q considered. We then obtain (X2)max N 10-3M&/JiT, in reasonable agreement with the numerical result ( 12). A 1/Jii (up to logarithms) estimate for (X2)maxappears in Refs. [ 2,8,11], although the way the estimate is obtained and the numerical coefficient are somewhat different from ours. We stress again that the l/& behavior may be an artifact of the Hartree approximation in this model. See also the discussion of the conformally invariant case below. Model 2: Massless injz’uton.In resealed variables, the action for the matter fields in this model is

-;-$gj$gx2- 2w2x2] ,

(15)

where Q and x are defined as before but 7 = J;&O)a(O)rl, 5 = G$(O)a(O)n. We again take 7 = 0 to be the first extremum of c,Q7). In this case we find & 0) = 0.35&l. There is no resonance for a massive enough X (see below). In that case, after the inflaton oscillations start, the Universe quickly becomes radiation dominated and we can use

(13) 1 ~0.517+

The second ingredient of our estimate is the observation [ 21 that the presence of (X2) changes the effective mass of the inflaton; the effective mass is given by m;tk = m2 + $(P), where (x2) is a slowly changing part of (X2). When (X2) becomes sufficiently large, m’& M 2(x2) N g2(X2), and the effective resonance parameter becomes qeR = g2$~~/4rn& N c$~/(X~). The resonance stops when q& 72) N 1; using IQ. (13), this gives4 4 Our usual condition qee 2 (Mx/rr~,~)~ is satisfied for qeE N I at sufficiently large q. due to the large rn,~. Notice also that because a large ,n,r lowers the Mx,/zn~n,~ ratio, the Hartree correction to

1.

(16)

All other formulas are obtained from the corresponding ones for Model 1 by replacing m2 with &p2(0). Accordingly, we now use m, = Mx/fi#( 0). For a massless field x the problem becomes conformally equivalent to the one in the Minkowski space-time: the expansion of the Universe decouples in the conformal coordinates. Parametric resonance is always effective. We have done a complete non-linear analysis of this problem in the range 1 5 q 5 1000, the mass of CJ actually helps the resonance for a while, as we have indeed observed in our numerical results.

86

S.Yu. Khlebnikou,

1.1. Tkachev /Physics

using the approach developed in Ref. [5]. We were able to trace the semiclassical thermalization [5] up to a slowly evolving state with a developed Kolmogorov turbulence. These results will be described in more detail elsewhere; for previous discussions of kinetic phenomena after preheating, see Refs. [ 5,163. For the set of large q for which we now have data for the full non-linear problem, the generic decrease in (X2) at the end of the resonance stage is faster than 1/q, in disagreement with the approximate l/J;i dependence obtained in the Hartree approximation. With a massive x field we have not found any parametric resonance for m, > 1 and q up to q = 106. Analytically, the condition for the resonance to exist at time r (neglecting back reaction) is q’/4 2 m,a( 7)/a(O). In the conformally invariant case (m, = 0), the fluctuations reach appreciable strength at r w 100. Using that value in the estimate for mx # 0, we would expect that for m, = 1, the resonance exists when q 2 107. We indeed observed a slow (physically insignificant) resonance at mx = 1 and q = lo*. In conclusion, we have studied the decay of massive and massless inflatons into massive excitations in the expanding Universe. We have found that a wide enough resonance can withstand expansion of the Universe and determined the range of parameters in which that happens. For a massive inflaton, the effective production of particles with mass ten times that of the inflaton requires very large values of the resonance parameter q, q 2 1O*.For these large q, the maximal size of produced fluctuations is significantly suppressed by back reaction, but at least within the Hartree approximation they are still not negligible. We do not know at present if this conclusion will persist in the fully nonlinear picture; the full account for rescattering in this model is thus an important problem for the future. For the massless inflaton with a A+4/4 potential, the Universe expansion completely prevents a resonant production of particles with masses larger than d&(O) for values of q smaller than q = 106.

Letters B 390 (1997) 80-86

We thank G. Dvali, L. Kofman, E. Kolb, A. Linde, A. Riotto, V. Rubakov, M. Shaposhnikov, and N. Tetradis for discussions and comments. S.K. thanks CERN Theory Division, where some of the discussions took place, for hospitality. The work of S.K. was supported in part by the US Department of Energy under Grant DE-FG02-91ER40681 (Task B), by the National Science Foundation under Grant PHY 95-01458, and by the Alfred P. Sloan Foundation. The work of LT. was supported by DOE Grant DE-AC02-76ER01545 at Ohio State. References [ I] A.D. Linde, Particle Physics and Inflationary Cosmology (Harwood, Chur, Switzerland, 1990); E.W. Kolb and M.S. Turner, The Early Universe (AddisonWesley, Redwood City, California, 1990). [2] L. Kofman, A. Linde and A.A. Starobinsky, Phys. Rev. Lett. 73 (1994) 3195. [3] J.H. Traschen and R.H. Brandenberger, Phys. Rev. D 42 (1990) 2491; Y. Shtanov, J. Traschen and R. Brandenberger, Phys. Rev. D 51 (1995) 5438. [4] D. Boyanovsky, H.J. de Vega, R. Holman, D.-S. Lee and A. Singh, Phys. Rev. D 51 ( 1995) 4419; D. Boyanovsky, M. D’Attanasio, H.J. de Vega, R. Holman and D.-S. Lee, Phys. Rev. D 52 ( 1995) 6805. 151 S. Khlebnikov and I. Tkachev, Phys. Rev. Lett. 77 ( 1996) 219. [6] M. Yoshimura, Prog. Theor. Phys. 94 ( 1995) 873; H. Fujisaki, K. Kumekawa, M. Yamaguchi and M. Yoshimun, Phys. Rev. D 53 ( 1996) 6805; hep-ph/95 1I38 1. [7] M. Yoshimum, hep_ph/9605246. [8] L. Kofman, A. Linde and A.A. Starobinsky, Phys. Rev. L&t. 76 (1996) 1011. [9] 1.1.Tkachev, Phys. L&t. B 376 (1996) 35. IO] E.W. Kolb, A.D. Linde and A. Riotto, hep-ph/9606260. I I l] G.W. Anderson, A. Linde and A. Riotto, hep-ph/9606416. 121 G. Dvali and A. Riotto, hep-ph/9606431. 131 J. Preskill, M.B. Wise and E Wilczek, Phys. Lett. B 120 (1983) 127; L.F. Abbott and P. Sikivie, Phys. Len. B 120 ( 1983) 133; 1.1. Tkachev, Sov. Astron. Lett. 12 ( 1986) 305; A.D. Dolgov and K. Freese, Phys. Rev. D 51 ( 1995) 2693. 141 See, e.g., M. Abramowitz and 1. Stegun, Handbook of Mathematical Functions (Dover, New York, 1972). 151 D.1. Kaiser, Phys. Rev. D 53 (1996) 1776. 161 D.T. Son, hep-ph/9604340.