text:
# 1 Introduction

## 1 Introduction

Suppose that during a de Sitter stage of expansion the coupling constant of an Abelian gauge field evolves in time. Thus the kinetic term of the gauge field can be written (in four space-time dimensions) as

$$S_{\mathrm{em}}=\frac{1}{4}d^4x\sqrt{G}f(\varphi )F_{\alpha \beta }F^{\alpha \beta },$$
(1.1)
where $`G`$ is the determinant of the space-time metric, $`\varphi `$ is a scalar field (which can depend upon space and time) and $`g(\varphi )=f(\varphi )^{1/2}`$ is the coupling 222The Heaviside electromagnetic system of units will be used throughout the investigation. The effective “electron” charge will then be given, in the present context, by $`e(\varphi )=e_1f(\varphi )^{1/2}`$.. This type of vertex is typical of scalar-tensor theories of gravity  and of the low energy string effective action . Early suggestions that the Abelian gauge coupling may change over cosmological times were originally made by Dirac  (and subsequently discussed in  and in ) mainly in the framework of ordinary electromagnetism.

The dynamics of the field $`\varphi `$ will be described by the action of a minimally coupled (massive) scalar. If the field is displaced from the minimum of its potential during inflation, there will be a phase where the field relaxes. Provided the scalar mass is much smaller than the curvature scale during inflation such a phase could be rather long. During the de Sitter stage, the specific form of the expanding background will dictate, through the equations of motion, the rate of suppression of the amplitude of $`\varphi `$.

While the field relaxes toward the minimum of its potential, energy is pumped from the homogeneous mode of $`\varphi `$ to the gauge field fluctuations. The function $`f(\varphi )`$ can be either an increasing function of $`\varphi `$ (leading to a decreasing coupling) or a decreasing function of $`\varphi `$ (leading to an increasing coupling). In both cases, depending upon the parameters of the model, the two-point correlation function of magnetic inhomogeneities increases during the inflationary stage. This implies that large scale magnetic fields can be potentially generated.

The implications of large scale magnetic fields for cosmology, astrophysics and high-energy physics have been recently reviewed in  where various useful references are reported. The interested student will find in  a morte physical introduction to the subject of large scale magnetic fields. The toy model discussed in the present paper is only an useful exercise suggesting a possible connection between the variation of gauge couplings and the production of large scale magnetic fields.

Recently  the variation in the fine structure constant has been suggested by observations of absorption lines at different frequencies. The results of these claims can be summarized by saying that $`\alpha _{\mathrm{em}}`$ was smaller in the past. Defining, indeed $`\mathrm{\Delta }\alpha _{\mathrm{em}}=\alpha _{\mathrm{today}}\alpha _{\mathrm{past}}`$, measurements suggest that

$$\frac{\mathrm{\Delta }\alpha _{\mathrm{em}}}{\alpha _{\mathrm{em}}}=(0.72\pm 0.18)10^5.$$
(1.2)
To avoid confusions, in the present scenario, gauge couplings are constant today and their variation occurs prior to big-bang nucleosynthesis (BBN).

The gauge coupling has to be (almost) constant by the time when the Universe is approximately old of one second, namely by the time of BBN. The abundances of light elements are very sensitive to any departure from the standard cosmological model. Hence, fluctuations in the baryon to photon ratio, matter–antimatter domains, anisotropies in the expansion of the four space-time dimensions, can all be successfully constrained by demanding that the abundances of the light elements are correctly reproduced. Following the same logic the variation in the gauge couplings can also be constrained from BBN .

The plan of the present contribution is then the following. In Section II the basic ideas concerning the model of evolution of the gauge coupling will be introduced. In Section III bounds coming both from the homogeneous and from the inhomogeneous evolution of $`\varphi `$ will be described. In Section IV the evolution of the magnetic inhomogeneities will be addressed along the various stages of the model with particular attention to the role of the two-point function. In Section V the large scale magnetic fields produced in the scenario will be estimated. Section VI contains some concluding remarks. The considerations reported in the present paper extend and complement the results of .

Finally, in concluding this Introduction, I wish to remind a couple of useful references of two other lecturers of the Chalonge school. In Ref.  an interesting discussion on the possible ambiguities arising in the observations of large scale fields in clusters is reported. In  another idea on the generation of magnetic fields is reported.

## 2 Basic Equations

Thanks to the high degree of isotropy and homogeneity of the observed Universe, the background geometry can be described using a (conformally flat) Friedmann-Robertson-Walker (FRW) line element

$$ds^2=G_{\mu \nu }dx^\mu dx^\nu =a^2(\eta )[d\eta ^2d\stackrel{}{x}^2],$$
(2.1)
where $`\eta `$ is the conformal time coordinate and $`G_{\mu \nu }`$ is the four-dimensional space-time metric. The cosmic time coordinate (often employed in this investigation) is related to $`\eta `$ as $`a(\eta )d\eta =dt`$.

From the anisotropies of the Cosmic Microwave Background (CMB) it is consistent to assume that the Universe underwent a period of inflationary expansion of de Sitter or quasi-de Sitter type. Therefore $`a(\eta )\eta _1/\eta `$ during a phase stopping, approximately, when the curvature scale was $`H_110^6M_\mathrm{P}`$. For $`\eta >\eta _1`$ (possibly after a transient period) the Universe gets dominated by radiation \[i.e. $`a(\eta )\eta `$\] and then, after decoupling, by dust matter \[i.e. $`a(\eta )\eta ^2`$\].

The action describing the dynamics of the (Abelian) gauge coupling in a given background geometry can be parametrized as

$$S=d^4x\sqrt{G}\left[\frac{1}{2}G^{\mu \nu }_\mu \varphi _\nu \varphi V(\varphi )\frac{1}{4g^2(\varphi )}F_{\mu \nu }F^{\mu \nu }\right].$$
(2.2)
In Eq. (2.2) $`\varphi `$ is not the inflaton but an extra (scalar) degree of freedom taking care of the evolution of the gauge coupling.

From Eq. (2.2) the equations of motion can be derived

$`{\displaystyle \frac{1}{\sqrt{G}}}_\mu \left[\sqrt{G}G^{\mu \nu }_\nu \varphi \right]+{\displaystyle \frac{V}{\varphi }}={\displaystyle \frac{1}{2g^3(\varphi )}}{\displaystyle \frac{g}{\varphi }}F_{\alpha \beta }F^{\alpha \beta },`$ (2.3)
$`{\displaystyle \frac{1}{\sqrt{G}}}_\alpha \left[{\displaystyle \frac{\sqrt{G}}{g^2(\varphi )}}F^{\alpha \beta }\right]=0.`$ (2.4)
If the gauge coupling does not change, the evolution of Abelian gauge fields is conformally invariant. Hence, using the the conformal time coordinate, the appropriately rescaled gauge field amplitudes obey a set of equations which is exactly the one they would obey Minkowski space. In order to simplify the explicit form of the equations of motion the electric and magnetic fields will be rescaled in such a way that the obtained system of equations will reproduce the usual (conformally invariant) system in the limit $`g`$ constant. The rescalings in the fields are

$`\stackrel{}{B}=a^2\stackrel{}{},\stackrel{}{E}=a^2\stackrel{}{},\stackrel{}{A}=a\stackrel{}{𝒜},`$ (2.5)
where $`\stackrel{}{}`$, $`\stackrel{}{}`$, $`\stackrel{}{𝒜}`$, are the flat-space quantities whereas $`\stackrel{}{B}`$, $`\stackrel{}{E}`$, $`\stackrel{}{A}`$, are the curved-space ones.

The explicit form of our system becomes, in the metric (2.1):

$`{\displaystyle \frac{^2\varphi }{\eta ^2}}+2{\displaystyle \frac{\varphi }{\eta }}+{\displaystyle \frac{V}{\varphi }}a^2^2\varphi ={\displaystyle \frac{1}{2g^3a^2}}{\displaystyle \frac{g}{\varphi }}\left[\stackrel{}{E}^2\stackrel{}{B}^2\right],`$ (2.6)
$`{\displaystyle \frac{\stackrel{}{B}}{\eta }}=\stackrel{}{}\times \stackrel{}{E},`$ (2.7)
$`{\displaystyle \frac{}{\eta }}\left[{\displaystyle \frac{1}{g^2(\varphi )}}\stackrel{}{E}\right]={\displaystyle \frac{1}{g^2(\varphi )}}\left[\stackrel{}{}\times \stackrel{}{B}{\displaystyle \frac{2}{g}}{\displaystyle \frac{g}{\varphi }}\stackrel{}{}\varphi \times \stackrel{}{B}\right],`$ (2.8)
$`\stackrel{}{}\stackrel{}{B}=0,`$ (2.9)
$`\stackrel{}{}\stackrel{}{E}={\displaystyle \frac{2}{g}}{\displaystyle \frac{g}{\varphi }}\stackrel{}{E}\stackrel{}{}\varphi ,`$ (2.10)
where $`\stackrel{}{v}`$ is the bulk velocity of the plasma. The spatial gradients used in Eqs. (2.6)–(2.10) are defined according to the metric (2.1). In Eq. (2.11) the quantity $`=\mathrm{ln}a/\eta `$ has also been introduced. $``$ is the Hubble factor in conformal time which is related to the Hubble factor in cosmic time as $`=H/a`$ where $`H=\mathrm{ln}a/t`$.

Once the background geometry is specified we are interested in the situation when the gauge field background is vanishing and the only fluctuations are the ones associate with the vacuum state of the Abelian gauge fields. Hence, Eqs. (2.6)–(2.10) allow to compute the evolution of $`\varphi `$ and the associated evolution of the two-point function of the gauge field fluctuations.

Suppose that $`\varphi `$ is originally displaced from the minimum of its potential. As far as the zero mode of $`\varphi `$ is concerned the system of equations can be further simplified:

$`{\displaystyle \frac{^2\varphi }{\eta ^2}}+2{\displaystyle \frac{\varphi }{\eta }}+{\displaystyle \frac{V}{\varphi }}a^2={\displaystyle \frac{1}{2g^3a^2}}{\displaystyle \frac{g}{\varphi }}\left[\stackrel{}{E}^2\stackrel{}{B}^2\right],`$ (2.11)
$`{\displaystyle \frac{\stackrel{}{B}}{\eta }}=\stackrel{}{}\times \stackrel{}{E},`$ (2.12)
$`{\displaystyle \frac{}{\eta }}\left[{\displaystyle \frac{1}{g^2(\varphi )}}\stackrel{}{E}\right]={\displaystyle \frac{1}{g^2(\varphi )}}\stackrel{}{}\times \stackrel{}{B},`$ (2.13)
$`\stackrel{}{}\stackrel{}{B}=0,\stackrel{}{}\stackrel{}{E}=0,`$ (2.14)
By now combining together the modified Maxwell’s equations we obtain the evolution of the magnetic fields

$$\stackrel{}{B}^{\prime \prime }2\frac{g^{}}{g}\stackrel{}{B}^{}^2\stackrel{}{B}=0,$$
(2.15)
where the prime denoted derivation with respect to the conformal time coordinate (the over-dot will denote, instead, derivation with respect to cosmic time).

Once the evolution of the metric is specified, Eq. (2.11) dictates a specific evolution for $`\varphi `$ and the evolution of $`\varphi `$ will determine, in its turn, the evolution of the gauge fields. The interesting initial conditions for the system are the ones where the classical gauge field background vanishes. Thus, when the homogeneous component of $`\varphi `$ starts its evolution during the de Sitter phase, quantum mechanical fluctuations will be postulated as initial conditions of gauge inhomogeneities.

When the background geometry evolves from the de Sitter phase to the subsequent epoch, massive quanta of $`\varphi `$ are produced. The amount of the produced inhomogeneous modes of $`\varphi `$ can be computed and it will be shown that the associated energy density will always be smaller than the one of the homogeneous mode. This analysis will be one of the subjects discussed in Section III.

## 3 Constraints on the evolution of the gauge coupling

Sub-millimiter tests of the Newton’s law show that no deviations are observed down to distances as small as $`0.1`$ mm . Therefore, $`\varphi `$ should be massive. Of course the potential of $`\varphi `$ may be much more complicated than the one provided by a simple mass term. However, for sake of simplicity, a massive scalar will be analyzed since, already in this case, interesting effects can be analyzed. In spite of the fact that this choice is apparently simple, various constraints on the scalar mass appear.

### 3.1 Evolution of the homogeneous mode

Suppose that the potential term driving the evolution of the gauge coupling is simply

$$V(\varphi )\frac{m^2}{2}\varphi ^2.$$
(3.1)
During the inflationary stage of expansion the scale factor evolves as $`a(\eta )=(\eta _1/\eta )`$ for $`\eta <\eta _1`$. The evolution of $`\varphi `$ is obtained by solving

$$\varphi ^{\prime \prime }+\frac{2}{\eta }\varphi ^{}+\frac{\mu ^2}{\eta ^2}\varphi =0,$$
(3.2)
where $`\mu =m/H_1`$ and where the relation $`\eta ^1aH`$ has been used. If $`\mu 1`$, for $`\eta <\eta _1`$ the solution of Eq. (3.2) can be written as

$$\varphi _\mathrm{i}(\eta )=\varphi _1\varphi _2\left(\frac{\eta }{\eta _1}\right)^3.$$
(3.3)
The end of the inflationary stage of expansion may not be directly followed by the radiation dominated phase. In the intermediate phase the scalar mass is still small than the curvature scale but the curvature decreases, in general, faster than during the inflationary phase since the background is neither of de Sitter nor of quasi-de Sitter type. The evolution of the scale factor can be parametrized as $`a(\eta )\eta ^\alpha `$ where, in order to fix the ideas, $`\alpha 2`$ could be assumed 333The case $`\alpha =2`$ corresponds to a matter-dominated intermediate stage.. The evolution of $`\varphi `$ will simply be

$$\varphi _{\mathrm{rh}}(\eta )=\varphi _1+\varphi _2\left[\frac{2\alpha 4}{2\alpha 1}+\frac{3}{2\alpha 1}\left(\frac{\eta }{\eta _1}\right)^{2\alpha 1}\right],\eta _1<\eta <\eta _\mathrm{r},$$
(3.4)
where the continuity between Eq. (3.3) and Eq. (3.4) has been required, so that $`\varphi _\mathrm{i}(\eta _1)=\varphi _{\mathrm{rh}}(\eta _1)`$ and $`\varphi _\mathrm{i}^{}(\eta _1)=\varphi _{\mathrm{rh}}^{}(\eta _1)`$.

After $`\eta _\mathrm{r}`$ the background enters a radiation dominated phase and the evolution of $`\varphi `$ can be explicitly solved in cosmic time. The equation for $`\varphi `$ , in this phase, is given by

$$\ddot{\varphi }+3H\dot{\varphi }+m^2\varphi =0,H=\frac{\dot{a}}{a},$$
(3.5)
which in terms $`\mathrm{\Phi }=a^{\frac{3}{2}}\varphi `$, becomes

$$\ddot{\mathrm{\Phi }}+\left[m^2\frac{3}{2}\dot{H}\frac{9}{4}H^2\right]\mathrm{\Phi }=0.$$
(3.6)
In the radiation-dominated stage of expansion Eq. (3.6) becomes

$$\ddot{\mathrm{\Phi }}+\left[m^2+\frac{3}{16t^2}\right]\mathrm{\Phi }=0,$$
(3.7)
whose solution can be written in terms of Bessel functions

$$\mathrm{\Phi }(mt)=\sqrt{mt}\left[AY_{\frac{1}{4}}(mt)+BJ_{\frac{1}{4}}(mt)\right].$$
(3.8)
For $`mt1`$, $`\varphi `$ has a constant mode and a solution as $`t^{1/2}`$. Recalling the relation between cosmic and conformal time and imposing the continuity of $`\varphi `$ and $`\varphi ^{}`$ (in $`\eta _r`$) with the solution of Eq. (3.4) the following form can be obtained:

$$\varphi _\mathrm{r}(\eta )=\varphi _1+\varphi _2\left[\left(\frac{2\alpha 4}{2\alpha 1}+\frac{6(1\alpha )}{2\alpha 1}\left(\frac{\eta _1}{\eta _\mathrm{r}}\right)^{2\alpha 1}\right)+3\left(\frac{\eta _1}{\eta _\mathrm{r}}\right)^{2\alpha 1}\frac{\eta _\mathrm{r}}{\eta }\right],$$
(3.9)
which is valid for $`\eta _r<\eta <\eta _\mathrm{m}`$. The time $`\eta _\mathrm{m}`$ marks the moment where $`Hm`$. When $`mt>1`$, the regime of coherent oscillations takes over and the solution (3.8) implies that the energy density stored in $`\varphi `$ decreases as $`a^3`$, meaning that $`\varphi _\mathrm{c}(\eta )\eta ^{3/2}`$. Since the coherent oscillations decrease as $`a^3`$ there will be a typical curvature scale $`H_\mathrm{c}`$ and a typical time $`\eta _\mathrm{c}`$ at which the coherent oscillations become dominant with respect to the radiation background. This moment is determined by demanding that

$$H_\mathrm{r}^2M_\mathrm{P}^2\left(\frac{a_\mathrm{r}}{a_\mathrm{c}}\right)^4m^2\varphi _1^2\left(\frac{a_\mathrm{m}}{a_\mathrm{c}}\right)^3,$$
(3.10)
which also implies that

$$H_\mathrm{c}m\phi ^4,$$
(3.11)
where $`\phi =\varphi _1/M_\mathrm{P}`$. Eq. (3.11) has been obtained without tuning the asymptotic value of $`\varphi `$ to the minimum of its potential. If such a tuning is made, the amplitude of oscillations at $`\eta _\mathrm{m}`$ will be smaller than $`\varphi _1`$ and it will be given, according to Eq. (3.9), by $`\varphi _2(\eta _1/\eta _\mathrm{r})^{2\alpha 1}`$. Thus, the scale $`H_c`$ will be defined by a different relation namely:

$$m^2\varphi _2^2\left(\frac{H_\mathrm{r}}{H_1}\right)^{\frac{2(2\alpha 1)}{\alpha +1}}\left(\frac{m}{H_\mathrm{r}}\right)\left(\frac{a_\mathrm{m}}{a_\mathrm{c}}\right)^3H_\mathrm{r}^2M_\mathrm{P}^2\left(\frac{a_\mathrm{r}}{a_\mathrm{c}}\right)^4,$$
(3.12)
leading, ultimately, to

$$H_\mathrm{c}=m\left(\frac{\varphi _2}{M_\mathrm{P}}\right)^4\left(\frac{H_\mathrm{r}}{H_1}\right)^{\frac{4(2\alpha 1)}{\alpha +1}}\left(\frac{m}{H_\mathrm{r}}\right)^2.$$
(3.13)
In the approximation of instantaneous reheating \[i.e. $`\eta _\mathrm{r}\eta _1`$\], $`H_\mathrm{r}H_1`$. Therefore, from Eq. (3.13) $`H_\mathrm{c}`$ is smaller than the value determined in Eq. (3.11) by a factor $`(m/H_1)^2`$. In the approximation of matter-dominated reheating (i.e. $`\alpha 2`$), the result of instantaneous reheating is further suppressed by a factor $`(H_\mathrm{r}/H_1)^4`$ as one can easily argue from Eq. (3.13). From Eqs. (3.5)–(3.6), the evolution of $`\varphi `$ will go as $`\eta ^3`$ when coherent oscillations start dominating.

In spite of the possible tunings made in the asymptotic values of $`\varphi `$, after $`\eta _\mathrm{c}`$ there will be a typical time at which the field $`\varphi `$ will decay. In order not to spoil the light elements produced at the epoch of BBN $`\varphi `$ has to decay at a scale larger than $`H_{\mathrm{ns}}T_{\mathrm{ns}}^2/M_\mathrm{P}`$ (where $`T_{\mathrm{ns}}\mathrm{MeV}`$. Since $`\varphi `$ is only coupled gravitationally the typical decay scale will be given by comparing the rate with the curvature scale giving that

$$H_\varphi \mathrm{\Gamma }\frac{m^3}{M_\mathrm{P}^2}>H_{\mathrm{ns}},$$
(3.14)
implying that $`m>10^4`$ GeV. This requirement also demands that the reheating temperature associated with the decay of $`\varphi `$ will be larger than the BBN temperature.

In order to illustrate some concrete examples of the various possibilities implied by our considerations, suppose that $`m10^3`$ TeV and suppose that the asymptotic value of $`\varphi `$ is fine-tuned to its minimum. Furthermore, suppose that the reheating is instantaneous. Then, according to the picture which has been presented, inflation stops at a scale $`H_110^{13}`$ GeV and $`\varphi `$ starts oscillating at a curvature scale $`H_\mathrm{m}10^3`$ TeV. The coherent oscillations will then become dominant at a curvature scale $`H_\mathrm{c}10^8`$ GeV (having assumed $`\varphi _2M_\mathrm{P}`$). The coherent oscillations of $`\varphi `$ will last down to $`H_\varphi 10^{20}`$ GeV. After this moment the Universe will be dominated by the radiation produced in the decay of $`\varphi `$. Notice, for comparison, that the BBN curvature scale is $`H_{\mathrm{ns}}10^{25}`$ GeV so that the decay occurs well before BBN (five orders of magnitude in curvature scale).

Another illustrative example is the one where $`m10^6`$ TeV. In this case the decay of $`\varphi `$ occurs prior to the EWPT epoch, namely

$$H_\varphi >H_{\mathrm{ew}}.$$
(3.15)
In fact $`H_{\mathrm{ew}}=\sqrt{N_{\mathrm{eff}}}T_{\mathrm{ew}}^2/M_\mathrm{P}10^{17}`$ GeV (with $`N_{\mathrm{eff}}=106.75`$ and $`T_{\mathrm{ew}}100`$ GeV) whereas, from Eq. (3.14), $`H_\varphi 10^9`$ GeV. In more general terms we can say that in order to have the $`\varphi `$ decay occurring prior to the EWPT epoch we have to demand that $`H_\varphi >H_{\mathrm{ew}}`$ which means that $`m>10^5`$ TeV.

In closing this section two general comments are in order. If no fine-tuning is made in the asymptotic amplitude of $`\varphi `$, the typical scale of the coherent oscillations will almost coincide with $`m`$. However, the possibility $`\phi 1`$ is still left if, for some reason, we want $`H_\mathrm{c}m`$.

The decay of $`\varphi `$ and the consequent freezing of the gauge coupling should occur prior to the EWPT epoch and the baryon number should be generated, in the present context, at the electroweak time. Suppose, for example, that this is not the case and that the BAU has been created prior to the electroweak scale. Suppose, moreover, that the decay of $`\varphi `$ occurs after baryogenesis. Then the temperature of the radiation gas before the decay of $`\varphi `$ will be $`T_\varphi T_\mathrm{m}(a_\mathrm{m}/a_\varphi )m(m/M_\mathrm{P})^5`$. Thus, the entropy increase due to the decay of $`\varphi `$ will be $`\mathrm{\Delta }S(T_{\mathrm{decay}}/T_\varphi )^3`$ where $`T_{\mathrm{decay}}\sqrt{H_\varphi M_\mathrm{P}}`$. This implies that $`\mathrm{\Delta }Sm/M_\mathrm{P}`$. It has been observed in different contexts that in order to preserve a pre-existing BAU one should have $`\mathrm{\Delta }S<10^5`$ . Thus, this bound would imply $`m>10^{14}`$ GeV. This is the reason why the present analysis will assume that the decay of $`\varphi `$ occurs prior to the electroweak time and that the BAU is generated at the EWPT or shortly after.

### 3.2 Evolution of the inhomogeneous modes

When the Universe passes from the inflationary stage to the subsequent radiation dominated expansion, inhomogeneities of the field $`\varphi `$ are generated. This may invalidate the original assumptions and introduce further complications by adding qualitatively new constraints on the scenario.

It is useful to recall that the inhomogeneities of $`\varphi `$ can be interpreted, in the framework of second quantization, as quanta of the field $`\varphi `$. Hence, the inhomogeneities produced because of the sudden change of the geometry from the de Sitter epoch to the radiation dominated epoch can be counted by estimating the number of quanta produced by the sudden change of the geometry according to the well known techniques of curved space-times .

Consider the first order fluctuations of the field $`\varphi `$

$$\varphi (\stackrel{}{x},\eta )=\varphi (\eta )+\delta \varphi (\stackrel{}{x},\eta ),$$
(3.16)
whose evolution equation is, in Fourier space,

$$\psi ^{\prime \prime }+2\psi ^{}+[k^2+m^2a^2]\psi =0,$$
(3.17)
where $`\psi (k,\eta )`$ is the Fourier component of $`\delta \varphi (\stackrel{}{x},\eta )`$. In order to count the number of quanta produced during the transition of the geometry from the inflationary to the radiation dominated stage of expansion the (canonically normalized) amplitude of fluctuations $`\mathrm{\Psi }=\psi a`$ should be defined so that Eq. (3.17) becomes:

$$\mathrm{\Psi }^{\prime \prime }+\left[k^2+m^2a^2\frac{a^{\prime \prime }}{a}\right]\mathrm{\Psi }=0.$$
(3.18)
In the de Sitter stage of expansion Eq. (3.18) reduces to

$$\mathrm{\Psi }_\mathrm{i}^{\prime \prime }+\left[k^2+\frac{\mu ^22}{\eta ^2}\right]\mathrm{\Psi }_\mathrm{i}=0,$$
(3.19)
whereas during the radiation dominated stage of expansion Eq. (3.19) takes the form

$$\mathrm{\Psi }_\mathrm{r}^{\prime \prime }+\left[k^2+\frac{\mu ^2(\eta +2\eta _1)^2}{\eta _1^4}\right]\mathrm{\Psi }_\mathrm{r}=0.$$
(3.20)
The solution of Eq. (3.19) (with the correct quantum-mechanical normalization for $`\eta \mathrm{}`$) can be written as

$$\mathrm{\Psi }_\mathrm{i}(\eta )=\frac{1}{\sqrt{2k}}p\sqrt{x}H_\rho ^{(1)}(x),$$
(3.21)
where $`x=k\eta `$ and $`H_\nu ^{(1)}`$ is the first order Hankel function . In the pure de Sitter case, $`\rho =3/2\sqrt{1(4/9)\mu ^2}`$ and since $`\mu 1`$, $`\rho 3/2`$; $`p`$ is a phase factor which has been chosen in such a way that

$$p=\sqrt{\frac{\pi }{2}}e^{i\frac{\pi }{4}(1+2\rho )}.$$
(3.22)
With this choice of $`p`$ we have that $`\mathrm{\Psi }_\mathrm{i}(\eta )e^{ik\eta }/\sqrt{2k}`$ for $`\eta \mathrm{}`$.

During the radiation dominated stage of expansion Eq. (3.20) is the equation of parabolic cylinder functions . The solutions turning into positive and negative frequencies for $`\eta +\mathrm{}`$ are then

$`f_\mathrm{r}(\eta )={\displaystyle \frac{1}{(2\gamma )^{1/4}}}e^{i\frac{\pi }{8}}D_{iq\frac{1}{2}}(ie^{i\frac{\pi }{4}}z),`$
$`f_\mathrm{r}^{}(\eta )={\displaystyle \frac{1}{(2\gamma )^{1/4}}}e^{i\frac{\pi }{8}}D_{iq\frac{1}{2}}(e^{i\frac{\pi }{4}}z),`$ (3.23)
where

$$z=\sqrt{2\gamma }(\eta +2\eta _1),q=\frac{k^2}{2\gamma },$$
(3.24)
and where $`D_\sigma `$ are the parabolic cylinder functions in the Whittaker’s notation. The solution of Eq. (3.20)

$$\mathrm{\Psi }_\mathrm{r}(\eta )=c_+(k)g_\mathrm{r}(\eta )+c_{}(k)g_\mathrm{r}^{}(\eta ),$$
(3.25)
is given in terms of $`c_+(k)`$ and $`c_{}(k)`$ which are the two (complex) Bogoliubov coefficients satisfying $`|c_+(k)|^2|c_{}|^2=1`$. In a second quantized approach $`|c_{}(k)|^2`$ is the mean number of created quanta, whereas in a semi-classical approach $`c_{}(k)`$ can be viewed as the coefficient parametrising the mixing between positive and negative frequency modes. In the case $`c_{}(k)0`$ no mixing takes place and no amplification is produced. In order to determine $`c_\pm (k)`$, $`\mathrm{\Psi }_\mathrm{i}(\eta )`$ and $`\mathrm{\Psi }_\mathrm{r}(\eta )`$ should be continuously matched in $`\eta =\eta _1`$, namely

$`\mathrm{\Psi }_\mathrm{i}(\eta _1)=\mathrm{\Psi }_\mathrm{r}(\eta _1),`$
$`\mathrm{\Psi }_\mathrm{i}^{}(\eta _1)=\mathrm{\Psi }_\mathrm{r}^{}(\eta _1),`$ (3.26)
By solving this system, an exact expression for the Bogoliubov coefficients is obtained which is, in general a function of two variables : $`\mu =m\eta _1`$ and $`x_1=k\eta _1`$. Since $`\mu 1`$ the exact result can be expanded, in this limit,

$`c_+(k)=\pi e^{i\frac{\pi }{8}}\left\{{\displaystyle \frac{i}{\sqrt{2}\mathrm{\Gamma }(\frac{3}{4})}}S_2(x_1,\rho )\mu ^{\frac{1}{4}}+{\displaystyle \frac{(1+i)}{2\mathrm{\Gamma }(\frac{1}{4})}}[S_1(x_1,\rho )+S_2(x_1,\rho )]\mu ^{\frac{1}{4}}\right\}+𝒪(\mu ^{\frac{5}{4}}),`$
$`c_{}(k)=\pi e^{i\frac{\pi }{8}}\left\{{\displaystyle \frac{i}{\sqrt{2}\mathrm{\Gamma }(\frac{3}{4})}}S_2(x_1,\rho )\mu ^{\frac{1}{4}}+{\displaystyle \frac{(i1)}{2\mathrm{\Gamma }(\frac{1}{4})}}[S_1(x_1,\rho )+S_2(x_1,\rho )]\mu ^{\frac{1}{4}}\right\}+𝒪(\mu ^{\frac{5}{4}}),`$ (3.27)
where $`S_1(x_1,\rho )`$ and $`S_2(x_1,\rho )`$ contain the explicit dependence upon the Hankel’s functions:

$`S_1(x_1,\rho )=e^{i\frac{\pi }{4}(1+2\rho )}H_\rho ^{(1)}(x_1),`$
$`S_2(x_1,\rho )=\sqrt{x_1}e^{i\frac{\pi }{4}(1+2\rho )}\left[\left(\rho +{\displaystyle \frac{1}{2}}\right){\displaystyle \frac{H_\rho ^{(1)}(x_1)}{\sqrt{x_1}}}\sqrt{x_1}H_{\rho +1}^{(1)}(x_1)\right].`$ (3.28)
If $`\rho 3/2`$ Eq. (3.28) gives

$`c_+(k)=e^{i\frac{\pi }{8}}\sqrt{\pi }\left\{{\displaystyle \frac{x_1^{\frac{3}{2}}}{2\mathrm{\Gamma }(\frac{3}{4})\mu ^{\frac{1}{4}}}}+{\displaystyle \frac{ix_1^{\frac{1}{2}}}{2\mathrm{\Gamma }(\frac{3}{4})\mu ^{\frac{1}{4}}}}+\left[{\displaystyle \frac{1}{2\mathrm{\Gamma }(\frac{3}{4})\mu ^{\frac{1}{4}}}}+{\displaystyle \frac{(i1)\mu ^{1/4}}{\sqrt{2}\mathrm{\Gamma }(\frac{1}{4})}}\right]\sqrt{x_1}\right\}+𝒪(\mu ^{\frac{5}{4}}),`$
$`c_{}(k)=e^{i\frac{\pi }{8}}\sqrt{\pi }\left\{{\displaystyle \frac{x_1^{\frac{3}{2}}}{2\mathrm{\Gamma }(\frac{3}{4})\mu ^{\frac{1}{4}}}}{\displaystyle \frac{ix_1^{\frac{1}{2}}}{2\mathrm{\Gamma }(\frac{3}{4})\mu ^{\frac{1}{4}}}}+\left[{\displaystyle \frac{1}{2\mathrm{\Gamma }(\frac{3}{4})\mu ^{\frac{1}{4}}}}+{\displaystyle \frac{(i+1)\mu ^{1/4}}{\sqrt{2}\mathrm{\Gamma }(\frac{1}{4})}}\right]\sqrt{x_1}\right\}+𝒪(\mu ^{\frac{5}{4}}).`$ (3.29)
In the limit $`x_1=k\eta _11`$ the mean number of created quanta can be finally approximated as

$$\overline{n}(k)|c_{}(k)|^2=q|k\eta _1|^{2\rho }\mu ^{1/2}$$
(3.30)
where $`q`$ is a numerical coefficient of the order of $`10^2`$. The energy density of the created (massive) quanta can be estimated from

$$d\rho _\psi =\frac{d^3\omega }{(2\pi )^3}m\overline{n}(k)$$
(3.31)
where $`\omega =k/a`$ is the physical momentum. In the case of a de Sitter phase ($`\rho =3/2`$) the typical energy density of the produced fluctuations is

$$\rho _\psi (\eta )qmH_1^3\left(\frac{m}{H_1}\right)^{1/2}\left(\frac{a_1}{a}\right)^3$$
(3.32)
The produced massive quanta may become dominant. If they become dominant after $`\varphi `$ already decayed they will not lead to further constraints on the scenario. If they become dominant prior to the decay of $`\varphi `$ further constraints may be envisaged. The scale at which the massive fluctuations become dominant with respect to the radiation background can be determined by requiring that $`\rho _\psi (\eta _{})\rho _\gamma (\eta _{})`$ implying that

$$qmH_1^3\left(\frac{m}{H_1}\right)^{1/2}\left(\frac{a_1}{a_{}}\right)^3H_1^2M_\mathrm{P}^2\left(\frac{a_1}{a_{}}\right)^4,$$
(3.33)
which translates into

$$H_{}q^2mϵ^4,$$
(3.34)
where $`ϵ=H_1/M_\mathrm{P}`$. In order to make sure that the non-relativistic modes will become dominant after $`\varphi `$ already decayed $`H_{}<H_\varphi `$ should be imposed, that is to say $`m>10^2`$ TeV for $`ϵ10^6`$.

The maximum tolerable amount of entropy, in order not to wash-out any preexisting BAU is model-dependent but, in general, $`\mathrm{\Delta }S<10^5`$ seems to be acceptable . Defining $`T_\varphi `$ as the radiation gas already present at the scale $`H_\varphi `$, the entropy increase from $`T_\varphi `$ to $`T_{\mathrm{decay}}\sqrt{H_\varphi M_\mathrm{P}}`$ is of the order of

$$\mathrm{\Delta }S=\left(\frac{T_{\mathrm{decay}}}{T_\varphi }\right)^3,$$
(3.35)
where

$$T_\varphi =T_{}\left(\frac{a_{}}{a_\varphi }\right)m\xi ^{1/6}ϵ^{1/2},$$
(3.36)
where $`\xi =m/M_\mathrm{P}`$. Demanding that $`\mathrm{\Delta }S<10^5`$ implies that

$$\xi >10^{10}ϵ^3.$$
(3.37)
Taking, as usual, $`ϵ10^6`$, $`m>10`$ GeV.

Hence, if the constraints pertaining to the homogeneous mode are enforced, the bounds coming from the inhomogeneous modes do not invalidate the conclusions of the analysis. According to the logic expressed in Eq. (3.15), an illustrative example is the case where $`m>10^5`$ TeV and the BAU is generated after EWPT. In this case the bounds obtained in the present section are satisfied and the analysis of the evolution of the inhomogeneous modes shows that the qualitatively new bounds introduced in the picture are less constraining than the ones obtained in the analysis of the dynamics of the homogeneous mode.

## 4 Evolution of the gauge field fluctuations

The evolution of the field $`\varphi `$ during and after the de Sitter stage implies, according to Eqs. (3.3)–(3.9), that the two-point function of the gauge field fluctuations may very well grow.

For $`\eta <\eta _1`$, a rough approximation suggests that the effect of the ohmic current is not present and the gauge field is in the vacuum state with $`k/2`$ energy in each of its modes. By promoting the classical fields to quantum mechanical operators we have that the physical polarizations of the magnetic field can be written as

$$\widehat{b}_i(\stackrel{}{x},\eta )=\frac{d^3k}{(2\pi )^{3/2}}\underset{\alpha }{}e_i^\alpha \left[\widehat{a}_{k,\alpha }b(k\eta )e^{i\stackrel{}{k}\stackrel{}{x}}+\widehat{a}_{k,\alpha }^{}b^{}(k\eta )e^{i\stackrel{}{k}\stackrel{}{x}}\right],$$
(4.1)
where $`b(k\eta )=B(k\eta )/g(\eta )`$ obey the equation

$$b^{\prime \prime }+\left[k^22\left(\frac{g^{}}{g}\right)^2+\frac{g^{\prime \prime }}{g}\right]b=0.$$
(4.2)
Notice that $`b(k\eta )`$ are the correct normal modes whose limit (for $`\eta \mathrm{}`$) should be normalized to $`\sqrt{k/2}e^{ik\eta }`$. The two-point correlation function of the magnetic fluctuations can then be expressed as

$$𝒢_{ij}(\stackrel{}{r},\eta )\widehat{b}_i(\stackrel{}{x},\eta )\widehat{b}_j(\stackrel{}{x}+\stackrel{}{r},\eta )=\frac{d^3k}{(2\pi )^3}P_{ij}(k)b(k,\eta )b^{}(k,\eta )e^{i\stackrel{}{k}\stackrel{}{r}},$$
(4.3)
where

$$P_{ij}(k)=\left(\delta _{ij}\frac{k_ik_j}{k^2}\right).$$
(4.4)
The magnetic energy density, derived from the energy-momentum tensor corresponding to the action of Eq. (2.2), is related to the trace of the correlation function reported in Eq. (4.1) over the physical polarizations:

$$\rho _B(r,\eta )=\rho _B(k,\eta )\frac{\mathrm{sin}kr}{kr}\frac{dk}{k}$$
(4.5)
where

$$\rho _B(k,\eta )=\frac{1}{\pi ^2}k^3|b(k,\eta )|^2.$$
(4.6)
A necessary condition in order to assess that gauge field fluctuations grow during the de Sitter phase is that the two-point function increases in the limit $`\eta \eta _1`$ .

Suppose that the gauge coupling decreases with monotonic dependence upon the field $`\varphi `$, namely

$$g(\eta )=\left(\frac{\varphi \varphi _1}{M_\mathrm{P}}\right)^{\frac{\lambda }{2}},\lambda >0.$$
(4.7)
This parametrisation is purely phenomenological, however, it allows to take into account, at once, some physically interesting cases like the one suggested by the low-energy string effective action where, in the limit of $`\varphi /M_\mathrm{P}<1`$, $`g^2(\varphi )\varphi `$.

Using Eq. (3.3) and Eq. (4.7) into Eq. (4.2) the time evolution of the normal modes of the magnetic field can be found analytically since the specific form of Eq. (4.2) falls in the same category of Eq. (3.19). Hence,

$$b(k,\eta )=N\sqrt{k\eta }H_\nu ^{(2)}(k\eta ),N=\frac{\sqrt{k\pi }}{2}e^{i\frac{\pi }{4}(1+2\nu )},$$
(4.8)
with $`\nu =(3\lambda +1)/2`$ and where $`H_\nu ^{(2)}(k\eta )`$ the Hankel function of second kind . Notice that the normalization $`N`$ has been chosen in such a way that for $`\eta \mathrm{}`$ the correct quantum mechanical normalization is reproduced. Consequently, following Eq. (4.3), the two-point function evolves as

$$\underset{\eta \eta _1}{lim}𝒢_{ij}(r,\eta )\left|\frac{\eta }{\eta _1}\right|^{3\lambda }.$$
(4.9)
Since the two-point function increases magnetic fluctuations are generated.

For sake of completeness the case of increasing gauge coupling will now be examined using the following phenomenological parameterisation

$$g(\eta )=\left(\frac{\varphi \varphi _1}{M_\mathrm{P}}\right)^{\frac{\delta }{2}},\delta >0.$$
(4.10)
Again different scenarios can be imagined. For instance, one coould argue in favour of scenarios where the gauge coupling depends upon $`\varphi `$ (or upon $`\eta `$) in a highly non-monotonic way. For the illustrative purposes of the present investigation it is however sufficient to focus the attention on the case of monotonic dependence.

Following now the same steps outlined in the case of decreasing gauge coupling the evolution of the two-point function can be obtained

$$\underset{\eta \eta _1}{lim}𝒢_{ij}(r,\eta )\left|\frac{\eta }{\eta _1}\right|^{23\delta }.$$
(4.11)
for $`\delta >1/3`$ and

$$\underset{\eta \eta _1}{lim}𝒢_{ij}(r,\eta )\left|\frac{\eta }{\eta _1}\right|^{3\delta }.$$
(4.12)
for $`\delta <1/3`$. If $`\delta <1/3`$ the correlation function decreases and this signals that large scale magnetic fields are not produced.

The back-reaction of the produced fluctuations can be safely neglected in de Sitter space. Looking at Eqs. (2.3) and (2.11) it can happen that if the magnetic fluctuations grow too much the term at the right hand side will become of the same order of the others. This is not the case. Using the conventions of this Section together with the explicit form of the scale factor in the de Sitter phase it can be shown that

$$\frac{1}{g^3a^2}\frac{g}{\varphi }\stackrel{}{B}^2\left|\frac{\eta }{\eta _1}\right|^{2\nu }|k\eta _1|^{52\nu }$$
(4.13)
where $`\nu `$ is determined from the specific power dependence of the coupling as a function of $`\varphi `$. Since we are interested in large scale modes we have $`k\eta _11`$. Therefore the back reaction effects are relevant towards the end of the de Sitter phase (i.e. $`\eta \eta _1`$ ) and for $`k\eta _1^1`$, namely exactly for the modes not relevant for the present investigation.

After the end of inflation the onset of the conductivity dominated regime may not be instantaneous. In this case after $`\eta _1`$ the presence of a reheating phase should be taken into account. Suppose, for instance, that $`g(\eta )`$ decreases according to Eq. (4.7). Suppose, moreover, that during reheating the background is dominated by the coherent oscillations of the inflaton. In this case the effective evolution of the geometry will be dominated by matter with $`a(\eta )\eta ^2`$. According to Eq. (3.4) (with $`\alpha 2`$), $`\varphi \eta ^3`$. The value of the magnetic inhomogeneities at $`\eta _\mathrm{r}`$ will be given by solving Eq. (4.2)

$$b(k,\eta _\mathrm{r})A_1\left(\frac{\eta _\mathrm{r}}{\eta _1}\right)^{\frac{3}{2}\lambda }+A_2\left(\frac{\eta _\mathrm{r}}{\eta _1}\right)^{1\frac{3}{2}\lambda },$$
(4.14)
where $`A_1`$ and $`A_2`$ are two arbitrary constants. Depending upon the value of $`\lambda `$ the fastest growing solution is selected in Eq. (4.14). This phase may induce further amplification on the two-point function since its main effect is to delay the conductivity-dominated regime.

### 4.1 Generalized MHD equations

For $`\eta >\eta _r`$ the role of the Ohmic diffusion becomes important and the evolution of the magnetic inhomogeneities will be described by the MHD equations generalized to the case of time varying gauge coupling.

The ordinary (i.e. fixed coupling) MHD treatment is an effective description valid for length scales larger than the Debye radius and for frequencies smaller than the iono-acoustic frequency. This means that MHD is accurate in reproducing the spectrum of plasma excitations obtained from the full kinetic (Vlasov-Landau) approach but only for sufficiently low frequencies and for sufficiently large scales.

Implicit in the ordinary MHD analysis is the assumption that the plasma has to be electrically neutral ($`\stackrel{}{}\stackrel{}{E}=0`$) over length scales larger than the Debye radius. Thus, this system of equations cannot be applied for distances shorter than the Debye radius and for frequencies larger than the plasma frequency where a kinetic description should be employed.

MHD equations can be derived from a microscopic (kinetic) approach and also from a macroscopic approach where the displacement current is neglected . If the displacement current is neglected the electric field can be expressed using the Ohm law and the magnetic diffusivity equation can be derived

$$\frac{\stackrel{}{B}}{\eta }=\stackrel{}{}\times (\stackrel{}{v}\times \stackrel{}{B})+\frac{1}{\sigma }^2\stackrel{}{B}.$$
(4.15)
The term containing the bulk velocity field is called dynamo term and it receives contribution provided parity is globally broken over the physical size of the plasma. In Eq. (4.15) the contribution containing the conductivity is usually called magnetic diffusivity term.

In the superconducting (or ideal) approximation the resistivity of the plasma goes to zero and the induced (Ohmic) electric field is orthogonal both to the bulk velocity of the plasma and to the magnetic field \[i. e. $`\stackrel{}{E}\stackrel{}{v}\times \stackrel{}{B}`$\]. In the real (or resistive) approximation the resistivity may be very small but it is always finite and the Ohmic field can be expressed as

$$\stackrel{}{E}\frac{\stackrel{}{}\times \stackrel{}{B}}{\sigma }\stackrel{}{v}\times \stackrel{}{B}.$$
(4.16)
If the gauge coupling changes with time the system of equations obtained by neglecting the displacement current receives new contributions and the relevant equations can be obtained, in the resistive approximation:

$`\stackrel{}{}\times \stackrel{}{E}={\displaystyle \frac{\stackrel{}{B}}{\eta }},`$ (4.17)
$`\stackrel{}{E}={\displaystyle \frac{\stackrel{}{J}}{\sigma }}\stackrel{}{v}\times \stackrel{}{B},`$ (4.18)
$`{\displaystyle \frac{1}{g^2}}\stackrel{}{}\times \stackrel{}{B}=\stackrel{}{J}2{\displaystyle \frac{g^{}}{2g^3}}\stackrel{}{E},`$ (4.19)
Using Eqs. (4.17)–(4.19) the generalized magnetic diffusivity equation can be obtained:

$$\left(1\frac{2}{\sigma g^2}\frac{g^{}}{g}\right)\frac{\stackrel{}{B}}{\eta }=\stackrel{}{}\times (\stackrel{}{v}\times \stackrel{}{B})+\frac{1}{\sigma g^2}^2\stackrel{}{B}.$$
(4.20)
Notice that Eq. (4.20) reproduces Eq. (4.15) if $`g^{}0`$.

Suppose now that the plasma, whose effective Ohmic description has been presented, is relativistic. In the case when the coupling is constant $`\sigma `$ is constant and it is given by

$$\sigma \sigma _c(\eta )a(\eta ),$$
(4.21)
where $`\sigma _cT/g^2`$ scales as the inverse of $`a(\eta )`$ if the evolution of the Universe is, to a good approximation, adiabatic.

If $`g`$ is not constant, $`\sigma `$ is not constant anymore but it decreases if the gauge coupling increases and, vice versa, it increases if the gauge coupling decreases. In spite of this, in the generalized MHD equations, the combination which appears is always $`\sigma g^2`$ which is roughly constant for an adiabatically expanding Universe.

In the approximation of instantaneous reheating, the solution of Eq. (4.20) is given by

$$B(k,\eta )=B(k,\eta _1)e^{{\scriptscriptstyle {\scriptscriptstyle \frac{k^2}{\sigma g^22\frac{g^{}}{g}}}𝑑\eta }}.$$
(4.22)
According to Eqs. (3.13)–(3.14), in order to get the coupling frozen prior to $`H_{\mathrm{ew}}10^{17}`$ GeV, $`m>10^5`$ TeV shall be required. If the gauge coupling is always decreasing as a function of $`\eta `$, it can be parametrized by Eq. (4.7). Hence Eq. (4.22) can be evaluated by using the explicit evolution of $`\varphi `$ as obtained from Eqs. (3.8)–(3.9) implying that $`\varphi _\mathrm{r}\eta ^1`$, $`\varphi _\mathrm{m}\eta ^{3/2}`$ and $`\varphi _\mathrm{c}\eta ^3`$. The result is

$$B(k,\eta _0)=(\eta _1,\eta _\mathrm{m},\eta _\varphi ,\eta _\mathrm{c})e^{\frac{k^2}{\sigma g^2}(\eta _1+\eta _0)}$$
(4.23)
where $`\eta _0`$ is the present time and

$$(\eta _1,\eta _\mathrm{m},\eta _\varphi ,\eta _\mathrm{c})=\left[\left(\frac{\lambda +\sigma g^2\eta _1}{\lambda +\sigma g^2\eta _\mathrm{m}}\right)\left(\frac{3\lambda +2\sigma g^2\eta _\mathrm{m}}{3\lambda +2\sigma g^2\eta _\mathrm{c}}\right)^{\frac{3}{2}}\left(\frac{3\lambda +\sigma g^2\eta _\mathrm{c}}{3\lambda +\sigma g^2\eta _\varphi }\right)^3\right]^{\lambda \left[\frac{k}{\sigma g^2}\right]^2}.$$
(4.24)
Concerning Eqs. (4.23)–(4.24) few comments are in order. From Eq. (4.23) all the modes

$$k^2>k_\sigma ^2\frac{\sigma g^2}{\eta _0}$$
(4.25)
are suppressed by the effect of the conductivity. The present value of $`\omega _\sigma (\eta _0)`$ 444 With $`\omega (\eta )k/a(\eta )`$ the physical momentum will be denoted. can be estimated by recalling that $`1/\eta _0H_0a_0`$ where $`H_010^{61}M_\mathrm{P}`$. Thus $`\omega _\sigma 10^3`$ Hz. Present modes of the magnetic fields are dissipated if $`\omega >\omega _\sigma `$.

As far as the problem of galactic magnetic fields is concerned, the relevant set of scales range around the Mpc corresponding to present modes of the magnetic field $`\omega _\mathrm{G}10^{14}`$ Hz, i.e. $`\omega _\mathrm{G}\omega _\sigma `$.

## 5 Estimates of large scale magnetic fields

In the approximation of instantaneous reheating the typical (present) frequency corresponding to the end of inflation can be computed and it turns out to be

$$\omega _1(\eta _0)z_{\mathrm{dec}}^1T_{\mathrm{dec}}ϵ^{1/2}\xi ^{1/3}\phi ^{\frac{2}{3}}.$$
(5.1)
Since $`T_{\mathrm{dec}}0.26`$ eV , $`z_{\mathrm{dec}}^1T_{\mathrm{dec}}100\mathrm{GHz}`$. Eq. (5.1) can be obtained by red-shifting the highest mode, i.e. $`\omega _1(\eta _1)H_1`$ through the different stages of the evolution of the model, namely, according to Eqs. (3.11) and (3.14), from $`\eta _1`$ down to $`\eta _\mathrm{m}`$ and from $`\eta _\mathrm{m}`$ down to $`\eta _{rmc}`$. Recall that from $`\eta _\mathrm{c}`$ to $`\eta _\varphi `$ the Universe is, effectively, matter dominated. The other typical frequencies appearing in the time evolution of the gauge coupling can be written, in units of $`\omega _1(\eta _0)`$, as

$`{\displaystyle \frac{\omega _\mathrm{m}(\eta _0)}{\omega _1(\eta _0)}}=ϵ^{1/2}\xi ^{1/2},`$
$`{\displaystyle \frac{\omega _\mathrm{c}(\eta _0)}{\omega _1(\eta _0)}}=\xi ^{1/2}ϵ^{1/2}\phi ,`$
$`{\displaystyle \frac{\omega _\varphi (\eta _0)}{\omega _1(\eta _0)}}=ϵ^{1/2}\xi ^{7/6}\phi ^{2/3},`$ (5.2)
where, as in the case of Eq. (5.1) all the frequencies are evaluated at the present time.

In the case of decreasing gauge coupling \[described by Eq. (4.7)\] the amount of generated large scale magnetic field can be estimated from Eqs. (4.2)–(4.8) together with Eqs. (4.20)–(4.23). Bearing in mind that the typical frequency scale corresponding to $`1`$ Mpc is $`10^{14}`$ Hz we have that the ratio of the magnetic to radiation energy density is:

$$r_B(\omega _\mathrm{G})=f(\lambda )ϵ^{\frac{3}{2}\lambda }\xi ^{\lambda \frac{4}{3}}\phi ^{2\lambda \frac{8}{3}}\mathrm{\hspace{0.17em}\hspace{0.17em}10}^{25(43\lambda )}𝒯(\omega _\mathrm{G}),$$
(5.3)
where

$$f(\lambda )=\frac{2^{3\lambda 1}}{\pi ^3}\mathrm{\Gamma }^2\left(\frac{3}{2}\lambda +\frac{1}{2}\right),$$
(5.4)
and

$$𝒯(\omega _\mathrm{G})e^{\frac{\omega _\mathrm{G}^2}{\omega _\sigma ^2}}\left[\left(\frac{\omega _\varphi }{\omega _1}\right)\left(\frac{\omega _\varphi }{\omega _\mathrm{m}}\right)^{\frac{1}{2}}\left(\frac{\omega _\varphi }{\omega _\mathrm{c}}\right)^{\frac{3}{2}}\right]^{\lambda \frac{\omega _\mathrm{G}^2}{T_0^2}},$$
(5.5)
which means, using Eqs. (5.1)–(5.2),

$$𝒯(\omega _\mathrm{G})=e^{\frac{\omega _\mathrm{G}^2}{\omega _\sigma ^2}}\left[ϵ^{1/2}\phi ^1\xi ^{5/2}\right]^{\lambda \frac{\omega _\mathrm{G}^2}{T_0^2}}.$$
(5.6)
Notice that in Eq. (5.5) is the present CMB temperature.

In order to illustrate the regions of the parameter space where magnetogenesis is possible $`\phi 1`$ will be assumed. In Fig. 1 and Fig. 2 the case of decreasing gauge coupling is discussed. In Fig. 1, $`\lambda `$ is fixed and the exclusion plot is given in terms of $`ϵ`$ and $`\xi `$. In particular, for illustration, $`\lambda =1`$ has been chosen. The choice of $`\lambda 1`$ implies that $`g^2(\varphi )\varphi `$. Such a case is favoured from the tree-level string effective action where the effective coupling can be approximated as $`g^2(\varphi )\varphi /M_\mathrm{P}`$ for $`\varphi <M_\mathrm{P}`$.

The two vertical lines lines mark, respectively, the bounds coming from BBN \[i.e. $`\xi >10^{15}`$, from Eq. (3.14)\] and from the electroweak epoch \[i.e. $`\xi >10^{11}`$, from Eq. (3.15)\]. For consistency with the assumptions of Eqs. (3.2)–(3.3) $`m/H_11`$, implying $`ϵ\xi `$. With the lower dashed line the bound $`ϵ>\xi `$ is reported. Notice that the requirement of Eq. (3.37) is not numerically relevant. The upper dashed line is obtained by requiring, according to Eq. (3.34), that $`H_\varphi >H_{}`$. Finally, the full (diagonal) line is derived by imposing on Eqs. (5.3)–(5.6) the dynamo requirement.

In Fig. 2 the value of $`\xi `$ has been fixed to $`10^{11}`$, as required by the considerations related to the electroweak epoch and the magnetogenesis region is described in terms of $`ϵ`$ and $`\lambda `$, both varying over their physical range. The two horizontal lines fix the bounds coming from $`ϵ10^6`$ and from $`ϵ>\xi `$. With the dashed line the curve $`r_B(\omega _\mathrm{G})10^{20}`$ is denoted. Hence, values larger than the dynamo requirement are allowed. This observation may be relevant in the context of magnetic fields associated with clusters.

In Fig. 2, $`\lambda `$ lies in the range $`0<\lambda <4/3`$. This choice guarantees the growth of the correlation function of the magnetic inhomogeneities during the de Sitter stage.

In order not to conflict with large scale bounds coming from the isotropy of the CMB the energy spectra of the produced gauge field fluctuations have to decay at large distance scales, implying that $`0<\lambda 4/3`$. To be compatible with CMB anisotropies $`r_\mathrm{B}(\omega _{\mathrm{dec}})10^{10}`$ should be imposed (where $`\omega _{\mathrm{dec}}10^{16}`$ Hz). If the condition $`0<\lambda <4/3`$ is enforced, the spectra increase with frequency and the possible bounds coming from the anisotropy of the CMB are satisfied.

An analogous estimate can be obtained in the case of increasing gauge coupling discussed in Eq. (4.10). In this case

$$r_\mathrm{B}(\omega _\mathrm{G})f(\delta )ϵ^{\frac{3}{2}\delta 1}\xi ^{\delta 2}\phi ^{42\delta }\mathrm{\hspace{0.17em}\hspace{0.17em}10}^{25(63\delta )}𝒯(\omega _\mathrm{G}),$$
(5.7)
where

$$f(\delta )=\frac{2^{3\delta 3}}{\pi ^3}\mathrm{\Gamma }^2\left(\frac{3}{2}\delta \frac{1}{2}\right),$$
(5.8)
and

$$𝒯(\omega _\mathrm{G})=e^{\frac{\omega _\mathrm{G}^2}{\omega _\sigma ^2}}\left[ϵ^{1/2}\phi ^1\xi ^{5/2}\right]^{\delta \frac{\omega _\mathrm{G}^2}{T_0^2}}.$$
(5.9)
In Fig. 3 and Fig. 4 the requirements coming from the dynamo mechanism as well as the other theoretical constraints are illustrated. For both plots $`ϵ<10^6`$ and $`\phi 1`$. In Fig. 3 the case $`\delta =2`$ is illustrated. The shaded area selects the region of the parameter space where the dynamo requirement is imposed on Eq. (5.7). The two vertical lines illustrate the conditions of Eqs. (3.14)–(3.15). As in Fig. 2 the two dashed lines correspond to the constraints coming from the inhomogeneous modes and from the condition $`ϵ>\xi `$.

The requirement that the spectra decrease at large distance scales implies, in this case that $`\delta 2`$. The condition on the growth of the correlation function obtained in Eqs. (4.11)–(4.12) imply $`\delta >1/3`$. Thus the interesting physical range of Eq. (4.10) and (5.7) will be $`1/3<\delta 2`$.

In Fig. 4 the parameter space is illustrated for fixed values of $`\xi `$, i.e. $`\xi 10^{11}`$. As in Fig. 3 the full and dashed curves correspond to the dynamo requirements imposed on Eq. (5.7). The shaded area selects the allowed region in the space of the parameters where magnetogenesis is possible for $`10^{11}<ϵ<10^6`$. As in the case of Fig. 2 there are regions in the shaded area where values much larger than the magnetogenesis requirement are possible (see the dashed line in Fig. 4).

As it has been pointed out in deriving the theoretical bounds on the scenario the requirement $`\xi 10^{11}`$ may be too restrictive since it excludes the variation of the gauge coupling at the electroweak time. To relax this assumption is possible and it would require a precise analysis of the dynamics of the EWPT in the presence of time varying gauge coupling. At the moment this kind of analysis is not available.

Summarizing this illustrative discussion, there are regions in the parameter space of the models wher all the theoretical constraints are satisfied and where magnetogenesis is possible. In particular, it is possible that the gauge coupling freezes prior to the electroweak epoch, leading still to magnetogenesis.

## 6 Concluding remarks

There are no reasons why the gauge couplings should be constant throughout all the history of the Universe. If they are allowed to change prior to the formation of the light elements they can lead to computable differences in the cosmological evolution.

In the present paper the interplay between inflationary magnetogenesis and the evolution of the Abelian gauge coupling has been addressed. In a phenomenologically reasonable model of inflationary and post-inflationary evolution the relaxation of the gauge coupling leads to a growth in the correlation function of magnetic inhomogeneities. Large scale magnetic fields are then generated. The evolution of the gauge coupling is driven, in the present context, by a massive scalar.

Since the gauge coupling evolves significant changes in the evolution of the plasma can be envisaged. In the present investigation the ordinary MHD equations have been generalized to the case of time evolving gauge coupling but other effects could be envisaged. In particular, the generlization of the full kinetic approach to the case of time evolving “electron” charge would be of related interest.

The value of large scale magnetic fields produced with this mechanism has been estimated. For a broad range of parameters the obtained values of the magnetic fields are much larger than the dynamo requirements.

In the present investigation the main assumption has been that the only gauge coupling free to evolve is the Abelian one. Furthermore, the parameters of the model have been chosen in such a way that the gauge coupling is not dynamical by the onset of the electroweak phase transition. It would be interesting to relax both assumptions since they may lead to potential differences with the standard scenarios.

Therefore, in many respects, the present investigation is not conclusive. At the same time it shows that acceptable models for the evolution of the gauge couplings can be obtained in a standard cosmological framework.

## Acknowledgements

The author wishes to thank N. Sanchez and H. de Vega for their kind invitation to the Chalonge School and for useful remarks on the subject of large scale magnetic fields. It is also a pleasure to thank interesting conversations with Yoel Raphaeli, Dan Boyanovsky, Yuri Pariiski, Sasha Burinskii and Peter Coles.


format: text

text:
# Is there a Glass Transition in Planar Vortex Systems?

## Abstract

The criteria for the existence of a glass transition in a planar vortex array with quenched disorder are studied. Applying a replica Bethe ansatz, we obtain for self-avoiding vortices the exact quenched average free energy and effective stiffness which is found to be in excellent agreement with recent numerical results for the related random bond dimer model Zeng et al. (1999). Including a repulsive vortex interaction and a finite vortex persistence length $`\xi `$, we find that for $`\xi 0`$ the system is at all temperatures in a glassy phase; a glass transition exists only for finite $`\xi `$. Our results indicate that planar vortex arrays in superconducting films are glassy at presumably all temperatures.

The existence of a vortex glass (VG) in type-II superconductors is of crucial importance for their performance. Quenched disorder affects both the translational order of the Abrikosov vortex lattice and the low current response Nattermann and Scheidl (2000). The thermodynamic response of a planar VG in a mesoscopic superconducting film Abrikosov (1964) has been studied in a recent experiment Bolle et al. (1999), providing ‘fingerprints’ of the sample’s disorder Emig and Kardar (2000). Glassy phases are expected to be generic to a wide class of disorder dominated systems. The planar VG is a famous starting point within this class since it is among the very few systems where analytical progress seems possible. The situation is here analogous to that of exactly solvable 1D quantum systems which have advanced the understanding of strongly correlated systems.

Thermal fluctuations compete against the VG, raising the question of the existence of a glass transition at a finite temperature $`T_g`$ Fisher (1989). However, the literature on this subject is rather contradictory Nattermann et al. (1991); Hwa and Fisher (1994); Zeng et al. (1999). Major support, both analytical and numerical, for a finite $`T_g`$ stems from the phenomenologically related random field (RF) XY model which has a finite $`T_g`$ transition Cardy and Ostlund (1982); Giamarchi and Le Doussal (1994). However, considering the entropic vortex interaction, the absence of a finite $`T_g`$ transition has been predicted under certain assumptions as well Nattermann et al. (1991). The theory underlying the vortex system can be applied to a broad class of other systems like surface steps Mullins and Sekerka (1964), polymers de Gennes (1968) and domain walls Pokrovsky and Talapov (1979). A more fundamental aspect of the planar VG is related to its often debated replica description, including the concept of replica symmetry breaking, and the non-trivial exchange of replica and thermodynamic limit. Due to the relation to the recently simulated random bond dimer model Zeng et al. (1999), the VG provides a unique system to test replica theory and could be useful in clarifying important concepts in the statistical mechanics of disordered systems.

In this paper, we analyze the conditions under which a planar vortex lattice undergoes a finite $`T_g`$ glass transition. For vanishing vortex persistence length $`\xi `$, we find the lattice to be glassy at all temperatures with an analytic free energy. Analyticity of the free energy is also obtained for the vortex-free RFXY model. Our replica symmetric Bethe ansatz results are found to be in excellent agreement with recent simulations Zeng et al. (1999). For finite $`\xi `$, there exists a finite $`T_g`$ which, however, is presumably too large to be observable in superconducting films.

Vortex system. — We consider a lattice of self-avoiding directed elastic lines of density $`\rho =1/a`$ in $`1+1`$ dimensions. Each line is characterized by its position $`x_j(z)`$ and line tension $`g`$. The lines interact via a short-ranged repulsive pair potential $`U(x)`$. Quenched disorder couples locally to the lines via the random potential $`V(𝐫)`$, $`𝐫=(x,z)`$, which is assumed to have short-ranged correlations $`\overline{V(𝐫)V(𝐫^{})}=\mathrm{\Delta }\delta (𝐫𝐫^{})`$. The total energy reads

$$H=𝑑z\underset{j}{}\left\{\frac{g}{2}\left(\frac{dx_j}{dz}\right)^2+\underset{i,ij}{}U(x_ix_j)+V(x_j,z)\right\}.$$
(1)
This model can be considered as the continuum version of lines placed on the bonds of a rectangular lattice Kardar and Nelson (1985). Below, we will also consider lines with a finite persistence length $`\xi `$, corresponding to correlated random walks in the lattice model.

Random field XY model. — The line lattice is usually treated within a 2D elastic theory for the displacement field $`u(𝐫)`$ so that the line positions are $`x_j(z)=ja+u(ja,z)`$. The effective Hamiltonian is then

$$H_{\mathrm{XY}}=d^2r\left\{\frac{c_{11}}{2}(_xu)^2+\frac{c_{44}}{2}(_zu)^2+\rho (𝐫)V(𝐫)\right\}$$
(2)
with compression modulus $`c_{11}=aU^{\prime \prime }(a)`$, tilt modulus $`c_{44}=g/a`$ and line density $`\rho (𝐫)=_j\delta (xx_j(z))`$. Rewriting $`\rho (𝐫)`$ as a Fourier series in $`u(𝐫)`$, the model becomes equivalent to the random field XY model without vortices Fisher (1989). It is well known that, for weak disorder, this model shows a phase transition to a glassy phase at a finite temperature $`T=T_g`$. A renormalization group (RG) analysis yields the scaling dimension of disorder Cardy and Ostlund (1982),

$$\lambda _\mathrm{\Delta }=2\left(1\frac{\pi T}{a^2\sqrt{c_{11}c_{44}}}\right),$$
(3)
giving $`T_g=a^2\sqrt{c_{11}c_{44}}/\pi `$. An important property of the XY model is that the elastic moduli remain unrenormalized due to a statistical symmetry Cardy and Ostlund (1982); Hwa and Fisher (1994). Thus, naively one would expect for the line lattice a transition at a $`T_g`$ determined by the elastic moduli on the length scale where $`H_{XY}`$ becomes applicable. However, as we will show below, there is, in fact, a renormalization of $`c_{11}`$ in the original model of Eq. (1) which can even destroy the transition.

What are the manifestations of the glass transition in the XY model? Anomalous variations in the response to an external tilt of the lines have been used to characterize the glassy phase Hwa and Fisher (1994). However, signatures of the transition in the free energy itself have not been studied so far. Using the known mapping of the XY model onto a vector Coulomb gas Cardy and Ostlund (1982) one can derive the RG flow of the free energy as shown in Kosterlitz (1974) for the pure XY model. From this analysis we find that the quenched average free energy is analytic at $`T_g`$. Weak disorder makes the contribution

$$f_\mathrm{\Delta }=\frac{\mathrm{\Delta }}{Ta^2}\left(1+\frac{\mathrm{\Delta }}{c_{11}^{3/2}c_{44}^{1/2}a^4}\frac{\pi ^3}{(1\tau )^2(12\tau )}\right)$$
(4)
to the free energy density for small $`\tau =1T/T_g`$, where we have neglected $`T`$-independent terms.

Another signature of the glassy phase are its correlations Nattermann and Scheidl (2000). For $`T>T_g`$, one finds $`[(u(𝐫)u(\mathrm{𝟎}))^2]\mathrm{log}(|𝐫|/a)`$, and in the glassy phase for small $`\tau >0`$ the slightly faster growth $`[(u(𝐫)u(\mathrm{𝟎}))^2]\mathrm{log}^2(|𝐫|/a)`$. Here $`\mathrm{}`$ and $`[\mathrm{}]`$ denote thermal and disorder average, respectively. The decay of the corresponding density-density correlations suggests that the correlation length is infinite on both sides of the transition which is consistent with the predicted analyticity of free energy. Next we will focus on the correlations of thermal fluctuations $`\delta u(𝐫)=u(𝐫)u(𝐫)`$ around the state pinned by disorder. On both sides of the transition one has

$$C(𝐫)=[(\delta u(𝐫)\delta u(\mathrm{𝟎}))^2]=\frac{a^2}{16\pi K}\mathrm{log}(|𝐫|/a),$$
(5)
with $`K=(\pi /16)T_g/T=(a^2/16)\sqrt{c_{11}c_{44}}/T`$. Below, we will calculate $`K`$ directly in terms of the parameters of the original model of Eq. (1).

Self-avoiding vortices. — A generic but simpler model is obtained if the only vortex interaction of Eq. (1) is the one representing the non-crossing condition, i.e., a hardcore repulsion $`U(x)=c\delta (x)`$ with $`c\mathrm{}`$. Using the replica method, with $`n`$ replicas, the model can be mapped to $`SU(n)`$ fermions and thus can be solved exactly by a replica Bethe ansatz (RBA) without replica-symmetry breaking Kardar and Nelson (1985). To date, only deep inside the glassy phase, e.g., for $`HH_{c1}`$, the quenched average free energy density $`[f]`$ has been studied by this technique Kardar and Nelson (1985); Emig and Kardar (2001). However, at arbitrary $`H`$ (or density $`\rho `$) the RBA yields

$$[f]=f_0\rho +\frac{T^2}{2g}_Q^Qq^2\varrho (q)𝑑q,$$
(6)
where $`f_0`$ is the single vortex free energy per unit length; its precise form is unimportant in what follows. $`\varrho (q)`$ is determined by

$$_Q^Q\left\{\frac{1}{l_d(qq^{})}+\pi \mathrm{coth}\left(\pi l_d(qq^{})\right)\right\}\varrho (q^{})𝑑q^{}=q$$
(7)
with length scale $`l_d=T^3/g\mathrm{\Delta }`$, and $`Q`$ is fixed by $`_Q^Q\varrho (q)𝑑q=\rho `$. The integral equation (7) can be solved perturbatively in $`Ql_d`$, leading to

$$\varrho (q)=\sqrt{1(q/Q)^2}\left[\frac{1}{2\pi }Ql_d\frac{\pi }{24}(Ql_d)^3+\mathrm{}\right].$$
(8)
Higher order terms can be easily calculated and will depend also on $`q`$, thus changing the functional form of the lowest order result. In the opposite limit, $`Ql_d\mathrm{}`$, i.e, for vanishing disorder, the solution is $`\varrho (q)=1/(2\pi )`$. With the result for $`\varrho (q)`$ at hand, both $`\rho l_d`$ and the integral of Eq. (6) can be obtained as a series in $`(Ql_d)^2`$. However, we are actually interested in eliminating $`Q`$ in order to calculate $`[f]`$ in terms of $`\rho `$ and $`l_d`$. At this stage, a particularly interesting property of $`\varrho (q)`$ is important. The exact result for the integral of Eq. (6) as a function of $`\rho l_d`$ is obtained from the two lowest order terms of Eq. (8). Indeed, it can be shown order by order in $`Ql_d`$ that the exact result is $`_Q^Qq^2\varrho (q)𝑑k=\frac{\pi ^2}{3}\rho ^3+\rho ^2l_d^1`$. This result we could also confirm numerically with high precision. Thus, the RBA gives for the quenched average free energy density the exact result

$$[f]=f_0\rho +\frac{\pi ^2}{6}\frac{T^2}{g}\rho ^3+\frac{\mathrm{\Delta }}{2T}\rho ^2.$$
(9)
For the interaction (last two) terms it is possible to give a simple explanation. The first is just the known Pokrovski-Talapov interaction Pokrovsky and Talapov (1979) from entropy loss due to collisions of thermally fluctuating lines in pure systems. The second term can be explained by a similar argument which gives the interaction in the dilute limit $`\rho l_d^1`$. It arises from an increase in energy due to collisions which prevent the lines from gaining their optimal energy from the random potential Kardar and Nelson (1985). It is rather unexpected that the total effective interaction is just the sum of these two contributions. As a consequence of that, the free energy is an analytic function. Based on our result for the XY model, this does not exclude the possibility of a finite-$`T`$ glass transition. However, its existence can be excluded by computing the scaling dimension $`\lambda _\mathrm{\Delta }`$, Eq. (3), in terms of the macroscopic compression modulus $`c_{11}^0`$. It is given by the effective interaction, $`c_{11}^0=aU_{\mathrm{eff}}^{\prime \prime }(a)`$, with $`U_{\mathrm{eff}}/a`$ corresponding to the two last terms of Eq. (9). Using $`c_{44}=g/a`$, and

$$c_{11}^0=\frac{\pi ^2T^2}{ga^3}+\frac{\mathrm{\Delta }}{Ta^2},$$
(10)
we find $`\lambda _\mathrm{\Delta }>0`$ at all finite temperatures, i.e., $`T_g=\mathrm{}`$. With the exact $`c_{11}^0`$ available, we are able to relate the phenomenological stiffness $`K`$ of Eq. (5) to the microscopic parameters of the model of Eq. (1), leading to

$$K=\frac{\pi }{16}\left(1+\frac{ag\mathrm{\Delta }}{\pi ^2T^3}\right)^{1/2}.$$
(11)
For vanishing disorder, $`\mathrm{\Delta }=0`$, the previously known result $`K=\pi /16`$ is recovered Henley (1997).

There are a few comments in order to clarify the apparent discrepancy between the exact RBA results for self-avoiding lines and the RG predictions for the random field XY model of Eq. (2). To begin with, on microscopic length scales there is no direct link between the two models since $`c_{11}=0`$ due to self-avoidance as only interaction. However, the XY model is generally believed to describe the vortex lattice at length scales larger than $`a`$ where a finite fluctuation induced $`c_{11}`$ exists Fisher (1989); Nattermann et al. (1991). That this point of view is incorrect for finding out about the existence of a glass transition in the vortex system shows the following argument. On large but finite length scales $`La^2g/T`$ the effective $`c_{11}`$ for $`\mathrm{\Delta }0`$ can be obtained from the analogy of the self-avoiding vortex lines to 1D fermions Kardar and Nelson (1985). Since a finite vortex length $`L`$ translates into a finite temperature for the fermions, the Sommerfeld expansion can be used to get the correction $`\delta c_{11}=ag/3L^2`$ to the result of Eq. (10). Thus, on any sufficiently large length scale, disorder seems to be irrelevant since then $`\lambda _\mathrm{\Delta }<0`$. In contrast, the RBA shows that disorder is always relevant. The contradicting results are explained by the fact that the statistical symmetry in the XY model is broken on all length scales by the non-crossing condition.

Dimer model. — The self-avoiding vortices can be simulated directly, i.e., not via the XY model, using a mapping to a discrete dimer model with quenched bond disorder Zeng et al. (1996). A recently developed polynomial algorithm allows to study large systems at finite temperature but fixed density $`\rho =1/2`$ (in units of the dimer lattice constant) Zeng et al. (1999). Comparison of our RBA predictions with the results of this numerical approach provides a critical test of the replica approach, including replica-symmetry and proper analytical continuation while taking the number of replicas to zero. The dimer model is defined as follows: Choose a subset (whose elements are called dimers) of the bonds on a square lattice, such that every lattice site is touched by exactly one of these dimers. The energy of one such dimer covering $`D`$ is defined by $`H_d=_{(i,j)D}ϵ_{ij}/T_d`$, where the sum is over all dimers of $`D`$, $`ϵ_{ij}`$ are random bond energies drawn from a Gaussian distribution with zero mean and a variance of one, and $`T_d`$ is the dimer temperature. Via the well-known mapping to the solid-on-solid model, it is possible to assign every dimer covering a height function which maps to the line displacement $`u(𝐫)`$ in the vortex system representation. This representation results from the superposition of an arbitrary covering $`D`$ and a fixed reference covering that specifies the direction of the lines, cf. Ref. Zeng et al. (1999). By measuring the height configurations, Zeng et al. could determine the correlation function of Eq. (5), and thus the stiffness $`K`$ as a function of $`T_d`$ Zeng et al. (1999). In order to compare their results to our prediction of Eq. (11), the parameters $`g/T`$ and $`\mathrm{\Delta }/T^2`$ have to be determined for the dimer representation of the vortex system. It is easy to see that the ratio of horizontal and vertical line steps on the lattice yields $`g/T=2`$. The precise relation between $`\mathrm{\Delta }`$ and $`T_d`$ can be obtained by comparing the annealed average free energy of the continuum model of Eq. (1) and the discrete dimer version of the vortex lattice. Using known methods for computing the partition function of dimer models Kasteleyn (1961), we get

$$\frac{\mathrm{\Delta }}{\xi _dT^2}=\frac{2}{T_d^2}\frac{4G}{\pi }+\frac{4}{\pi }_0^{e^{T_d^2/2}}𝑑x\frac{\mathrm{arctan}(x)}{x},$$
(12)
with $`\xi _d`$ a regularization length for the disorder correlator not , and Catalan’s constant $`G=0.915966`$. With the above relations and $`a=1/\rho =2`$ in Eq. (11), we obtain $`K(T_d)`$ as plotted in Fig. 1, together with the data provided by the dimer algorithm of Ref. Zeng et al. (1999). Obviously, there is excellent agreement between our RBA results and simulations.

Finite interaction and persistence length. — So far, we have shown that the model of Eq. (1) with a vanishing single line persistence length $`\xi `$ and $`U(x)`$ representing only a non-crossing condition is always in the glassy phase. We now consider the effects of an arbitrary repulsive interaction $`U(x)`$ and a finite persistence length $`\xi `$. For the vortex system, $`U(x)`$ diverges for $`x0`$ and decays exponentially for large $`x`$; $`\xi `$ can be identified with the coherence length of the superconductor.

Let us first discuss a general interaction while keeping $`\xi =0`$. The relevance of disorder now depends on the interaction induced change in $`c_{11}`$ for $`\mathrm{\Delta }0`$, cf. Eq. (3). This change can be obtained by employing the known mapping Pokrovsky and Talapov (1979) of elastic lines to 1D spin-less fermions which however interact via $`U(x)`$, denoting now the finite distance interaction on top of the self-avoidance. The compressibility $`\kappa `$ of the fermions is related to the compression modulus of the vortex system according to $`1/\kappa =a^2c_{11}`$. For weak interaction, $`\kappa `$ is easily calculated, e.g., by bosonization. Using the standard notation of Luttinger liquid theory Voit (1995), one gets $`\kappa =1/\pi \mathrm{}v_N`$ with $`v_N=v_F+(\stackrel{~}{U}(0)+\stackrel{~}{U}(2\pi /a)/2)/\mathrm{}\pi `$, where $`\stackrel{~}{U}(q)`$ is the Fourier transform of $`U(x)`$. Thus, any repulsive interaction increases the effective compression modulus $`c_{11}=c_{11}^{\mathrm{PT}}+c_{11}^{\mathrm{int}}`$ with

$$c_{11}^{\mathrm{int}}=\frac{1}{a^2}\left(\stackrel{~}{U}(0)+\frac{1}{2}\stackrel{~}{U}(2\pi /a)\right),$$
(13)
and for $`U(x)0`$, the known Pokrovsky-Talapov result $`c_{11}=c_{11}^{\mathrm{PT}}=\pi ^2T^2/ga^3`$ is recovered. Thus according to Eq. (3), the scaling dimension of disorder is increased by any finite vortex interaction, and the system remains always in the glassy phase.

Next, we neglect interactions for the time being and now consider the thermally fluctuating vortex lines as correlated random walks, i.e., we assume a finite persistence length $`\xi `$. These correlations increase the fluctuations of a single line. Denoting the so decreased macroscopic line tension again by $`g`$, entropy reduction due to line collisions leads not only to the Pokrovsky-Talapov result but also to an additional renormalization of the modulus $`c_{11}=c_{11}^{\mathrm{PT}}c_{11}^\xi `$ where, for $`\xi a^2g/T`$,

$$c_{11}^\xi =\frac{10}{3}\frac{T\xi }{a^2g}c_{11}^{\mathrm{PT}}.$$
(14)
For larger $`\xi `$, $`c_{11}`$ increases relative to this first order result. Including the effect of disorder on $`c_{11}`$, cf. Eq. (10), one obtains from Eq. (3) the finite transition temperature

$$T_g=0.42\left(\frac{g^2\mathrm{\Delta }}{\rho ^3\xi }\right)^{1/4}.$$
(15)
Thus, a finite persistence length yields a finite temperature glass transition.

Finally, we study the combined effect of interactions and a finite persistence length. If disorder is strong compared to the interaction, $`\mathrm{\Delta }T(\stackrel{~}{U}(0)+\stackrel{~}{U}(2\pi /a)/2)`$, the result of Eq. (15) remains applicable. In the opposite limit, $`T_g`$ is independent of disorder,

$$T_g=0.31\left(\frac{g^2(\stackrel{~}{U}(0)+\stackrel{~}{U}(2\pi /a)/2)}{\rho ^3\xi }\right)^{1/3}.$$
(16)
For vanishing persistence length, $`T_g`$ diverges, and the system is glassy at all temperatures.

Discussion. — The model of Eq. (1) applies to vortex systems in films of thickness $`da,\lambda `$ with $`\lambda `$ the penetration depth. Then the crucial questions is if the glass transition temperature $`T_g`$ can be smaller than the critical temperature $`T_c`$ of the superconductor. Since the only mechanism leading to a finite $`T_g`$ is a finite persistence length, we maximize the effect from $`c_{11}^\xi `$ in Eq. (14) by considering small $`a\lambda `$. Then Ginzburg-Landau theory yields a $`\stackrel{~}{U}(q)`$ Abrikosov (1964) which is too large for the bosonization result of Eq. (13) to be valid. However, a lower bound is $`c_{11}^{\mathrm{int}}>bc_{11}^{\mathrm{PT}}`$ with $`b`$ of order one. Using Eq. (14), we find $`T_g>3b/10\times a^2g/\xi `$. With the stiffness $`g=\mathrm{\Phi }_0^2/16\pi \times d/a\lambda ^2`$, $`\mathrm{\Phi }_0=hc/2e`$, Abrikosov (1964), we obtain in the limit of infinitesimal disorder $`T_g\mathrm{\Phi }_0^2/16\pi ^2\times ad/\xi \lambda ^2`$ which becomes, with $`a`$, $`d>\xi `$, $`T_g/T_c>1510^4`$ for $`T_c=100`$K and with the length $`\xi _0=10`$Å, $`\lambda _0=100`$Å at $`T=0`$. Thus, $`T_g`$ is extremely close to $`T_c`$, and presumably any small amount of quenched disorder will raise $`T_g`$ above $`T_c`$, making it impossible to observe a glass transition in planar vortex systems.

We thank M. Kardar, T. Nattermann, S. Scheidl and C. Zeng for useful discussions, and especially C. Zeng for sending the data of Fig. 1. This work was supported by DFG through Emmy-Noether grant No. EM70/2-1.


text:
# Indications of Neutrino Oscillation in a 250 km Long-baseline Experiment

The K2K Collaboration

## Abstract

The K2K experiment observes indications of neutrino oscillation: a reduction of $`\nu _\mu `$ flux together with a distortion of the energy spectrum. Fifty-six beam neutrino events are observed in Super-Kamiokande (SK), 250 km from the neutrino production point, with an expectation of $`80.1_{5.4}^{+6.2}`$. Twenty-nine one ring $`\mu `$-like events are used to reconstruct the neutrino energy spectrum, which is better matched to the expected spectrum with neutrino oscillation than without. The probability that the observed flux at SK is explained by statistical fluctuation without neutrino oscillation is less than 1%.

Recent atmospheric and solar neutrino data indicate the existence of neutrino oscillation and therefore the existence of neutrino mass Fukuda et al. (1998a, 2002); Ahmad et al. (2002). The zenith angle distribution of atmospheric neutrinos observed by Super-Kamiokande (SK) shows a clear deficit of upward-going muon neutrinos, which is well explained by two-flavor $`\nu _\mu `$-$`\nu _\tau `$ oscillations with $`\mathrm{\Delta }m^2`$ around $`3\times 10^3\mathrm{eV}^2`$, and $`\mathrm{sin}^22\theta `$ close to or equal to unity.

The KEK to Kamioka long-baseline neutrino oscillation experiment (K2K) Ahn et al. (2001) uses an accelerator-produced neutrino beam with a neutrino flight distance of 250 km to probe the same $`\mathrm{\Delta }m^2`$ region as that explored with atmospheric neutrinos. The neutrino beam is produced by a 12 GeV proton beam from the KEK proton synchrotron. After the protons hitting an aluminum target, the produced positively charged particles, mainly pions, are focused by a pair of pulsed magnetic horns Ieiri et al. (1997). The neutrinos produced from the decays of these particles are 98% pure muon neutrinos with a mean energy of 1.3 GeV.

This analysis is based on data taken from June 1999 to July 2001, corresponding to $`4.8\times 10^{19}`$ protons on target (POT). The pion momentum and angular distributions downstream of the second horn are occasionally measured with a gas-Cherenkov detector (PIMON) Maruyama (2000) in order to verify the beam Monte Carlo (MC) simulation and to estimate the errors on the flux prediction at SK. The direction of the beam is monitored on a spill-by-spill basis by observing the profile of the muons from the pion decays with a set of ionization chambers and silicon pad detectors located just after the beam dump. The neutrino beam itself is measured in a set of near neutrino detectors (ND) located 300 m from the proton target. The measurements made at the ND are used to verify the stability and the direction of the beam, and to determine the flux normalization and the energy spectrum before the neutrinos travel the 250 km to SK. The flux at SK is estimated from the flux of the ND by multiplying the Far/Near (F/N) ratio, the ratio of fluxes between the far detector (SK) and ND, to that of the ND.

Since both a suppression in the number of events and a distortion of the spectrum are expected for neutrinos which travel a fixed path length in the presence of oscillations, both the number of observed events and the spectral shape information at SK are compared with expectation. All of the beam-induced neutrino events observed within the fiducial volume of SK are used to measure the overall suppression of flux. In order to study the spectral distortion, 1 ring $`\mu `$-like events ($`1R\mu `$) are selected to enhance the fraction of charged-current (CC) quasi-elastic (QE) interactions ($`\nu _\mu +n\mu +p`$). Only the muon is visible in these reactions since the proton momentum is typically below Cherenkov threshold. The energy of the parent neutrino can be calculated by using the observed momentum of the muon, assuming QE interactions, and neglecting Fermi momentum:

$`E_\nu ^{\mathrm{rec}}={\displaystyle \frac{m_NE_\mu m_\mu ^2/2}{m_NE_\mu +P_\mu \mathrm{cos}\theta _\mu }},`$ (1)
where $`m_N`$, $`E_\mu `$, $`m_\mu `$, $`P_\mu `$ and $`\theta _\mu `$ are the nucleon mass, muon energy, the muon mass, the muon momentum and the scattering angle relative to the neutrino beam direction, respectively.

The ND is comprised of two detector systems: a 1 kiloton water Cherenkov detector (1KT) and a fine-grained detector (FGD) system. The flux normalization is measured by the 1KT to estimate the expected number of events at SK. Since the 1KT has the same detector technology as SK, most of systematic uncertainties on the measurement are canceled. The energy spectrum is measured by analyzing the muon momentum and angular distributions in both detector systems. The 1KT has a high efficiency for reconstructing the momentum of muons below 1 $`\mathrm{GeV}/c`$, and full 4$`\pi `$ coverage in solid angle. However, the 1KT has little efficiency for reconstructing muons with momentum above 1.5 $`\mathrm{GeV}/c`$ since they exit the detector. The FGD, on the other hand, has high efficiency for measuring muons above 1 $`\mathrm{GeV}/c`$, and the two complementary detectors are thus able to completely cover the relevant energy range.

In the 1KT analysis, a cylindrical fiducial volume of 25 tons oriented along the beam direction is used. Event selection criteria for the flux normalization are the same as those in reference Ahn et al. (2001). Events which deposit more than $`100`$ MeV of energy are used for the measurement of the integrated flux. The measurement has a 5% systematic uncertainty, of which the largest contribution comes from the vertex reconstruction Ahn et al. (2001). For the spectrum measurement, further cuts are imposed in order to select $`1R\mu `$ events which stop in the detector. Among the events selected for the flux normalization measurement, 53% of the events have 1 ring. The events with a muon exiting from the detector are effectively eliminated by requiring the maximum charge of any PMT to be less than 200 photo-electrons; 68% of the 1 ring events remain after this requirement. The largest systematic uncertainty for the spectrum measurement is an uncertainty on the energy scale. The energy scale is understood to within $`{}_{3}{}^{}{}_{}{}^{+2}`$ %, which is confirmed with both cosmic-ray muons and beam-induced $`\pi ^0`$s.

The FGD is comprised of a scintillating fiber and water detector (SciFi) Suzuki et al. (2000), a lead-glass calorimeter (LG), and a muon range detector (MRD) Ishii et al. (2002). In the FGD analysis, events containing one or two tracks with vertex within the 5.9 ton fiducial volume of the SciFi are used. The track finding efficiency is 70% for a track passing through three layers of scintillating fiber and close to 100% for more than 5 layers Kim et al. (2002). Three layers is the minimum track length required in this analysis. Events which have at least one track passing into the MRD are chosen in order to select $`\nu _\mu `$-induced CC interactions. The momentum of the track is measured by its range through the SciFi, LG, and MRD, with accuracy of 3%.

If the proton produced in the QE interaction has a momentum greater than 600 MeV/$`c`$, its track may also be reconstructed. In the case where a second track is visible, the kinematic information is used to enhance the fraction of QE events in the sample. Assuming QE interaction, the direction of the proton can be predicted from the muon momentum. The QE enhanced sample is selected by requiring that the direction of the second track agrees with the prediction to within 25 degrees. Events where the direction of the second track differs from the prediction by more than 30 degrees are put into a non-QE enhanced sample. In the QE enhanced sample, 62% of the events are estimated to be QE events from the MC simulation. In the non-QE enhanced sample, $`82`$% of events are estimated to come from interactions other than QE. The SciFi events are divided into three event categories: 1-track, 2-track QE enhanced, and 2-track non-QE enhanced samples.

The 2-dimensional distributions of the muon momentum versus angle with respect to the beam direction of four event categories (the 1KT event sample and the three SciFi event samples) are used to measure the neutrino spectrum. A $`\chi ^2`$-fitting method is used to compare these data against the MC expectation. The neutrino spectrum is divided into 8 energy bins as defined in Table 1. During the fit, the flux in each energy bin is re-weighted relative to the values in the beam MC. These weights are normalized relative to 1.00 for the $`E_\nu =1.0`$$``$1.5 GeV bin. An overall normalization is included as a free parameter in the fit. The parameter $`R_{nqe}`$ is used to re-weight the ratio between the QE and non-QE cross-section relative to the MC simulation. The systematic uncertainties of the ND are incorporated into the fitting parameters. They are the energy scales, the track finding efficiencies, and the detector thresholds. In addition, the spectrum measurement by PIMON is used as a constraint on the re-weighting factors.

The value of $`\chi ^2`$ is 227.2/197 d.o.f. at the best fit point. All the parameters including the detector systematics are found to lie within their expected errors. The best fit values of the flux re-weighting factors are shown in Table 1. The muon momentum and angular distributions of $`1R\mu `$ events in the 1KT, and the muon momentum distributions of the 2-track QE enhanced and non-QE enhanced events in SciFi are overlaid with the re-weighted MC in Figure 1. The fit result agrees well with the data. The errors of the measurement are provided in the form of an error matrix. Correlations between the parameters of the fit are taken into account in the oscillation analysis using this matrix. The diagonal elements in the error matrix are shown in Table 1.

The uncertainty due to neutrino interaction models is separately studied. In QE scattering, the axial vector mass in the dipole formula is set to a central value of 1.1 $`\mathrm{GeV}/c^2`$ and is varied by $`\pm `$10%. The axial mass for single pion production is set to a central value of 1.2 $`\mathrm{GeV}/c^2`$ and is varied by $`\pm `$20% Bernard et al. (2002). This affects both the $`q^2`$ dependence of the cross-section and the total cross-section. For coherent pion production, two different models are compared: one is the Rein and Sehgal model Rein and Sehgal (1983), and the other is a model by Marteau Marteau et al. (2000). For deep inelastic scattering, GRV94 Gluck et al. (1995) and the corrected structure function by Bodek and Yang Bodek and Yang (2002) are studied. For the oscillation analysis the Marteau model and Bodek and Yang structure functions are employed. Varying the choice of models causes the fitted value of $`R_{nqe}`$ ($`=0.93`$) to change by $`20\%`$. In order to account for this, an additional systematic error of $`\pm 20\%`$ on $`R_{nqe}`$ is added to the error matrix. The choice of models and axial mass does not affect the spectrum measurement itself beyond the size of the fitted errors.

The F/N ratio from the beam simulation is used to extrapolate the measurements at the ND to those at SK. The errors including correlations above 1 GeV, where the PIMON is sensitive, are estimated based on the PIMON measurements. The errors on the ratio for $`E_\nu `$ below 1 GeV are estimated based on the uncertainties in the hadron production models used in the K2K beam MC Ahn et al. (2001). The diagonal elements in the error matrix for the F/N ratio are summarized in Table 1.

The events in SK are selected using the timing information provided by the global positioning system. Events detected in SK must occur within an expected beam arrival time window of $`1.5`$ $`\mu \mathrm{sec}`$. In addition, the detected events must have no activity in outer detector, and have an energy deposit greater than 30 MeV with a vertex reconstructed within the 22.5 kiloton fiducial volume Ahn et al. (2001). This sample of events is referred to as the fully contained (FC) sample. The efficiency of this selection is 93% for CC interactions. Fifty-six events satisfy the criteria. With the timing cut, the expected number of atmospheric neutrino background is approximately $`10^3`$ events.

The expected number of FC events at SK without oscillation is estimated to be $`80.1_{5.4}^{+6.2}`$. The correlations between energy bins from the spectrum measurement at the ND and the F/N ratio are taken into account in the estimation of the systematic errors. The major contributions to the errors come from the uncertainties in the F/N ratio ( $`{}_{5.0\%}{}^{}{}_{}{}^{+4.9\%}`$) and the normalization (5.0%), dominated by uncertainties of the fiducial volumes due to vertex reconstruction both at the 1KT and SK.

A two flavor neutrino oscillation analysis, with $`\nu _\mu `$ disappearance, is performed by the maximum-likelihood method. In the analysis, both the number of FC events and the energy spectrum shape for $`1R\mu `$ events are used. The likelihood is defined as $`=_{norm}\times _{shape}`$. The normalization term $`_{norm}(N_{obs},N_{exp})`$ is the Poisson probability to observe $`N_{obs}`$ events when the expected number of events is $`N_{exp}(\mathrm{\Delta }m^2,\mathrm{sin}^22\theta ,f)`$. The symbol $`f`$ represents a set of parameters constrained by the systematic errors. These parameters are described in detail later. The shape term, $`_{shape}=_{i=1}^{N_{1R\mu }}P(E_i;\mathrm{\Delta }m^2,\mathrm{sin}^22\theta ,f),`$ is the product of the probability for each $`1R\mu `$ event to be observed at $`E_\nu ^{\mathrm{rec}}=E_i`$, where $`P`$ is the normalized $`E_\nu ^{\mathrm{rec}}`$ distribution estimated by MC simulation and $`N_{1R\mu }`$ is the number of $`1R\mu `$ events.

In the oscillation analysis, the whole data sample is used for $`_{norm}`$, i.e. $`N_{obs}=56`$. The spectrum shape in June 1999 was different from that for the rest of the running period because the target radius and horn current were different. The estimation of energy correlations in the spectrum at the ND and in the far/near ratio has not been completed for this period. Thus, data taken in June 1999 are discarded for $`_{shape}`$. The discarded data correspond to 6.5% of total POT. The number of $`1R\mu `$ events observed excluding the data of June 1999 is 29, and the corresponding number of $`1R\mu `$ events expected from MC simulation in the case of no oscillation is 44.

The parameters $`f`$ consist of the re-weighted neutrino spectrum measured at the ND ($`\mathrm{\Phi }_{ND}`$), the F/N ratio, the reconstruction efficiency ($`ϵ_{SK}`$) of SK for $`1R\mu `$ events, the re-weighting factor for the QE/non-QE ratio $`R_{nqe}`$, the SK energy scale and the overall normalization. The errors on the first 3 items depend on the energy and have correlations between each energy bin. The diagonal parts of their error matrices are summarized in Table 1 as described earlier. The error on the SK energy scale is 3% Fukuda et al. (1998b). Two different approaches are taken for the treatment of systematic errors in the likelihood. The first is to treat the parameters $`f`$ as fitting parameters with an additional constraint term in the likelihood (method 1) Fukuda et al. (1998a). The other approach is to average the $`(f)`$ sampled over many random trials weighted according to the probability density distribution of the systematic parameters $`f`$ (method 2) Swain and Taylor (1998).

The likelihood is calculated at each point in the $`\mathrm{\Delta }m^2`$ and $`\mathrm{sin}^22\theta `$ space to search for the point where the likelihood is maximized. The best fit point in the physical region of oscillation parameter space is found to be at ($`\mathrm{sin}^22\theta `$, $`\mathrm{\Delta }m^2`$)=$`(1.0,2.8\times 10^3\mathrm{eV}^2)`$ in method 1 and at $`(1.0,2.7\times 10^3\mathrm{eV}^2)`$ in method 2. If the whole space including the unphysical region is considered the values are $`(1.03,2.8\times 10^3\mathrm{eV}^2)`$ in method 1 and $`(1.05,2.7\times 10^3\mathrm{eV}^2)`$ in method 2. The results from two methods are consistent with each other. At the best fit point in the physical region the total number of predicted events is 54.2, which agrees with the observation of 56 within statistical error. The observed $`E_\nu ^{\mathrm{rec}}`$ distribution of the $`1R\mu `$ sample is shown in Figure 2 together with the expected distributions for the best fit oscillation parameters, and the expectation without oscillations. Consistency between the observed and best-fit $`E_\nu ^{\mathrm{rec}}`$ spectrum is checked by using Kolgomorov-Smirnov(KS) test. A KS probability of 79% is obtained. The best fit spectrum shape agrees with the observations.

The probability that the observations are due to a statistical fluctuation instead of neutrino oscillation is estimated by computing the likelihood ratio of the no-oscillation case to the best fit point. The no-oscillation probabilities are calculated to be 0.7% and 0.4% for method 1 and 2 respectively. When only normalization (shape) information is used, the probabilities are 1.3% (16%) and 0.7% (14.3%) for the two methods. Allowed regions of oscillation parameters are evaluated by calculating the likelihood ratio of each point to the best fit point, and are drawn in Figure 3. Both methods give essentially the same results. In order to be conservative, the result from method 1 is shown in the figure as it gives a slightly larger allowed region at the 99% C.L. The 90% C.L. contour crosses the $`\mathrm{sin}^22\theta =1`$ axis at $`1.5`$ and $`3.9\times 10^3`$ eV2 for $`\mathrm{\Delta }m^2`$. The oscillation parameters preferred by the total flux suppression and the energy distortions alone also agree well. Finally, the uncertainties of neutrino interactions are studied using the same procedure as the spectrum measurement at the ND. It is found that the effects of the interaction model difference on all the results are negligible due to the cancellation caused by using the same models in both the ND and SK.

In conclusion, both the number of observed neutrino events and the observed energy spectrum at SK are consistent with neutrino oscillation. The probability that the measurements at SK are explained by statistical fluctuation is less than 1%. The measured oscillation parameters are consistent with the ones suggested by atmospheric neutrinos. At the time of this letter the K2K experiment has collected approximately one-half of its planned $`10^{20}`$ protons on target.

###### Acknowledgements.

We thank the KEK and ICRR Directorates for their strong support and encouragement. K2K is made possible by the inventiveness and the diligent efforts of the KEK-PS machine and beam channel groups. We gratefully acknowledge the cooperation of the Kamioka Mining and Smelting Company. This work has been supported by the Ministry of Education, Culture, Sports, Science and Technology, Government of Japan and its grants for Scientific Research, the Japan Society for Promotion of Science, the U.S. Department of Energy, the Korea Research Foundation, and the Korea Science and Engineering Foundation.