text:
# Faster than Lyapunov decays of classical Loschmidt echo

## Abstract

We show that in the classical interaction picture the echo-dynamics, namely the composition of perturbed forward and unperturbed backward hamiltonian evolution, can be treated as a time-dependent hamiltonian system. For strongly chaotic (Anosov) systems we derive a cascade of exponential decays for the classical Loschmidt echo, starting with the leading Lyapunov exponent, followed by a sum of two largest exponents, etc. In the loxodromic case a decay starts with the rate given as twice the largest Lyapunov exponent. For a class of perturbations of symplectic maps the echo-dynamics exhibits a drift resulting in a super-exponential decay of the Loschmidt echo.

Analyzing the parametric stability of quantum dynamics through the so-called fidelity or Quantum Loschmidt Echo (QLE) has become increasingly popular and useful tool, either in the context of experiments with many-body systems Usaj , quantum computation Qcomp , or simple dynamical systems Peres ; Jalabert ; Prosen ; PZ ; Tomsovic ; Beenakker . The corresponding Classical Loschmidt Echo (CLE) has been defined PZ ; BC1 ; Eckhardt ; Veble , as

$$F(t)=\mathrm{d}^N𝒙\rho (𝒙)\rho _\mathrm{E}(𝒙,t)$$
(1)
where $`\rho (𝒙)`$ is an $`L^2`$ normalized non-negative initial density in $`N=2d`$-dimensional classical phase space with coordinates $`𝒙`$, and $`\rho _\mathrm{E}(𝒙,t)`$ is its image after the echo-dynamics, i.e. a composition of a hamiltonian flow with a slightly perturbed time-reversed hamiltonian flow. It has been found that for classically chaotic systems CLE follows QLE only for a short time that scales as $`\mathrm{log}\mathrm{}`$. For longer times, CLE of chaotic systems has been shown to follow time correlation functions Veble . However, CLE is in itself an interesting quantity in classical statistical mechanics as it provides a way to quantify the old Loschmidt-Boltzmann controversy.

For shorter times, before the relaxation of the initially localized density under echo-dynamics takes place, it has been found numericallyBC1  that CLE decays exponentially with the rate given by the Lyapunov exponent. No classical mechanism for this phenomenon has been given, apart from the necessary correspondence with QLE for which a semiclassical theory of Lyapunov decay exists Jalabert . In this letter we report several surprising analytical results for the case of localized initial densities and for sufficiently weak perturbations: (i) In many-dimensional systems, a cascade of Lyapunov decays is predicted, with the exponents which are given as consecutive sums of largest few Lyapunov exponents. Hence precise conditions for previously observed simple Lyapunov decay are understood. (ii) In loxodromic case of degenerate largest Lyapunov exponent $`\lambda _1=\lambda _2`$ we find initial decay of CLE with exponent $`2\lambda _1`$. (iii) For maps under special conditions we find non-zero average drift of the echo-dynamics resulting in faster than exponential decay of CLE. The same results should apply for quantum fidelity (QLE) up to the log-time, namely until the Wigner function follows the classical density.

The propagation of classical densities in phase space is governed by the unitary Liouville evolution $`\widehat{U}_\delta (t)`$

$$\frac{\mathrm{d}}{\mathrm{d}t}\widehat{U}_\delta (t)=\widehat{}_{H_\delta (𝒙,t)}\widehat{U}_\delta (t)$$
(2)
where $`\widehat{}_{A(𝒙,t)}=\left(\mathbf{}A(𝒙,t)\right)𝐉\mathbf{},`$ $`A`$ is any observable, and

$$H_\delta (𝒙,t)=H_0(𝒙,t)+\delta V(𝒙,t),$$
(3)
is a generally time-dependent family of Hamiltonians with perturbation parameter $`\delta `$. Matrix $`𝐉`$ is the usual symplectic unit. Similarly, $`\mathrm{d}\widehat{U}_\delta ^{}(t)/\mathrm{d}t=\widehat{U}_\delta ^{}(t)\widehat{}_{H_\delta (𝒙,t)}.`$ The classical echo propagator that composes perturbed forward evolution with the unperturbed backward evolution is also unitary and is given by

$$\widehat{U}_\mathrm{E}(t)=\widehat{U}_0^{}(t)\widehat{U}_\delta (t).$$
(4)
Using eqs. (2,3,4) and writing $`\widehat{U}_\delta (t)=\widehat{U}_0(t)\widehat{U}_\mathrm{E}(t)`$ we get

$$\frac{\mathrm{d}}{\mathrm{d}t}\widehat{U}_\mathrm{E}(t)=\left\{\widehat{U}_0^{}(t)\widehat{}_{\delta V(𝒙,t)}\widehat{U}_0(t)\right\}\widehat{U}_\mathrm{E}(t).$$
(5)
The classical dynamics has the nice property that the evolution is governed by characteristics that are simply the classical phase space trajectories, so the action of the evolution operator on any phase space density is given as $`\widehat{U}_0(t)\rho =\rho \mathit{\varphi }_t^1,`$ where $`\varphi _t^1`$ denotes the backward (unperturbed) phase space flow from time $`t`$ to $`0`$. Similarly, $`\widehat{U}_0^{}(t)\rho =\rho \mathit{\varphi }_t,`$ where $`\varphi _t`$ represents the forward phase space flow from time $`0`$ to time $`t`$. Here and in the following we assume the dynamics to start at time $`0`$.

We note that echo-dynamics (4) can be treated as Liouvillian dynamics in interaction picture, since

$`\left\{\widehat{U}_0^{}(t)\widehat{}_{A(𝒙,t)}\widehat{U}_0(t)\rho \right\}(𝒙)=`$ (6)
$`=\widehat{U}_0^{}(t)\left(\mathbf{}_𝒙A(𝒙,t)\right)𝐉\mathbf{}_𝒙\rho \left(\mathit{\varphi }_t^1(𝒙)\right)=`$
$`=\left(\mathbf{}_{\mathit{\varphi }_t(𝒙)}A(\mathit{\varphi }_t(𝒙),t)\right)𝐉\mathbf{}_{\mathit{\varphi }_t(𝒙)}\rho \left(\mathit{\varphi }_t^1\left(\mathit{\varphi }_t(𝒙)\right)\right)=`$
$`=\left(\mathbf{}_𝒙A(\mathit{\varphi }_t(𝒙),t)\right)𝐉\mathbf{}_𝒙\rho (𝒙)=\left\{\widehat{}_{A(\mathit{\varphi }_t(𝒙),t)}\rho \right\}(𝒙).`$
In the last line we used the invariance of the Poisson bracket under the flow. This extends eq. (5) to form (2)

$$\frac{\mathrm{d}}{\mathrm{d}t}\widehat{U}_\mathrm{E}(t)=\widehat{}_{H_\mathrm{E}(𝒙,t)}\widehat{U}_\mathrm{E}(t)$$
(7)
where the echo Hamiltonian is given by

$$H_\mathrm{E}(𝒙,t)=\delta V(\mathit{\varphi }_t(𝒙),t).$$
(8)
The function $`H_\mathrm{E}`$ is nothing but the perturbation part $`\delta V`$ of the original Hamiltonian, which, however, is evaluated at the point that is obtained by forward propagation with the unperturbed original Hamiltonian. It is important to stress that this is not a perturbative result but an exact expression. Also, even if the original Hamiltionian system was time independent, the echo dynamics obtains an explicitly time dependent form. Trajectories of the echo-flow are given by Hamilton equations

$$\dot{𝒙}=𝐉\mathbf{}H_\mathrm{E}(𝒙,t).$$
(9)
At this point we limit our discussion only to time independent original Hamiltonians and perturbations. Slightly more general case of periodically driven systems reducible to symplectic maps shall also be discussed later.

Inserting (8) into eq. (9) yields

$$\dot{𝒙}=\delta 𝐉\mathbf{}_𝒙V(\mathit{\varphi }_t(𝒙))=\delta \mathrm{𝐉𝐌}_t^T(𝒙)(\mathbf{}V)(\mathit{\varphi }_t(𝒙)).$$
(10)
Here we have introduced the stability matrix $`𝐌_t(𝒙)`$, $`[𝐌_t(𝒙)]_{i,j}=_j[\mathit{\varphi }_t(𝒙)]_i`$. From now on we assume that the flow $`\mathit{\varphi }_t`$ is Anosov. To understand the dynamics (10) we need to explore the properties of $`𝐌_t`$. We start by writing the matrix $`𝐌_t^T(𝒙)𝐌_t(𝒙)=_je^{2\lambda _jt}d_j^2(𝒙,t)𝒗_j(𝒙,t)𝒗_j(𝒙,t)`$ expressed in terms of orthonormal eigenvectors $`𝒗_j(𝒙,t)`$ and eigenvalues $`d_j^2(𝒙,t)\mathrm{exp}(2\lambda _jt)`$. After the ergodic time $`t_\mathrm{e}`$ necessary for the echo trajectory to explore the available region of phase space, Osledec theoremOsledec  guarantees that the eigenvectors of this matrix converge to Lyapunov eigenvectors being independent of time, while $`d_j(𝒙,t)`$ grow slower than exponentially, so the leading exponential growth defines the Lyapunov exponents $`\lambda _j`$. Similarly, the matrix $`𝐌_t(𝒙)𝐌_t^T(𝒙)=_je^{2\lambda _jt}c_j^2(𝒙_t,t)𝒖_j(𝒙_t,t)𝒖_j(𝒙_t,t),`$ where $`𝒙_t=\mathit{\varphi }_t(𝒙)`$, has the same eigenvalues \[$`c_j^2(𝒙_t,t)d_j^2(𝒙,t)`$\], and its eigenvectors depend on the final point $`𝒙_t`$ only, as the matrix in question can be related to the backward evolution. The vectors $`\{𝒖_j(𝒙_t)\}`$, $`\{𝒗_j(𝒙)\}`$, constitute left, right, part, respectively, of the singular value decomposition of $`𝐌_t(𝒙)`$, so we write for $`tt_\mathrm{e}`$

$$𝐌_t(𝒙)=\underset{j=1}{\overset{N}{}}\mathrm{exp}(\lambda _jt)𝒆_j(\mathit{\varphi }_t(𝒙))𝒇_j(𝒙)$$
(11)
assuming that the limits $`𝒆_j(𝒙)=lim_t\mathrm{}c_j(𝒙,t)𝒖_j(𝒙,t)`$, $`𝒇_j(𝒙)=lim_t\mathrm{}d_j(𝒙,t)𝒗_j(𝒙,t)`$ exist. Rewriting eq. (10) by means of eq. (11) we obtain

$$\dot{𝒙}=\delta \underset{j=1}{\overset{N}{}}\mathrm{exp}(\lambda _jt)W_j(\mathit{\varphi }_t(𝒙))𝒉_j(𝒙)$$
(12)
where $`𝒉_j(𝒙)=𝐉𝒇_j(𝒙)`$, and introducing new observables

$$W_j(𝒙)=𝒆_j(𝒙)\mathbf{}V(𝒙).$$
(13)
At this point it is perhaps necessary to discuss the nature of the vector fields $`𝒆_j(𝒙)`$, $`𝒇_j(𝒙)`$ remark . While the theorem Osledec  guarantees the existence of these directions, the actual sizes of these fields as given by $`c_j`$, $`d_j`$ are not yet well understood. Numerical data suggest that these quantities do converge as $`t\mathrm{}`$ for an Anosov system.

For small perturbations the echo trajectories remain close to initial point $`𝒙(0)`$ for times large in comparison to the internal dynamics of the system ($`t_\mathrm{e}`$, Lyapunov times, decay of correlations, etc), and in this regime the echo evolution can be linearly decomposed along different independent directions $`𝒉_j(𝒙(0))`$

$$𝒙(t)=𝒙(0)+\underset{j=1}{\overset{N}{}}y_j(t)𝒉_j(𝒙(0)).$$
(14)
For longer times, the point $`𝒙(t)`$ moves away from the initial point, but the dynamics is still governed by the local unstable vectors at the evolved point. Therefore the decay of fidelity is governed by the spreading of the densities along the conjugated unstable manifolds defined by the vector field $`𝒉(𝒙)`$.

Inserting (14) into (12) we obtain for each direction $`𝒉_j`$

$$\dot{y}_j=\delta \mathrm{exp}(\lambda _jt)W_j(\mathit{\varphi }_t(𝒙)).$$
(15)
For stable directions with $`\lambda _j<0`$, clearly after a certain time the variable $`y_j`$ becomes a constant of the order $`\delta `$.

For unstable directions with $`\lambda _j>0`$, we introduce $`z_j`$ as $`y_j=\delta \mathrm{exp}(\lambda _jt)z_j`$ and rewrite the above equation as

$$\dot{z}_j+\lambda _jz_j=W_j(\mathit{\varphi }_t(𝒙)).$$
(16)
The right hand side of this equation is simply the evolution of the observable $`W_j`$ starting from a point in phase space $`𝒙=𝒙(0)`$. Due to assumed ergodicity of the flow $`\mathit{\varphi }_t`$, $`W_j(\mathit{\varphi }_t(𝒙))`$ has well defined and stationary statistical properties such as averages and correlation functions. Thus the solution $`z_j(t)`$ of the linear damped equation (16) has also stationary statistics and well defined time- and $`\delta `$-independent probability distribution $`P_j(z_j)`$. Its moments can be expressed in terms of moments and correlation functions of $`W_j`$, in particular $`\overline{W_j}=0`$ firstmoment . The analysis remains valid in a general case of explicitly time-dependent $`W_j`$ remark .

Going back to the original coordinate $`y_j`$ we obtain its distribution as $`K_j(y_j)=P_j(z_j)\mathrm{d}z_j/\mathrm{d}y_j`$, or $`K_j(y_j)=P_j\left(\mathrm{exp}(\lambda _jt)y_j/\delta \right)\mathrm{exp}(\lambda _jt)/\delta .`$ This probability distribution tells us how, on average, points within some initial (small) phase space set of characteristic diameter $`\nu `$ spread along locally well defined unstable Lyapunov direction $`j`$ and therefore represents an averaged kernel of the evolution of such densities along this direction. Starting from the initial localized density $`\rho _0`$, of small width $`\nu `$ such that the decomposition (14) does not change appreciably along $`\rho _0`$, the echo dynamics for densities solves as $`\rho _t(𝒚)=\mathrm{d}^N𝒚^{}\rho _0(𝒚^{})_jK_j(y_jy_j^{}).`$ For stable directions $`j`$ we set $`K_j(y_j)=\delta (y_j)`$, as the shift of $`y_j`$ (of order $`\delta `$) can be neglected as compared to unstable directions. This also implies that the assumption $`\delta \nu `$ is necessary in order to get any echo at all after not too short times. CLE (1) can now be written as $`F(t)=\mathrm{d}^N𝒚\rho _0(𝒚)\rho _t(𝒚)`$. As long as the width $`\nu _j`$ of $`\rho _0`$ along the unstable direction $`j`$ is much larger than the width of the kernel $`K_j`$, there is no appreciable contribution to the fidelity decay in that direction. At time

$$t_j=(1/\lambda _j)\mathrm{log}(\nu _j/(\delta \gamma _j)),$$
(17)
where $`\gamma _j`$ is a typical width of the distribution $`P_j`$, the width of the kernel is of the order of the width of the distribution along the chosen direction. After that time, the overlap between the two distributions along the chosen direction starts to decay with the same rate as the value of the kernel in the neighborhood of $`y_j=0`$, which is $`\mathrm{exp}(\lambda _jt)`$. The total overlap decays as

$$F(t)\underset{j;t_j<t}{}\mathrm{exp}\left[\lambda _j(tt_j)\right],$$
(18)
where only those unstable directions contribute to the decay for which $`t_j>t`$. As the time $`t_j`$ is shorter the higher the corresponding Lyapunov exponent $`\lambda _j`$, fidelity will initially decay with the largest Lyapunov exponent $`\lambda _1`$. In chaotic systems with more than two degrees of freedom we, however, expect to observe an increase of decay rate after the time $`t_2`$, etc. Eq. (18) provides good description for CLE as long as $`F(t)`$ does not approach the saturation value $`F_{\mathrm{}}\nu ^N`$ where the asymptotic decay of CLE is then given by leading Perron-Frobenius eigenvalue Veble .

Though the above theory has been developed for smooth flows, the generalization to ergodic symplectic maps on bounded phase space $``$ is straightforward. We adopt notation of Ref. PZ , sect. 4: $`\mathit{\varphi }\mathit{\varphi }_1`$, discrete time $`t`$ is an integer and a general small perturbation is given by a composition $`\mathit{\varphi }_\delta =\mathit{\varphi }𝒈_\delta `$ with a near-identity map $`𝒈_\delta `$ generated by a vector field $`𝒂(𝒙)`$, $`\mathrm{d}𝒈_\delta /\mathrm{d}\delta =𝒂(𝒈_\delta )`$ with initial condition $`𝒈_0(𝒙)=𝒙`$. For the perturbed map to remain symplectic we write $`𝒂(𝒙)=𝐉\mathbf{}V(𝒙)`$ for some potential $`V(𝒙)`$. We note that on compact phase space $`V(𝒙)`$ does not need to be unique and continuous, e.g. for a unit 2-torus $`V(q+1,p)=V(q,p)+\alpha ,V(q,p+1)=V(q,p)+\beta `$ where $`\alpha ,\beta `$ are arbitrary constants. Provided hyperbolic orbits with inversion do not exist, $`\eta 0`$ remark , one finds a non-vanishing drift of echo-dynamcis, $`\overline{W}_j0`$, resulting in a possible super-exponential decay of CLE if the mean $`\overline{W_j}/\lambda _j`$ of the distribution $`P_j(z_j)`$ is larger than its width $`\gamma _j`$.

In order to illustrate super-exponential versus exponential decay of CLE we consider the perturbed cat map

$`\overline{p}`$ $`=`$ $`p+q{\displaystyle \frac{\mu }{2\pi }}\mathrm{sin}(2\pi q)+\delta s(q),\mu [0,1)`$
$`\overline{q}`$ $`=`$ $`q+\overline{p}.`$ (19)
The perturbation was chosen either as: (a) $`s(q)=\mathrm{sin}(2\pi q)/(2\pi )`$, or (b) $`s(q)=1/(2\pi )`$ where the perturbation is a shift in phase space. The main difference between the two is that the case (a) corresponds to zero drift since a unique smooth potential exists ($`\alpha =\beta =0`$), while in case (b) $`𝒂=\mathrm{const}0`$ and the Lyapunov fields have a predominant direction in phase space for this system, e.g. for the unperturbed cat map ($`\mu =0`$) $`𝒆(𝒙),𝒇(𝒙),𝒉(𝒙)`$ are constant. Since the perturbed cat has no orbits with inversionremark  the corresponding phase space observable $`W_j`$ for the case (b) has a distinct nonzero average value, causing the kernel $`K_j`$ to drift exponentially in time. The difference in the qualitative nature of the two decays is shown in figure 1. In figure 2 we show the behaviour of fidelity as a function of time for the two cases, where a super-exponential decay is observed in the case of drift.

Another result, which applies only to systems with two or more unstable directions, is the occurrence of decays which are exponential but faster than Lyapunov. In the case of well separated individual Lyapunov exponents the decay is expected to go through a cascade of increasing decay rates given by (18), whereas in the loxodromic case $`\lambda _1=\lambda _2`$ the rate is $`2\lambda _1`$. We illustrate this numerically for 4D cat maps Rivas : $`𝒙^{}=𝐂𝒙(mod1),`$ $`𝒙[0,1)^4`$, and

$$𝐂_{\mathrm{d}\mathrm{h}}=\left[\begin{array}{cccc}\hfill 2& \hfill 2& \hfill 1& \hfill 0\\ \hfill 2& \hfill 3& \hfill 1& \hfill 0\\ \hfill 1& \hfill 2& \hfill 2& \hfill 1\\ \hfill 2& \hfill 2& \hfill 0& \hfill 1\end{array}\right],𝐂_{\mathrm{lox}}=\left[\begin{array}{cccc}\hfill 0& \hfill 1& \hfill 0& \hfill 0\\ \hfill 0& \hfill 1& \hfill 1& \hfill 0\\ \hfill 1& \hfill 1& \hfill 1& \hfill 1\\ \hfill 1& \hfill 1& \hfill 2& \hfill 0\end{array}\right]$$

are two examples representing the doubly-hyperbolic and loxodromic case. Matrix $`𝐂_{\mathrm{d}\mathrm{h}}`$ has the unstable eigenvalues $`5.22,2.11`$, while the large eigenvalues of $`𝐂_{\mathrm{lox}}`$ are $`1.70\mathrm{exp}(\pm \mathrm{i1}.12)`$. The perturbation for both cases was done by performing an additional mapping at each timestep $`\overline{x}_1=x_1^{}+\delta \mathrm{sin}\left(2\pi x_3\right)(mod1)`$, $`\overline{x}_{2,3,4}=x_{2,3,4}^{}`$. In figure 3 we show the two types of decay which agree with theoretical predictions.

In conclusion, we have developed a theory for short-time decay of CLE based on classical interaction picture. Our theory predicts several new phenomena, in particular a cascade of exponential decays in systems with more than one unstable direction and doubly-Lyapunov decay for the particular case of loxodromic stability. Besides being related to quantum computation, our results for systems with many degrees of freedom provide a way to understand macroscopic irreversibility in classical statistical mechanics.

We acknowledge useful discussion with T.H.Seligman and financial support by the Ministry of Education, Science and Sport of Slovenia, and in part by the U.S. ARO grant DAAD19-02-1-0086.


format: text

text:
# ITEP-LAT/2003-20 September, 2003 Self-tuning of the P-vortices Talk presented by A.V.K. at Lattice 2003, Tsukuba

## 1 INTRODUCTION

The Abelian monopoles and central vortices seem to be most promising effective degrees of freedom to explain the confinement, for review and references see, e.g., . Straightforward interpretation of the lattice data is made difficult, however, by the fact that monopoles and vortices are defined in terms of projected fields. One still can, upon detection of the vortices or of monopoles by means of projections, measure their entropy and non-Abelian action. In case of the monopoles, it turns out that there is a non-trivial behavior of the action in the ultraviolet, see  and references therein. In other words, the monopoles can be observed only as a result of cancellation between the mass and entropy, for discussion see . The point-like facet of the monopoles is manifested also in properties of the monopole clusters .

In this note we summarize the first steps made to realize a similar program of studying the anatomy of the P-vortices . Our main result is that the non-Abelian action of the P-vortices is comparable to the action of the zero-point fluctuations and is ultraviolet divergent in the limit $`a0`$:

$$S_{vortex}0.54A_{vort}a^2,$$
(1)
where $`A_{vort}`$ is the area of the vortex, the measurements refer to the $`SU(2)`$ case and the P-vortices defined in the so called indirect central gauge, for details see . Note that naively one would assume that for non-perturbative fluctuations $`S_{npert}A\mathrm{\Lambda }_{QCD}^2`$.

In Sect. 2 we present the results of the measurements while in Sect. 3 we discuss in more detail structure of the P-vortices and some implications of our measurements.

## 2 RESULTS

We have performed our calculations in pure $`SU(2)`$ lattice gauge theory for $`2.35\beta 2.6`$. The lattice spacing $`a`$ is fixed using $`\sqrt{\sigma }=440MeV`$ where $`\sigma `$ is the string tension. At each value of $`\beta `$ we have considered 20 statistically independent configurations generated on symmetric $`L^4`$ lattices. The lattice size was $`L`$=16 for $`\beta `$=2.35, $`L`$=24 for $`\beta `$= 2.4,2.45,2.5 and $`L`$=28 at $`\beta `$=2.55,2.6. The indirect maximal center gauge  was employed to define the P-vortices. For further details see .

In Fig. 1 we summarize the results of the measurements of the P-vortex density, $`\rho _{vort}N_{PV}/(6L^4a^2)`$, where $`N_{PV}`$ is the number of plaquettes occupied by the P-vortex. Alternatively, $`A_{vort}=6\rho _{vort}V_4`$ where $`V_4`$ is the volume of the lattice. We see that the vortex density is approximately a constant in physical units:

$$\rho _{vort}4(fm)^2.$$
(2)
Note that the scaling of the vortex area was observed earlier . We confirm this result on better statistics and for smaller values of $`a`$.

In Fig. 2 we summarize measurements of the non-Abelian action associated with P-vortices. In more detail, we measure the average action density, $`S_{PV}`$, on the plaquettes dual to those forming P-vortices. We subtract from it the average plaquette action, $`S_{vac}=\beta (1TrU_P/2)`$. Both $`S_{PV}`$ and $`S_{vac}`$ are measured in the lattice units. The difference $`(S_{PV}S_{vac})=S_{vortex}`$ is shown on the Fig.2 by circles and is approximately a constant in the lattice units. Note that for particular value of $`\beta =2.4`$ this difference was measured first in Ref.  and we agree with the results obtained therein.

To further probe the structure of the vortices we have measured also the average action density next to the P-vortex world sheet. Geometrically, there are in fact two various types of plaquettes neighboring the plaquettes belonging to the P-vortex, for details see . We have measured the average action for these two types of plaquettes separately (up and down triangles on the Fig. 2). In the both cases, the average action is very close to that of vacuum. Thus, we conclude that the vortex thickness is smaller than the resolution available:

$$R_{vort}0.06fm,$$
(3)
where $`0.06fm`$ is the smallest lattice spacing used in our simulations.

## 3 SELF TUNING

Thus, the vortices represent quite a unique type of the vacuum fluctuations which exhibits both the ultraviolet and infrared scales. Let us remind the reader that the most common examples demonstrate sensitivity to one of the two scales. Thus, zero-point fluctuations are controlled by the ultraviolet scale $`1/a`$. In case of instantons both the effective size and action are determined by $`\mathrm{\Lambda }_{QCD}`$. Much less common are examples of fine tuned fluctuations when one can tune a parameter to ensure that, say, inverse size is much larger than the mass. For example, for the monopole mass in case of the compact $`U(1)`$ one has

$$m_{mon}^2const(e^2e_{crit}^2)a^2,$$
(4)
where $`e`$ is the electric charge and $`e_{crit}^2`$ is fixed (for a review and references see, e.g., ). Tuning $`e^2`$ to $`e_{crit}^2`$ one ensures that $`m_{mon}^2a^2`$.

In case of a non-Abelian theory there is no parameter to tune since the coupling is running. Nevertheless, the radiative mass of the monopole is largely cancelled by the entropy, for discussion see . Thus, the underlying gauge theory ensures a kind of self tuning for the monopoles.

Now we observe a similar phenomenon for the P-vortices: their area scales in the physical units while the tension is a constant in the lattice units. Suppression of the vortices due to their action is approximately

$$\mathrm{exp}(S_{vortex})\mathrm{exp}(0.54A/a^2).$$
(5)
For the area of the vortices not to shrink this factor is to be cancelled by the entropy. P-vortices represent closed surfaces and their entropy can be estimated knowing the position of the phase transition to the vortex percolation in case of the lattice $`Z(2)`$:

$$(Entropy)_{vort}\mathrm{exp}(+0.88A/a^2),$$
(6)
(for discussion of the $`Z(2)`$ theories and further references see, e.g., ).

Apparently, the action and entropy factors do not cancel each other to the accuracy needed. In other words, the self tuned vortices cannot be described in terms of a local action density alone but are to possess further structure. And, indeed, it was observed, for a particular $`\beta `$ that the vortices are populated by the monopoles . We confirm this picture for all the values of $`\beta `$ tested. Moreover, the plaquettes shared by the monopoles and vortices have even higher action than the P-vortices on average, see Fig. 2. Accounting for the monopoles might bridge the gap between (6) and (5) although no explicit calculation is available.

Our final remark is that we defined the vortex size in terms of the distribution of the non-Abelian action. Another manifestation of the small physical size of the vortices is the linearity of the heavy quark potential associated with vortices . However, one could define thickness of the vortices it terms of their flux, for discussion see . Then the vortex seems to be not localized to the ultraviolet cut off. Indeed, the action which we observe is still considerably smaller than that of the negative plaquettes.

In conclusion let us mention that further observations on the interplay between the monopoles and vortices were presented in another talk at this conference .


text:
# Alpha-Relaxation Processes in Binary Hard-Sphere Mixtures

## I Introduction

If one compresses or cools a liquid, there appear slow dynamical processes which are referred to as structural relaxation. These processes are precursors of the liquid-glass transition. The study of these phenomena has been a very active field of research in recent years. Several new experimental techniques were introduced to measure the evolution of structural relaxation spectra within the GHz band. Molecular-dynamics simulation techniques have been improved considerably so that correlation functions of liquids in equilibrium can be obtained for time intervals covering five to six orders of magnitude. The wealth of information obtained on glassy dynamics is a challenge for the theory of amorphous condensed matter. However, consensus on the understanding of the slow dynamics in glass-forming liquids has not yet been achieved cre (2002); pis (2003).

Simple monoatomic liquids crystallize before structural relaxation dynamics is fully developed. Therefore, studies of the glassy dynamics have to be performed on simple molecular systems or suitable mixtures. Recently, for example, a four-component mixture was studied by neutron-scattering spectroscopy. This system transforms to a metallic glass at low temperatures, but it exhibits the same scenario for the evolution of structural relaxation as known for molecules Meyer (2002). The first molecular-dynamics studies of structural relaxation in an equilibrium liquid were performed for a binary mixture of particles interacting by purely repulsive potentials Bernu et al. (1985, 1987); Roux et al. (1989). A binary Lennard-Jones system has been introduced Kob and Andersen (1994), whose interaction potentials are similar to the ones proposed for the description of the glass-forming Ni-P mixture Stillinger and Weber (1982). This system has been used extensively to analyze all facets of glassy dynamics in the equilibrium liquid and also for the quenched non-equilibrium system Kob (2003); Sciortino and Tartaglia (2001). In the mentioned previous studies, mixing was merely introduced as a means of suppressing crystallization. In the present paper, we analyze the influence of mixing on the structural relaxation.

In order to identify the effect of mixing on the glassy dynamics, we have performed molecular-dynamics simulations for four binary hard-sphere mixtures differing in the size ratio of the constituents and in the composition. By increasing the total packing fraction up to $`0.605`$, the evolution of structural relaxation was detected for a time interval up to five orders of magnitude. As reported earlier Foffi et al. (2003), two scenarios for mixing effects have been identified. For a mixture with a small size disparity of the constituents, the increase of the mixing percentage of the small particles for a fixed total packing fraction leads to a slowing down of the long-time dynamics. In this case, mixing stabilizes the glass state. However, upon mixing particles with a large size disparity, the increase of the percentage of the small particles at fixed packing fraction speeds up the structural relaxation. In this case, mixing stabilizes the liquid. The present paper reports a detailed analysis of the long-time relaxation processes, traditionally referred to as alpha processes, for the two scenarios mentioned.

Glassy dynamics and a liquid-glass transition can also be observed experimentally in colloidal suspensions. In particular, one can prepare glass-forming colloids where the interaction potential is a very good approximation to a hard-sphere repulsion Pusey (1991). Informative light-scattering studies of structural relaxation for such a hard-sphere suspension have been reported van Megen and Underwood (1993). To suppress crystallization, a narrow distribution of particle sizes was chosen. Strictly, such a system is a multi-component mixture. But ignoring the small polydispersity, it can be viewed as a one-component system. In the same sense, one can consider a colloid studied by Henderson et al. Henderson et al. (1996) as an approximation for a binary hard-sphere mixture with a size ratio $`0.8`$ between the two groups of particles, and a mixture studied by Williams and van Megen Williams and van Megen (2001) as one with size ratio $`0.6`$. For the first mixture, a dramatic effect of mixing on the nucleation ratio was observed, but no effect on the glassy dynamics has been reported. For the second mixture, it was shown that the time scale for the alpha relaxation decreased upon mixing. Our simulation results Foffi et al. (2003) suggest that the cited experiments do not deal with colloid-specific features. Rather, they exemplify the two scenarios for mixing effects on structural relaxation. Therefore, the present paper provides a detailed list of quantitative predictions for correlation functions of glassy colloids, which can be measured by photon-correlation spectroscopy.

The mode-coupling theory for ideal liquid-glass transitions provides a physical explanation for the evolution of structural relaxation in simple systems and allows for the first-principle evaluation of the density-correlation functions Götze and Sjögren (1992). The results of this theory for the hard-sphere system have been used for a detailed analysis of the light-scattering data obtained for hard-sphere colloids with a small polydispersity van Megen (1995). The theory has been extended recently to a discussion of binary hard-sphere mixtures. In particular, the above mentioned two mixing scenarios had been obtained Götze and Th. Voigtmann (2003). It was possible to describe a major part of the scattering data for a mixture Williams and van Megen (2001) quantitatively by the theoretical results Th. Voigtmann (2003a). These findings provide a motivation to use our simulation results also for a detailed quantitative test of the mode-coupling theory for the alpha-relaxation process.

The paper is organized as follows. In Sec. II, the simulation details are described, representative results for the two mixtures considered are exhibited, and the mode-coupling-theory formulas for the alpha-relaxation process are listed. Then, in Sec. III, it will be explained how the alpha-relaxation processes are parameterized. The results for the parameters are compared with the corresponding ones obtained from the mode-coupling-theory findings. Section IV presents a comparison of the alpha-relaxation master functions for the density-fluctuation correlation functions of the simulation data with the corresponding theoretical results. In Sec. V, the findings are summarized.

## II Basic concepts and results

### II.1 Specification of the systems

Let $`N`$ and $`\varrho `$ denote the total number of particles and the total number density, respectively, of the binary hard-sphere mixture (HSM) to be studied. Further numbers specifying the system are the particle diameters $`d_\alpha `$, the particle masses $`m_\alpha `$, their thermal velocities $`v_\alpha `$, the partial number densities $`\varrho _\alpha `$ and number concentrations $`x_\alpha =\varrho _\alpha /\varrho =N_\alpha /N`$, as well as the partial packing fractions $`\phi _\alpha =(\pi /6)\varrho _\alpha d_\alpha ^3`$. Here, $`\alpha =\text{A}`$ and B labels the big and small particles, respectively. In the present work, the size ratio $`\delta =d_\text{B}/d_\text{A}`$, the total packing fraction $`\phi =\phi _\text{A}+\phi _\text{B}`$, and the relative packing fraction of the smaller species, $`x=\phi _\text{B}/\phi `$, will be used as convenient control parameters to characterize the thermodynamic state.

Density fluctuations of species $`\alpha `$ for wave vector $`\stackrel{}{q}`$ are constructed from the positions $`\stackrel{}{r}_k^{(\alpha )}`$, $`k=1,\mathrm{},N_\alpha `$, of particles of type $`\alpha `$: $`\varrho _\alpha (\stackrel{}{q})=_k\mathrm{exp}[i\stackrel{}{q}\stackrel{}{r}_k^{(\alpha )}]/\sqrt{N}`$. The partial structure factors $`S_{\alpha \beta }(q)=\varrho _\alpha (\stackrel{}{q})^{}\varrho _\beta (\stackrel{}{q})`$ provide the simplest statistical information on the equilibrium distribution of particles, i.e., on the structure. Here, $``$ denotes canonical averaging. The structure factors depend on the wave vector only via $`q=|\stackrel{}{q}|`$. They can be written as $`S_{\alpha \beta }(q)=x_\alpha \delta _{\alpha \beta }+\varrho x_\alpha h_{\alpha \beta }(q)x_\beta `$, with $`h_{\alpha \beta }(q)`$ denoting the Fourier transform of the pair-correlation function. The latter can be expressed in terms of the direct correlation functions $`c_{\alpha \beta }(q)`$ via the Ornstein-Zernike equation: $`h_{\alpha \beta }(q)=c_{\alpha \beta }(q)+\varrho _\gamma c_{\alpha \gamma }(q)x_\gamma h_{\gamma \beta }(q)`$ Hansen and McDonald (1986). The $`S_{\alpha \beta }(q)`$, $`h_{\alpha \beta }(q)`$, and $`c_{\alpha \beta }(q)`$ are elements of real symmetric two-by-two matrices. The discussions will be restricted to such stable and metastable states where $`S_{\alpha \beta }(q)`$ and $`c_{\alpha \beta }(q)`$ are smooth functions of $`q`$ and of the control parameters $`\phi `$, $`\delta `$, and $`x`$.

The main quantities of interest in this paper are the density correlators $`\mathrm{\Phi }_{\alpha \beta }(q,t)=\varrho _\alpha (\stackrel{}{q},t)^{}\varrho _\beta (\stackrel{}{q})`$. These correlation functions provide the simplest statistical characterization of the structural dynamics. They are real even functions of time $`t`$, and they form the elements of a symmetric two-by-two matrix. In principle, the correlators can be measured as intermediate coherent scattering functions by neutron-scattering experiments for conventional liquids, or by photon-correlation spectroscopy for colloidal suspensions. A short-time expansion yields $`\mathrm{\Phi }_{\alpha \beta }(q,t)=S_{\alpha \beta }(q)(1/2)(qv_\alpha t)^2x_\alpha \delta _{\alpha \beta }+𝒪(t^3)`$ Hansen and McDonald (1986). Within the regime of normal-liquid states, say, $`\phi <0.4`$, the short-time dynamics varies on a $`30\%`$ level upon changes of the control parameters. There are no structural relaxation phenomena apparent in the transient dynamics.

Two mixtures shall be considered in the following. A system with $`\delta =0.60`$ and $`x=0.20`$, referred to as the $`\delta =0.60`$-system, is representative for a mixture with large size disparity. In this case, $`54\%`$ of the particles are of species B. A system with $`\delta =0.83`$ and $`x=0.37`$, referred to as the $`\delta =0.83`$-system, contains $`50\%`$ of small particles. It is representative for a mixture with small size disparity. These mixtures have been used before, together with systems of smaller percentages $`x`$, in order to demonstrate the evolution of mixing anomalies with variations of $`x`$ Foffi et al. (2003). Since we are not interested in details of the short-time dynamics, the masses of the particles are chosen equal, i.e. $`v_\text{A}=v_\text{B}`$. The units of length and time are chosen such that $`d_\text{A}=1`$ and $`v_\text{A}=v_\text{B}=1`$. With these units, the natural time scale for the microscopic motion is $`t_{\text{mic}}=d_\text{A}/v_\text{A}=1`$.

### II.2 Results from molecular dynamic simulations

We perform standard molecular dynamics simulations for binary mixtures of $`N=1237`$ and $`N=700`$ hard-sphere particles with size ratios $`\delta =0.60`$ and $`\delta =0.83`$, respectively. The algorithm follows the usual event-driven scheme for the simulation of hard-sphere particles Rapaport (1997), where the trajectory of the system is propagated from one collision to the next one. To generate dense enough initial configurations without particle overlaps, we applied the same procedure as described earlier Zaccarelli et al. (2002): Starting from a random distribution of points, particles were separated growing their diameters in successive steps until the desired size was reached. From the initial configuration, each simulation proceeds by an equilibration run, followed by a production run during which positions and velocities are saved for subsequent analysis. In all cases, the equilibration time was larger than the time it takes for the particles’ average displacement to reach one diameter of the large species, $`d_\text{A}`$. Up to four independent runs per state point have been performed to reduce statistical errors. Density correlation functions $`\mathrm{\Phi }_{\alpha \beta }(q,t)`$ and static structure factors $`S_{\alpha \beta }(q)`$ have been calculated by averaging over the independent runs and over $`300`$ different wave vectors $`\stackrel{}{q}`$ of the same modulus $`q`$. The longest simulation run requested about three weeks, i.e. the largest density studied took about three months of CPU time on a fast AMD Athlon processor to be completed.

We checked that no crystallization occurred during the production runs by monitoring the time evolution of the pressure of the system, and by visual inspection of the configurations. We also evaluated the wave-vector resolved structure factor without angular averaging to make sure that no crystalline peaks have developed. Other mixture compositions than the ones presented below have been tried; for $`\delta =0.83`$ and $`x=0.276`$ as well as for $`\delta =0.60`$ and $`x=0.10`$, it was also possible to study the glassy dynamics in the liquid phase Foffi et al. (2003), despite a stronger tendency to crystallization. A system with $`\delta =0.60`$ and $`x=0.05`$ did not stay in the homogeneous liquid phase long enough for a study of structural relaxation.

Figures 1 and 2 exhibit a typical set of structure factors and the corresponding pair distribution functions. The results refer to the $`\delta =0.60`$ system, and the lines are calculated using the Percus-Yevick theory for HSM Lebowitz and Rowlinson (1964); Baxter (1970). Obviously, this approximation theory accounts for the data rather well, even though the packing fraction $`\phi =0.60`$ considered is rather large. But there are small systematic deviations of the kind known from the discussion of hard-sphere mixtures at smaller packing fractions Malijevský et al. (1997). For example, the theory overestimates the height of the first and second peak of $`S_{\text{A}\text{A}}(q)`$ by about $`10\%`$. The contact values for the radial distribution function are $`11.6`$, $`6.87`$, and $`8.76`$ for the AA, BB and AB function, respectively, while the Percus-Yevick theory yields $`8.88`$, $`6.32`$, and $`7.28`$, respectively. It will be discussed below that these discrepancies have to be acknowledged if one intends to consider the results of the mode-coupling theory (MCT) quantitatively. Even though the results of Percus-Yevick theory are well-known, a side remark on its qualitative features might be in order. Increasing the size disparity, i.e. decreasing $`\delta `$ below unity, the height of the first diffraction peak in $`S_{\text{A}\text{A}}(q)`$ decreases. Simultaneously, the wing of the peak at $`qd_\text{A}9`$ increases. Within MCT, the first trend stabilizes the liquid state, while the second trend stabilizes the glass. The first trend dominates at the glass transition for small $`\delta `$, the second one that for larger $`\delta `$ Götze and Th. Voigtmann (2003). For the radial distribution functions, these trends correspond to a reduction of the contact values and a decrease of the averaged radius of the first neighbor shell with decreasing $`\delta `$. Moreover, at larger size disparities, a splitting of the first-neighbor-shell peak into a double peak is observed, as is evident in Fig. 2.

Figures 3 and 4 exhibit representative examples for the evolution of the glassy dynamics upon compressing the $`\delta =0.60`$ mixture. The wave vectors $`qd_\text{A}=5`$ and $`qd_\text{A}=10`$ have been chosen since they exhibit the characteristic differences in mixing effects for small and large wave vectors that have been discussed before Foffi et al. (2003); Götze and Th. Voigtmann (2003). Note that the values of the partial correlation functions can be quite different; in particular the small values of the AB correlator at $`qd_\text{A}=10`$ are responsible for the worse signal to noise ratio observed there.

The correlators for $`\phi =0.40`$ are close to exponentials whose characteristic decay time is near the natural time scale $`t_{\text{mic}}=1`$ for the microscopic dynamics. They are typical for normal-liquid behavior which can be described on a $`30\%`$ accuracy level by Enskog’s theory for dense gases Boon and Yip (1980). If the packing fraction increases to $`\phi =0.45`$, the time scale for the decay of the correlators increases by about $`40\%`$, so that for $`t5`$ all correlators have decayed to below the $`5\%`$ level of their initial values. Increasing $`\phi `$ above $`0.45`$, a new relaxation pattern evolves for the dynamics outside the transient regime, say for $`t>5`$. The $`\mathrm{\Phi }_{\alpha \beta }(q,t)`$-versus-$`\mathrm{log}t`$ diagrams exhibit a two-step relaxation scenario that has repeatedly been observed before in simulation studies and experiments. First, the correlators decrease towards some plateau. The curves become flatter and the plateau lengths increase if $`\phi `$ increases. Then, the correlators decrease from the plateau to zero. The dynamics for $`tt_{\text{mic}}`$ is called structural relaxation. Our simulations document this process, which is characteristic for glass-forming liquids, for a time interval extending over nearly five orders of magnitude.

The second step of the structural relaxation, i.e., the decay below the mentioned plateau is conventionally referred to as the alpha process. The figures demonstrate that the time scale for the alpha process increases the faster with increasing $`\phi `$ the larger the packing fraction. The decay cannot be described by an exponential function, rather it is stretched over wider time intervals. Obviously, the whole structural relaxation pattern, in particular the alpha process, shows a subtle dependence on the wave vector $`q`$. It is the goal of this paper to characterize the alpha process, in particular its $`q`$ dependence, quantitatively.

### II.3 Some mode-coupling-theory equations

Within MCT, the concept of a plateau and of an alpha-relaxation process can be defined precisely in the sense of asymptotic laws describing the dynamics near an idealized liquid-to-glass transition. These laws provide a motivation for the parameterization of the data. In addition, our data shall be used to test quantitatively the results of the theory. In this section, the required formulas are compiled.

Let us introduce an obvious matrix notation to get the following equations in a transparent form. $`𝑺(q)`$, $`𝚽(q,t)`$, etc., shall denote two-by-two matrices with elements $`S_{\alpha \beta }(q)`$, $`\mathrm{\Phi }_{\alpha \beta }(q,t)`$, etc. The Zwanzig-Mori formalism Hansen and McDonald (1986) can be used to derive the exact equation of motion for the density correlators,

$$𝑰(q)_t^2𝚽(q,t)+𝑺^1(q)𝚽(q,t)+_0^t𝑴(q,tt^{})_t^{}𝚽(q,t^{})dt^{}=\mathrm{𝟎}.$$
(1a)
Here, $`𝑰(q)`$ is a matrix of inertia parameters, $`I_{\alpha \beta }(q)=1/(q^2v_\alpha ^2x_\alpha )\delta _{\alpha \beta }`$. The kernels $`M_{\alpha \beta }(q,t)`$ are fluctuating-force correlation functions, and they reflect the complicated many-body interaction effects. The equation of motion has to be solved with the initial conditions $`𝚽(q,t=0)=𝑺(q)`$ and $`_t𝚽(q,t=0)=\mathrm{𝟎}`$. The essential step in the theory is the application of Kawasaki’s factorization approximation in oder to identify the kernel contribution $`𝑴^{\text{MCT}}(q,t)`$, which expresses the coupling of the forces to the density fluctuations. The remainder of the kernel, $`𝑴^{\text{reg}}(q,t)`$, is assumed to describe the normal-liquid-state dynamics. It is anticipated to vary regularly in the control parameters and to decay on the scale $`t_{\text{mic}}`$ for the transient motion. One gets $`𝑴(q,t)=𝑴^{\text{reg}}(q,t)+𝑴^{\text{MCT}}(q,t)`$, where the mode-coupling kernel is given by the mode-coupling functional $``$,
$$𝑴^{\text{MCT}}(q,t)=[𝚽(t)](q),$$
(1b)
which reads
$$_{\alpha \beta }[𝚽](q)=\frac{1}{2q^2}\frac{\varrho }{x_\alpha x_\beta }\underset{\alpha ^{}\beta ^{}\alpha ^{\prime \prime }\beta ^{\prime \prime }}{}\underset{\stackrel{}{k}}{}V_{\alpha \alpha ^{}\alpha ^{\prime \prime }}(\stackrel{}{q},\stackrel{}{k},\stackrel{}{p})\mathrm{\Phi }_{\alpha ^{}\beta ^{}}(k)\mathrm{\Phi }_{\alpha ^{\prime \prime }\beta ^{\prime \prime }}(p)V_{\beta \beta ^{}\beta ^{\prime \prime }}(\stackrel{}{q},\stackrel{}{k},\stackrel{}{p}).$$
(1c)
In the sum, the abbreviation $`\stackrel{}{p}=\stackrel{}{q}\stackrel{}{k}`$ is used. The vertices $`V_{\alpha \alpha ^{}\alpha ^{\prime \prime }}(\stackrel{}{q},\stackrel{}{k},\stackrel{}{p})`$ are expressed in terms of the direct correlation functions $`c_{\alpha \beta }(k)`$ and a triple average $`\varrho _\alpha (\stackrel{}{q})^{}\varrho _\beta (\stackrel{}{k})\varrho _\gamma (\stackrel{}{p})`$. Simplifying the latter by the convolution approximation, one gets
$$V_{\alpha \alpha ^{}\alpha ^{\prime \prime }}(\stackrel{}{q},\stackrel{}{k},\stackrel{}{p})=(\stackrel{}{q}\stackrel{}{k}/q)c_{\alpha \alpha ^{}}(k)\delta _{\alpha \alpha ^{\prime \prime }}+(\stackrel{}{q}\stackrel{}{p}/q)c_{\alpha \alpha ^{\prime \prime }}(p)\delta _{\alpha \alpha ^{}}.$$
(1d)
Note that the mode-coupling functional is determined by the equilibrium structure alone. Specifying the structure factors and a model for $`𝑴^{\text{reg}}(q,t)`$, the preceding equations are closed. To proceed towards a numerical solution, one introduces a grid of equally spaced wave numbers extending up to a cutoff $`q^{}`$. In this paper, we use $`140`$ wave numbers and $`qd_\text{A}^{}=56.8`$ for the calculations based on the simulated structure factors. Some results based on the Percus-Yevick approximation shall also be shown, and for those we used $`200`$ wave numbers up to $`qd_\text{A}^{}=79.8`$. We refer to Ref. Götze and Th. Voigtmann (2003) and the papers quoted there for further details.
The MCT equations exhibit bifurcations for the long-time limits of the solutions. For packing fractions $`\phi `$ below some critical value $`\phi ^c=\phi ^c(\delta ,x)`$, one gets $`𝚽(q,t\mathrm{})=\mathrm{𝟎}`$. In this parameter regime, the solutions describe ergodic liquid dynamics. For $`\phi \phi ^c`$, the long-time limits are non-degenerate symmetric positive definite matrices $`𝚽(q,t\mathrm{})=𝑭(q)`$. For these states, the solutions describe amorphous solids, i.e., ideal glasses. The long-time limits obey the implicit equations

$$𝑭(q)=𝑺(q)\left[𝑺(q)^1+[𝑭](q)\right]^1.$$
(2)
These equations are defined by the equilibrium structure alone; neither the inertia matrix $`𝑰(q)`$ nor the regular memory kernel $`𝑴^{\text{reg}}(q,t)`$ enter. The above equation for $`𝑭(q)`$ can be solved by a standard iteration procedure Franosch and Th. Voigtmann (2002).

Let $`𝑭^c(q)`$ denote the non-degenerate positive definite matrix of long-time limits at the transition point $`\phi =\phi ^c`$. For reasons of continuity, $`𝚽(q,t)`$ has to be close to $`𝑭^c(q)`$ for a large time interval if $`|\phi \phi ^c|`$ is small. The correlators are the closer to $`𝑭^c(q)`$ the smaller $`|\phi \phi ^c|`$, and the time interval of this close approach extends simultaneously. Thereby, the evolution of the plateaus, which were discussed above in connection with Figs. 3 and 4, is explained by MCT, and the $`𝑭^c(q)`$ are the MCT expressions for the plateaus.

The decay of the correlators from the plateaus to zero for small negative $`\phi \phi ^c`$ shall be characterized by some time scale $`\tau (\phi )`$. Obviously, $`lim_{\phi \phi ^c}\tau (\phi )=\mathrm{}`$. Let us consider the dynamics of the relaxation from the plateau on the time scale $`\tau `$ by writing $`t=\stackrel{~}{t}\tau `$ with $`\stackrel{~}{t}`$ fixed but positive. There holds

$$\underset{\phi \phi ^c}{lim}𝚽(q,\stackrel{~}{t}\tau )=𝚿(q,\stackrel{~}{t}),$$
(3)
where $`𝚿(q,\stackrel{~}{t})`$ obeys the equation Fuchs and Latz (1993)

$$𝚿(q,\stackrel{~}{t})=𝑺(q)\stackrel{~}{𝑴}(q,\stackrel{~}{t})𝑺(q)𝑺(q)\frac{d}{d\stackrel{~}{t}}_0^{\stackrel{~}{t}}\stackrel{~}{𝑴}(q,\stackrel{~}{t}\stackrel{~}{t}^{})𝚿(q,\stackrel{~}{t}^{})𝑑\stackrel{~}{t}^{},$$
(4a)
to be solved with the initial condition $`𝚿(q,\stackrel{~}{t}=0)=𝑭^c(q)`$. Here, $`\stackrel{~}{𝑴}(q,\stackrel{~}{t})`$ is determined by the mode-coupling functional for the critical point,
$$\stackrel{~}{𝑴}(q,\stackrel{~}{t})=^c[𝚿(\stackrel{~}{t})](q).$$
(4b)
The numerical solution of the equation for $`𝚿(q,\stackrel{~}{t})`$ is done similar to that for the full equations of motion.
Equation (3) implies the following conclusion. Given some $`\stackrel{~}{t}_{}`$ and some error margin,

$$𝚽(q,t)=𝚿(q,t/\tau )$$
(5)
is valid within the margin for $`\stackrel{~}{t}=t/\tau >\stackrel{~}{t}_{}`$, provided $`\phi \phi ^c`$ is small enough. This is the superposition principle for the MCT alpha process. It describes the correlators in terms of a $`\phi `$-independent master function $`𝚿(q,\stackrel{~}{t})`$, and attributes the strong $`\phi `$-dependence to that of the scale $`\tau `$. Presenting the correlators as functions of $`\stackrel{~}{t}`$, the interval for $`\stackrel{~}{t}`$ where they coincide expands to lower values of $`\stackrel{~}{t}`$ if $`|\phi \phi ^c|`$ decreases. The master functions $`𝚿(q,\stackrel{~}{t})`$ depend only on the equilibrium structure. Neither the inertia parameters in $`𝑰(q)`$ nor the regular kernel $`𝑴^{\text{reg}}(q,t)`$ have any influence on $`𝚿(q,\stackrel{~}{t})`$. These quantities enter the time scale $`\tau `$ only.

There are complicated but straight-forward formulas to evaluate from $`^c`$ the so-called von Schweidler exponent $`b`$, $`0<b<1`$, a critical amplitude $`𝑯^{(1)}(q)`$, which is a positive definite matrix, and a correction amplitude $`𝑯^{(2)}(q)`$ Th. Voigtmann (2003b). These quantities determine the von Schweidler expansion of the master functions,

$$𝚿(q,\stackrel{~}{t})=𝑭^c(q)𝑯^{(1)}(q)\stackrel{~}{t}^b+𝑯^{(2)}(q)\stackrel{~}{t}^{2b}.$$
(6)
Here, terms of order $`\stackrel{~}{t}^{3b}`$ have been dropped. Thereby an equation is obtained for the beginning of the alpha process.

The MCT equations have been studied before for binary HSM using the Percus-Yevick approximation for the structure factors Götze and Th. Voigtmann (2003). In the present paper, results will be presented using the structure factors $`𝑺(q)`$ obtained from the simulation work. For the two mixtures we have calculated

$`\phi _{\text{MCT}}^c=0.548,b=0.44`$ ($`\delta =0.60`$, $`x=0.20`$); (7a)
$`\phi _{\text{MCT}}^c=0.545,b=0.43`$ ($`\delta =0.83`$, $`x=0.37`$). (7b)
From the simulation data, the critical packing fractions $`\phi ^c`$ for the liquid-glass transitions of the two mixtures have been determined from the alpha-relaxation times of the density auto-correlation functions, see Sec. III.5, yielding values of $`0.606`$ and $`0.586`$, respectively. The errors of $`10\%`$ and $`7\%`$, respectively, exhibited by the values noted in Eqs. (7a) and (7b), indicate the uncertainty one should expect for MCT results. It is worth to stress that if one bases MCT on Percus-Yevick approximation for the structure factors, one gets as critical values for the two mixtures $`0.520`$ and $`0.515`$, respectively. Hence, the use of correct instead of approximated structure factor input to the theory improves the results for the critical points. It is remarkable that the modest errors of the Percus-Yevick theory, which are exhibited in Fig. 1, lead to noteworthy changes in the MCT results for the critical points.

## III Parameterization of the alpha-relaxation processes

### III.1 Evolution of the alpha process

The evolution of the alpha-relaxation scaling law is examined in Fig. 5. The upper panel is a typical example for the majority of correlators obtained in our simulations. For every packing fraction $`\phi `$, some time scale $`\tau \tau (\phi )`$ can be defined so that the long-time parts of the correlators coincide if these are considered as functions of $`\stackrel{~}{t}=t/\tau `$. This coinciding part provides the master functions $`\mathrm{\Psi }_{\alpha \beta }(q,\stackrel{~}{t})`$ for the alpha process of the fluctuation considered. Upon increasing $`\phi `$, the $`\stackrel{~}{t}`$ interval where Eq. (3) holds, expands to lower values of the rescaled time $`\stackrel{~}{t}`$. Thus, the observed scenario confirms the MCT prediction. However, our data also exhibit violations of the above described scenario, which cannot be understood in the framework of MCT. These occur only for the BB correlators of the $`\delta =0.60`$ mixture for wave vectors around the structure-factor peak position, $`qd_\text{A}7`$; and this only for the two largest packing fractions examined, $`\phi =0.600`$ and $`0.605`$. The lower panel of Fig. 5 shows a representative example. The following analysis of $`𝚿(q,\stackrel{~}{t})`$ shall therefore be based on those data sets that do not exhibit the described phenomenon; i.e. for $`\phi <0.60`$.

The stretched exponential,

$$\varphi _\alpha ^\text{K}(q,\stackrel{~}{t})=A_\alpha (q)\mathrm{exp}\left[(\stackrel{~}{t}/\stackrel{~}{\tau }_\alpha (q))^{\beta _\alpha (q)}\right],$$
(8)
is an often used empirical function for the description of alpha processes. It was introduced by Kohlrausch for the description of dielectric relaxation data. The description of the alpha-process master function by an amplitude—also called plateau value—$`A_\alpha (q)`$, a time scale $`\stackrel{~}{\tau }_\alpha (q)`$, and a Kohlrausch exponent $`\beta _\alpha (q)`$ shall be used here as well. The dashed lines in Fig. 6 exhibit representative examples for an analysis of the normalized auto-correlation functions $`\varphi _\text{A}(q,t)=\mathrm{\Phi }_{\text{A}\text{A}}(q,t)/S_{\text{A}\text{A}}(q)`$. The figures contain rescaled data for $`\phi =0.595`$ and $`0.590`$, in oder to identify a major part for the interval of rescaled times $`\stackrel{~}{t}`$, for which the superposition principle is valid. Similarly, Fig. 7 exhibits examples for a fit of the stretched exponentials to the numerical solutions of the MCT equations. The $`q`$ vectors are as close to the ones of Fig. 6 as permitted by the use of discrete wave vector grids. Note that the choice of the overall time scale is irrelevant for the discussion of the master functions $`𝚿(q,\stackrel{~}{t})`$.

Figures 6 and 7 confirm an observation often made earlier: the fits by Eq. (8) provide a good description of a major part of the alpha process. However, there are also systematic deviations between the fit function $`\varphi _\alpha ^\text{K}(q,\stackrel{~}{t})`$ and the master functions $`\mathrm{\Psi }_{\alpha \alpha }(q,\stackrel{~}{t})`$. This holds in a similar manner for the fits to the data and to the MCT results. The fit contains unavoidable systematic errors, because the fit parameters depend somewhat on the time interval chosen for the fit optimization. In our analysis, the fit was done so that the large-$`\stackrel{~}{t}`$ part is described best. Thereby, the errors of the fit appear solely for the small-$`\stackrel{~}{t}`$ part of the master functions.

Equation (6) suggests another fit formula, which is valid for the small-$`\stackrel{~}{t}`$ part of the master functions. But the small-$`\stackrel{~}{t}`$ part can be identified only to the limited extend to which the scaling interval can be established. Consequently, also the fit using the von Schweidler series contains unavoidable uncertainties. The dotted lines in Fig. 6 exhibit representative examples for such fits. Fits could be achieved with von Schweidler exponents $`b`$ chosen between $`0.40`$ and $`0.50`$. Therefore, the predicted exponents, Eq. (7), are confirmed within an uncertainty of $`\pm 0.05`$. All the results shown are obtained with the cited theoretical exponents. The remaining fit parameters shall be discussed below.

Note that these results can depend somewhat on the time window chosen for the fit, which is $`3.5<t<200`$ for the fits discussed here. Similarly, the dotted lines in Fig. 7 exhibit the results of Eq. (6) for the MCT results. But here, the functions $`𝑭^c(q)`$, $`𝑯^{(1)}(q)`$, and $`𝑯^{(2)}(q)`$, as well as the scale $`\tau `$ are calculated from the MCT equations. In this sense, the dotted lines are not fit results. The discrepancies between the dotted and the full lines represent the ones between the full solution of the MCT equations of motion and the specified second-order asymptotic description of the solution. It is reassuring that the discrepancies between the full and the dotted lines in Fig. 6 exhibit similar trends as the ones shown in Fig. 7. A quantitative account of the differences between the results shown in the two figures is included in the discussions of the following sections.

### III.2 The plateau values

The circles in Figs. 8 and 9 exhibit the plateau values of the two-step relaxation process, $`A_\alpha (q)`$, obtained by fitting Eq. (8) to the simulation data for the normalized auto-correlation functions $`\varphi _\alpha (q,t)`$. For $`qd_\text{A}<6`$, the $`A_\alpha (q)`$ are very large; they almost reach their upper limit unity for $`q`$ tending to zero. This is a typical mixing phenomenon. For one-component systems, $`A(q)0.4`$ is expected for small $`q`$ Bengtzelius et al. (1984); Fuchs and Latz (1993). The width $`q_{1/2}`$ of the $`A_\alpha (q)`$-versus-$`q`$ curves, defined by $`A_\alpha (q_{1/2})=1/2`$, is about $`7\%`$ larger for $`\alpha =\text{A}`$ in Fig. 8 than in Fig. 9. This means that the large particles are better localized in the mixture with the larger size disparity. For the smaller B particles, the opposite trend is observed. Similarly, the plateaus for the A correlators exhibit some small peak near the position $`qd_\text{A}=7.5`$ of the structure-factor peak, while the B-correlator plateau values exhibit a pronounced minimum there. For the $`\delta =0.60`$ mixture, $`A_\text{A}(q)`$ exhibits a minimum near $`qd_\text{A}=11`$, which is accompanied by a maximum near $`qd_\text{A}=12`$. Instead, the $`\delta =0.83`$ system has a shoulder for $`A_\text{A}(q)`$ for $`9<q<12`$. All these details are reproduced semi-quantitatively by the results obtained from the MCT values, which are shown as squares in the figures. There are only a few cases where the plateau values deduced from the data differ by up to $`10\%`$ from the ones deduced from the MCT results: the B plateaus for the $`\delta =0.60`$ system for $`qd_\text{A}16`$ or the A plateaus for the $`\delta =0.83`$ system for $`qd_\text{A}17`$, for example.

The symbols in Figs. 10 and 11 show the results for the plateau values $`f_{\alpha \beta }(q)`$ obtained from fitting Eq. (6) to the simulation data, normalized by $`\mathrm{\Phi }_{\alpha \beta }(q,t)/\sqrt{S_{\alpha \alpha }(q)S_{\beta \beta }(q)}`$. The fit results obtained for the $`\delta =0.60`$ mixture for $`\phi =0.590`$ and those for $`\phi =0.595`$ differ by about $`5\%`$. This difference thus appears as the inherent uncertainty of the data analysis. No such difference could be identified for the $`\delta =0.83`$ mixture. The full lines in the figures exhibit the MCT plateau values $`f_{\alpha \beta }^c(q)=F_{\alpha \beta }^c(q)/\sqrt{S_{\alpha \alpha }(q)S_{\beta \beta }(q)}`$. For the diagonal functions, $`\alpha \beta =\text{A}\text{A}`$ and $`\alpha \beta =\text{B}\text{B}`$, these lines are identical to the full lines shown in Figs. 8 and 9. Those lines should help to compare the plateau results $`A_\alpha (q)`$ with the corresponding fit results for the $`f_{\alpha \alpha }(q)`$. Obviously, all qualitative features of the plateau fits based on the Kohlrausch function agree with the ones based on the von Schweidler series, both for the fits to the simulation data and for those to the MCT results. The MCT results for $`f_{\text{A}\text{A}}^c(q)`$ of the $`\delta =0.83`$ mixture are in perfect agreement with the simulation data, while $`f_{\text{B}\text{B}}(q)`$ is underestimated systematically by MCT. But the difference is only about $`5\%`$, except for $`qd_\text{A}14`$, where the discrepancy reaches about $`10\%`$. The deviations for $`f_{\text{A}\text{B}}(q)`$ are of similar size. For the system with larger size disparity, the discrepancies between data and MCT result is somewhat larger, but it is not seriously larger than the inherent uncertainty of the fits.

Figure 7 exhibits also the leading term of the von Schweidler expansion. In general, accounting for the next-to-leading term of $`𝒪(\stackrel{~}{t}^{2b})`$ increases the range of validity of the von Schweidler expansion dramatically Fuchs and Latz (1993). Indeed, a data analysis with a $`q`$-independent exponent $`b`$ is possible only if the $`𝒪(\stackrel{~}{t}^{2b})`$ term is included Sciortino et al. (1996). However, if $`\phi `$ is not close enough to $`\phi ^c`$, it may happen that the $`𝒪(\stackrel{~}{t}^{3b})`$ terms cancel $`𝒪(\stackrel{~}{t}^{2b})`$ contributions. In such case, the fit range may shrink upon inclusion of the $`𝒪(\stackrel{~}{t}^{2b})`$ terms Franosch et al. (1997). This accident is demonstrated by Fig. 7 for $`qd_\text{A}=10.2`$. Such phenomenon cannot be foreseen in an unbiased data analysis, which then, necessarily, must lead to errors in the fit amplitudes. This explains why the sign of the correction amplitude identified in the lower panel of Fig. 6 differs from the one shown in the lower panel of Fig. 7.

An obvious source of errors in the MCT results is due to using incorrect equilibrium-structure information in the mode-coupling functional. It was mentioned in connection with Eq. (7) that replacing the structure factors by their Percus-Yevick approximations increases remarkably the difference between the MCT results for the critical packing fractions and the results derived from the simulation data. The dashed lines in Figs. 8 and 9 exemplify the same phenomenon for the plateau values. These are the MCT results for the $`f_{\alpha \alpha }^c(q)`$ based on the use of the Percus-Yevick structure factors. This approximation for the equilibrium structure leads to underestimations of the plateau values by more than $`10\%`$. Figures 8 and 9 demonstrate that MCT is so sensitive to the small deviations between $`S_{\alpha \beta }(q)`$ and the Percus-Yevick results which are exhibited in Fig. 1, that they lead to serious flaws in the quantitative predictions for the dynamics.

The difference between the MCT results for the $`f_{\alpha \beta }^c(q)`$ based on the Percus-Yevick approximation and those using $`S_{\alpha \beta }(q)`$ as obtained from numerical simulations is of similar size for both the $`\delta =0.60`$ and the $`\delta =0.83`$ system. It is not a specialty of mixtures, as is demonstrated by Fig. 12. There, the comparison is repeated for the one-component hard-sphere system (HSS). The MCT result for the critical point based on the Percus-Yevick structure factor is $`\phi _{\text{PY}}^c=0.516`$, and the corresponding plateau values $`f_{\text{PY}}^c(q)`$ are shown as a dashed line. The Verlet-Weis modification of the Percus-Yevick structure factor Verlet and Weis (1972) is an empirical improvement of $`S(q)`$, in particular for the contact values. Using this as input for the MCT, one gets as the critical point $`\phi _{\text{VW}}^c=0.525`$. The corresponding $`f_{\text{VW}}^c(q)`$, which are shown as the dash-dotted line, are systematically larger than the ones based on the Percus-Yevick results for $`S(q)`$. We could obtain simulation data for the structure factor of the metastable HSS for packing fractions up to $`\phi =0.54`$. Beyond this packing fraction, crystallization was always taking place before particles could diffuse one nearest neighbor distance, making it impossible to generate data meeting our equilibration criteria. Extrapolating the smoothly $`\phi `$-dependent $`S(q)`$ up to $`\phi =0.55`$, we calculated as the critical point $`\phi _{\text{HSS}}^c=0.546`$ and obtained the full line as the plateau values. It is remarkable that the differences between $`f^c(q)`$ and $`f_{\text{VW}}^c(q)`$ for $`qd_\text{A}>10`$ are larger than those between $`f_{\text{VW}}^c(q)`$ and $`f_{\text{PY}}^c(q)`$. Note that we find, using the simulation results for the structure factors, $`\phi ^c(\delta =0.83)<\phi _{\text{HSS}}^c<\phi ^c(\delta =0.60)`$, i.e., the system with small size disparity shows a change of $`\phi ^c`$ upon mixing that is qualitatively different from the one seen for large size disparity. This MCT result is qualitatively the same as predicted originally on the basis of the Percus-Yevick approximation Götze and Th. Voigtmann (2003) and is confirmed by our simulations Foffi et al. (2003).

### III.3 Von Schweidler-expansion amplitudes

The amplitudes $`H_{\alpha \beta }^{(1)}(q)`$ in Eq. (6) are the most important parameters quantifying the dynamics in a time interval where the correlators are close to their plateaus. The upper panel of Fig. 13 exhibits a set of representative results. It compares the amplitudes $`h_\alpha ^{(1)}(q)`$ of the normalized master functions $`\mathrm{\Psi }_{\alpha \alpha }(q,\stackrel{~}{t})/S_{\alpha \alpha }(q)`$ obtained from the analysis of the simulation data for the $`\delta =0.60`$ mixture with the corresponding quantities calculated within MCT. For the quantitative comparison, a scale factor $`S`$ has to be adjusted since the arbitrariness of the time scale $`\tau `$ implies an arbitrariness in the prefactor $`\tau ^b`$ of $`h_\alpha ^{(1)}(q)=H_{\alpha \alpha }^{(1)}(q)/S_{\alpha \alpha }(q)`$. The $`h_\alpha ^{(1)}`$-versus-$`q`$ curves exhibit a subtle structure. While $`h_\text{A}^{(1)}(q)`$ has a maximum for $`qd_\text{A}5.7`$ followed by a minimum for $`qd_\text{A}7.4`$, i.e. for a wave vector near the structure-factor peak position, $`h_\text{B}^{(1)}(q)`$ increases monotonically to a maximum for $`qd_\text{A}7.4`$. While $`h_\text{A}^{(1)}(q)`$ increases for $`qd_\text{A}>7.4`$ monotonically to a maximum for $`qd_\text{A}11`$, $`h_\text{B}^{(1)}(q)`$ exhibits a sharp minimum for $`q`$ slightly above $`7.4`$ before it also reaches a maximum for $`qd_\text{A}11`$. For $`qd_\text{A}>11`$, $`h_\text{B}^{(1)}(q)`$ decreases monotonically, while $`h_\text{A}^{(1)}(q)`$ has a minimum for $`qd_\text{A}12.7`$ and then exhibits a broad maximum. These features are reproduced by the MCT results. The MCT results agree with the data on a $`10\%`$ level, except for the B amplitudes for $`qd_\text{A}`$ near $`10`$, where there are $`20\%`$ discrepancies.

The amplitudes $`h_\alpha ^{(2)}(q)`$, which describe corrections to von Schweidler’s law $`\mathrm{\Psi }_{\alpha \alpha }(q,\stackrel{~}{t})/S_{\alpha \alpha }(q)=f_{\alpha \alpha }^c(q)h_\alpha ^{(1)}(q)\stackrel{~}{t}^b`$, exhibit a zero at some wave vector $`q_\alpha ^{}`$. For $`q<q_\alpha ^{}`$, the amplitudes are negative, and for $`q>q_\alpha ^{}`$ they are positive. These features and also the value $`q_\text{B}^{}6`$ are reproduced by the MCT results. Notably, the MCT results for $`qd_\text{A}12`$ still share some qualitative features with the results obtained from the fit to the simulation data, e.g. the sharp peak followed by a sharp minimum in $`h_\text{B}^{(2)}(q)`$ at $`qd_\text{A}7`$, and the peak in $`h_\text{A}^{(2)}(q)`$ at $`qd_\text{A}11`$. Otherwise, one notices serious discrepancies between data and MCT results. For example, MCT predicts a particularly large range of validity of von Schweidler’s law for the density fluctuations of the large particles for a wave vector $`qd_\text{A}15`$. But the data analysis is done best for this case by using a correction amplitude near $`0.2`$. The identified discrepancies signalize the limitations in the determination of a correction amplitude for an asymptotic law from data which cannot be chosen sufficiently close to the singularity.

### III.4 The stretching exponents

The exponent $`\beta `$ in Eq. (8) provides a convenient overall measure for the alpha-relaxation stretching. It quantifies in an averaged manner deviations of the alpha-relaxation process from a Debye law, $`\varphi (t)\mathrm{exp}(\stackrel{~}{t}/\stackrel{~}{\tau })`$. The latter is the universal result for the dynamics of a variable coupled to a white-noise field. Figure 14 shows that for $`2<q<6`$, the exponent decreases considerably with increasing $`q`$. There is no difference between the fluctuations for the large and the small particles in this wave-vector interval; and the stretching is larger for the system with smaller size disparity. For the wave vectors near the structure-factor-peak position, $`qd_\text{A}7`$, the stretching of the B-fluctuations of the $`\delta =0.60`$ mixture is much bigger than the one of the A-fluctuations: $`\beta _\text{B}0.5`$ versus $`\beta _\text{A}0.75`$. There is an indication of the same phenomenon for the $`\delta =0.83`$ mixture. For larger wave vectors, $`\beta _\text{A}`$ is somewhat larger than $`\beta _\text{B}`$ for the mixture with large size disparity. For the $`\delta =0.83`$ mixture, $`\beta _\text{A}`$ equals $`\beta _\text{B}`$ for $`qd_\text{A}>8`$ within the noise of the data.

For $`qd_\text{A}<6`$, the MCT results overestimate $`\beta _\alpha (q)`$ by about $`10\%`$. The $`qd_\text{A}7`$ anomaly and also the large-$`q`$ variation for the $`\delta =0.60`$ mixture are described well by the theory. For the system with small size disparity, MCT overestimates the $`qd_\text{A}7`$ anomaly, and there is a slight trend to underestimate $`\beta _\alpha (q)`$ for large $`q`$.

For large wave vectors, the MCT alpha-process-master functions approach the Kohlrausch law. In this asymptotic regime, $`\beta _\alpha (q)=b`$ Fuchs (1994). Figure 14 illustrates that this theoretical result is consistent with the simulation data.

### III.5 The alpha-relaxation time scale

The fit of Eq. (8) to the long-time parts of the correlators yields the time scales $`\stackrel{~}{\tau }_\alpha (q)`$ for the alpha processes up to an overall scale $`\tau `$. The results for the fits allow for a comparison of the scales for the fluctuations of the large particles with those for the small ones. They also allow to discuss the $`q`$ dependence of the relaxation times. In order to compare the scales from the data with those calculated from MCT, one has to fit an overall scale factor $`s`$.

Figure 15 exhibits the results for the two mixtures. For $`qd_\text{A}6`$, $`\stackrel{~}{\tau }_\text{A}(q)=\stackrel{~}{\tau }_\text{B}(q)`$, and both scales decrease with increasing wave vector $`q`$. These features are reproduced by MCT, but the theory overestimates the time scales seriously. For wave vectors near the structure-factor-peak position, $`\stackrel{~}{\tau }_\text{A}`$ exhibits a pronounced maximum, and $`\stackrel{~}{\tau }_\text{B}`$ has a sharp minimum. The ratios of the scales is about $`5`$ and $`3`$ for the systems with large and small size disparity, respectively. This feature is reproduced well by MCT. For the $`\delta =0.60`$ mixture, $`\stackrel{~}{\tau }_\text{A}(q)`$ exhibits a maximum for $`qd_\text{A}13`$, while $`\stackrel{~}{\tau }_\text{B}(q)`$ has a shoulder there. For $`qd_\text{A}15`$, the relaxation times decrease with increasing $`q`$. The time scales for the large particles are somewhat larger than those for the small particles. These features are reproduced qualitatively by MCT, but the theory underestimates the time scales by a factor $`2`$ to $`3`$ for wave vectors above the structure-factor-peak position. For the $`\delta =0.83`$ mixture, the scale $`\stackrel{~}{\tau }_\text{A}(q)`$ exhibits a shoulder for $`10<q<13`$ in accord with MCT. The time $`\stackrel{~}{\tau }_\text{B}(q)`$ exhibits a minimum for $`qd_\text{A}12`$, while MCT shows a kink there. Again, MCT underestimates the $`\stackrel{~}{\tau }_\alpha (q)`$ for $`qd_\text{A}10`$.

Let us consider the variation of the alpha-relaxation time scales as functions of the packing fraction $`\phi `$. To this end, we have determined a time scale $`\tau ^\alpha (q,\phi )`$ for this process by arbitrarily choosing $`\varphi _{\alpha \alpha }(q,t=\tau ^\alpha (q,\phi ))=0.1`$, for those values of $`q`$ where the plateau values are still appreciably larger than $`0.1`$. Equation (3) formulates the scale coupling of the MCT results. While the time scales of two alpha processes, say $`\tau _1`$ and $`\tau _2`$, diverge for vanishing distance parameter $`\epsilon =(\phi ^c\phi )/\phi ^c`$, the ratio $`\tau _1/\tau _2`$ is a smooth function of $`\epsilon `$. For example, let $`\tau _1`$ and $`\tau _2`$ refer to the alpha processes of the correlator $`\varphi _{\text{A}\text{A}}(q,t)`$ of the $`\delta =0.60`$ system for $`qd_\text{A}=10`$ and $`qd_\text{A}=5`$, respectively. If $`\epsilon `$ increases up to about $`0.05`$, $`\tau _1/\tau _2`$ decreases linearly with $`\epsilon `$ by about $`4\%`$. The simulation results behave similarly, as is shown in the upper panel of Fig. 16. If the packing fraction increases from $`0.50`$ to $`0.605`$, the alpha-relaxation time scales of the density-fluctuation correlators increase by more than three orders of magnitude. Still, the shown three representative ratios for the scales vary by less than a factor two. Hence, the scale-coupling prediction is verified on a $`10\%`$ level for the three ratios shown in the upper panel of Fig. 16 by the open symbols. These examples are representative for density-fluctuation scales with intermediate and large wave numbers. If one of the wave numbers decreases to small values, the violation of the scale coupling becomes larger, as is demonstrated by the full symbols in Fig. 16. The diffusivity $`D`$ is proportional to the inverse of the alpha-relaxation scale $`\tau _D`$ of the mean-squared displacement. Hence, a $`\tau D`$-versus-$`\phi `$ diagram demonstrates the coupling of the scales for the processes described by $`\tau `$ and that for the diffusivity. The lower panel of Fig. 16 shows the results for the $`\delta =0.60`$ system. Here, $`\tau `$ refers to the scale for the density fluctuations of the A particles, and the $`D`$ are the simulation data for the A-particle diffusivities Foffi et al. (2003). The scale coupling holds for $`\epsilon 0.05`$. However, approaching the transition point more closely, contrary to the MCT results, the scale for the diffusivity decouples from that for the density fluctuations. The diffusivity does not decrease with decreasing $`\epsilon `$ as strongly as $`1/\tau `$. The results of our simulations for the $`\delta =0.83`$ system behave similarly. The described decoupling is in qualitative agreement with the behavior found earlier for the simulation results of a binary Lennard-Jones system Kob and Andersen (1994) and for a model for water Starr et al. (1999).

MCT predicts a power-law divergence in the asymptotic limit of vanishing $`\epsilon `$ for the common scale $`\tau `$ in Eq. (3): $`\tau \epsilon ^\gamma `$. The exponent $`\gamma `$ is determined by the von Schweidler exponent $`b`$ Götze and Sjögren (1992). The value $`b0.44`$ used throughout the preceding discussions implies $`\gamma 3.0`$. Figure 17 demonstrates this property for the MCT results for $`\tau _\text{A}(q)`$ for three representative wave numbers $`q`$ in form of a rectification diagram. In agreement with typical results for the simple HSS Franosch et al. (1997), the asymptotic description holds well for $`\epsilon `$ up to about $`0.05`$, and there appear deviations if $`\phi `$ differs from $`\phi ^c`$ by more than $`5\%`$. Analogous rectification diagrams for the simulation data are shown in Fig. 18. Linear extrapolations of the data to large $`\phi `$ for three wave numbers $`q`$ yield estimates for the critical packing fractions: $`\phi ^c(\delta =0.60)=0.606\pm 0.001`$ and $`\phi ^c(\delta =0.83)=0.585\pm 0.001`$. The rectification curves for the diffusivity relaxation times for both types of particles are included in the figure as $`D^{1/\gamma }`$-versus-$`\phi `$ plots. The data seem to follow the power-law predictions, but lead to slightly different estimates of $`\phi ^c`$, also depending on the species $`\alpha =\text{A},\text{B}`$. For the $`\delta =0.60`$ system, we get $`\phi _{D_\text{A}}^c0.609`$ and $`\phi _{D_\text{B}}^c0.619`$, while for the $`\delta =0.83`$ system, $`\phi _{D_\text{A}}^c0.588`$ and $`\phi _{D_\text{B}}^c0.590`$. The decoupling of the diffusivity scales from the ones for the density fluctuations mentioned above yields this overestimation of the critical packing fractions. Let us emphasize that the described estimations of $`\phi ^c`$ have been done with the bias of a given exponent $`\gamma =3.0`$. An unbiased three-parameter fit of the scale as a function of $`\phi `$ by the formula $`\tau =C(\phi ^c\phi )^\gamma `$ suffers from correlations between the fit parameters $`\phi ^c`$ and $`\gamma `$. Indeed, such fits to the diffusivities of the two species lead to differing exponents $`\gamma `$ Foffi et al. (2003), in disagreement with MCT. A similar result has also been found in a simulation of a binary Lennard-Jones mixture Kob and Andersen (1994). Fig. 18 however shows that on the basis of our simulation data, one can, given the restricted range of validity of the asymptotic law, not distinguish between these different $`\gamma `$ values. Note that the largest discrepancy, either in $`\gamma `$ or in $`\phi ^c`$, emerges for the B particles in the $`\delta =0.60`$ system. Unbiased fits lead to larger uncertainties for $`\phi ^c`$, since a decrease of the fit parameter $`\phi ^c`$ can partly be compensated by a decrease of the fit parameter $`\gamma `$. One could get the crossing points of the $`D^{1/\gamma ^{}}`$-versus-$`\phi `$ curves closer to that of the $`\tau ^{1/\gamma }`$-versus-$`\phi `$ curves in Fig. 18, if one would use some $`\gamma ^{}<\gamma `$. Such formulation of the decoupling phenomenon is suggested by Fig. 16, since the increase in $`D\tau `$ for $`\phi `$ increasing above $`0.58`$ can be fitted by $`(\phi ^c\phi )^x`$, $`x=\gamma ^{}\gamma <0`$.

A comparison of Figs. 17 and 18 leads to two questions, which we cannot answer. Why is the range of the distances $`\epsilon `$, where the power-law asymptote describes the alpha scales, so much larger for the simulation results than for the MCT solutions? Why do the deviations of $`\tau ^{1/\gamma }`$ for larger $`\epsilon `$ from the small-$`\epsilon `$ asymptote have a different sign for the simulation results than for the MCT solutions?

## IV The alpha-process shape functions

The preceding parameterization of the correlation functions shows that there is no universal master function for the alpha-relaxation processes. For example, the stretching of the density-fluctuation auto-correlation functions for the large particles generally differs from the one for the small particles, and it depends on the wave vector of the fluctuations. It is a challenge for a microscopic theory to describe the alpha-process shape functions for different probing variables. Figures 19 and 20 present our simulation results for the two mixtures under consideration for four representative wave vectors in comparison with the MCT curves. All results are rescaled by the plateau values and by the time scales discussed in Sec. III. The data are shown for two packing fractions in order to document the asymptotic regime of validity of the superposition principle. The MCT master function is complemented by a curve for some small $`\phi ^c\phi `$, in order to indicate the effect of preasymptotic corrections to the alpha-process asymptote.

Figure 19 shows that, typically, the decay of the correlators of the $`\delta =0.60`$ system from $`90\%`$ to $`5\%`$ of their respective plateau values is stretched on a time interval of about three decades. With two exceptions, this decay is described well by the MCT master functions for the alpha process. The first exception is the AA correlator for $`qd_\text{A}=7.5`$ for rescaled times $`t4\times 10^3`$. Here, the data decay less rapidly as the exponential long-time tail obtained by MCT. The second exception concerns the BB correlator for $`qd_\text{A}=10`$ and $`qd_\text{A}=12.5`$. Here, the data fall on the master functions only for such long times, where the correlator is below $`70\%`$ of the plateau value. But, the MCT results for $`(\phi \phi ^c)/\phi ^c0.015`$, which are added in rescaled presentation as full lines in Fig. 19, exhibit the same phenomenon. The reason for this is the large size of the critical amplitude $`h_\text{B}^{(1)}(q)`$ for these wave vectors, cf. Fig. 13. They cause particularly large preasymptotic corrections to the alpha-scaling law near the plateau. Thus, the specified exceptions are not a defect of MCT, but a confirmation of a subtlety of that theory.

The test of the alpha-process shape functions for $`\delta =0.83`$, Fig. 20, exhibits a series of problems. Deviations of the data from exponential decay for the very long rescaled times $`t5\times 10^5`$ occur for $`qd_\text{A}=5`$ and $`qd_\text{A}=7.5`$ for all three correlators. The AB correlator for $`qd_\text{A}=12.5`$ shows more stretching than the MCT solutions. Furthermore, the BB correlators for $`qd_\text{A}=10`$ and $`qd_\text{A}=12.5`$ miss the plateau. However, the latter is an obvious mistake of the data analysis, which could be eliminated by correcting the plateau value. We did not carry out this correction in order to emphasize that $`5\%`$ errors in the determination of the plateau values are almost unavoidable in an unbiased data analysis. The reduction of the scaling regime for the BB correlator for $`qd_\text{A}=10`$ and $`12.5`$ occurs for the $`\delta =0.83`$ system as discussed above for the $`\delta =0.60`$ system; and it can be explained in the same manner.

The most severe problem exhibited by Fig. 20 is the following. Even for $`\epsilon 0.02`$, in the simulation, several correlators stay within a $`5\%`$ interval around the plateau for time intervals as large as $`2.5`$ decades or more. This holds, e.g., for the AA correlator for $`qd_\text{A}=7.5`$, for the BB correlator for $`qd_\text{A}=5`$, and for the AB correlator for $`qd_\text{A}=7.5`$ and $`10`$. The full lines show that this feature cannot be explained by MCT. Even if the distance parameter $`\epsilon =(\phi ^c\phi )/\phi ^c`$ is as small as $`10^2`$, the calculated correlators cross their plateaus much steeper than exhibited by the simulation data.

## V Conclusions

Molecular-dynamics simulations have been presented for two dense binary mixtures of hard spheres. One mixture deals with a size ratio $`\delta =0.60`$ for the two particle species and the other with $`\delta =0.83`$. The first system is representative for the situation where mixing stabilizes the liquid state, and the other for the one where mixing stabilizes the glass Foffi et al. (2003). The data demonstrate the evolution of a two-step relaxation scenario with increasing packing fraction $`\phi `$, Figs. 3 and 4, which is similar to the one detected previously for other systems. In this paper, a comprehensive analysis of the second relaxation step, usually referred to as the alpha process, is presented. It deals with the decay of the correlators from some plateau to zero. The process was identified as that part of the correlators exhibiting the superposition principle predicted originally by mode-coupling theory (MCT), Fig. 5. This pattern is exhibited by all our simulation data except for the correlators for the small particles of the $`\delta =0.60`$ system for wave numbers near the structure-factor peak position, and this for the two largest densities $`\phi =0.600`$ and $`0.605`$, Fig. 5. This violation of the superposition principle might indicate that MCT ignores relaxation processes which become important close to the liquid-glass transition point in this mixture. A possibility that we cannot exclude, however, are precursors of crystallization or phase separation.

The simulation data for the alpha process have been fitted by Kohlrausch functions and by the von Schweidler expansion, Eqs. (8) and (6), respectively. These fits provide two estimates of the plateau values. Usually, these estimates agree on a $`5\%`$ accuracy level, Fig. 6. Analyzing the MCT results with Kohlrausch functions, one gets plateaus which agree with the correct values within $`5\%`$, Fig. 7. This means, that both mentioned fit procedures yield reliable estimations for the plateaus, within the indicated uncertainty level. These plateau values exhibit a remarkable structure as functions of the wave number $`q`$. There are qualitative differences in the structure of the plateau functions referring to the A particles and the B particles. And there are quantitative differences between the results for the two mixtures. All these results are described quantitatively by the MCT results except for some rare cases, where data and theory differ by up to $`10\%`$, Figs. 811.

MCT requires the structure factors as input for the equations of motion. The MCT results reported in this paper are based on the simulation results for the studied systems, Fig. 1. Replacing these structure factors by their Percus-Yevick approximation results, MCT still reproduces all qualitative features of the mentioned plateau functions. However, one systematically underestimates the data. The error of the MCT results caused by the specified erroneous structure information can be as large as $`20\%`$, Figs. 8 and 9. These findings do not deal with mixture-specific effects. They apply also for the simple hard-sphere system, Fig. 12.

The stretching of the alpha processes is parameterized by the Kohlrausch exponents $`\beta _\alpha (q)`$, and the relative time scales are quantified by the scales $`\tau _\alpha (q)`$. These quantities vary with the wave vector and they depend on the index $`\alpha `$ for the species, $`\alpha =\text{A},\text{B}`$. MCT reproduces the $`\beta _\alpha (q)`$ reasonably, but there occur discrepancies up to $`20\%`$. The trends for the $`\tau _\alpha (q)`$ are reproduced by the theory, but there occur large quantitative errors, Figs. 14 and 15.

Writing the von Schweidler expansion, Eq. (6), for the diagonal correlators in the form $`\mathrm{\Psi }_{\alpha \alpha }(q,\stackrel{~}{t})=F_{\alpha \alpha }^c(q)\left[1[h_\alpha ^{(1)}(q)/F_{\alpha \alpha }^c(q)]\stackrel{~}{t}^b+[h_\alpha ^{(2)}(q)/F_{\alpha \alpha }^c(q)]\stackrel{~}{t}^{2b}\right]`$, one notices that $`h_\alpha ^{(1)}(q)`$ defines a relative time scale $`\stackrel{~}{\tau }_\alpha (q)=[h_\alpha ^{(1)}(q)/F_{\alpha \alpha }^c(q)]^{1/b}`$ while $`h_\alpha ^{(2)}(q)`$ specifies the shape. The amplitudes $`h_\alpha ^{(1)}(q)`$ are reproduced reasonably by MCT, but there occur errors up to $`20\%`$. MCT describes the trend of the $`q`$\- and species-dependence of the correction amplitudes $`h_\alpha ^{(2)}(q)`$, but there are large discrepancies between data and theory, Fig. 13.

The contradicting conclusions concerning the description of the alpha process parameters arrived at in the preceding two paragraphs indicate that the described problems are ones of the fitting procedures. Indeed, Figs. 19 and 20 show that MCT describes the alpha-process master functions well.

Qualitative discrepancies between MCT and simulation data concern the ratio of the time scales $`\tau (q)`$ for the alpha-relaxation processes for the density fluctuations of intermediate wave numbers $`q`$ and the time scale $`\tau _D`$ determining the strong variation of the particle diffusivity $`D`$, $`D1/\tau _D`$. If the packing fraction of the $`\delta =0.60`$ system increases from $`0.51`$ to $`0.605`$, $`\mathrm{log}\tau (q)`$ increases by $`3.5`$, but $`\mathrm{log}\tau _D`$ increases only by $`3.0`$. The ratio $`\tau (q)/\tau _D\tau (q)D`$ is practically constant if $`\phi `$ increases from $`0.51`$ to $`0.58`$, i.e. there is perfect scale coupling within this density interval. But increasing $`\phi `$ further, the ratio increases by up to a factor $`4`$, Fig. 16. MCT overestimates the trend to particle localization near the glass transition point. This overestimation of $`\tau _D`$ seems to be the reason why the calculated relaxation times $`\tau (q)`$ exceed the data for small $`q`$, as is demonstrated in Fig. 15 for $`qd_\text{A}5`$. The longest time scale for $`q0`$ is that for the collective diffusion process: $`\tau ^{\text{coll}}(q)1/(q^2D^{\text{coll}}(q))`$. The coupling of this mode to the density fluctuations causes the divergence of $`\tau (q)`$ for $`q0`$. One should expect that an underestimation of the tagged-particle diffusivity $`D`$ implies the same mistake for the collective diffusivity $`D^{\text{coll}}`$. This explains, why the decoupling of $`\tau (q)`$ and $`\tau (q^{})`$ increases if $`q^{}`$ decreases to small values, as shown by the filled symbols in Fig. 16.

The increase of the time scales $`\tau (q)`$ with increasing packing fraction is described well by the asymptotic power law for the MCT results, Fig. 18. However, the mentioned scale decoupling implies that the extrapolation to zero of the $`(1/\tau _D)^{1/\gamma }`$-versus-$`\phi `$ graphs leads to an estimation of the critical packing fraction $`\phi ^c`$, which exceeds the value obtained from the $`(1/\tau (q))^{1/\gamma }`$-versus-$`\phi `$ extrapolation by $`0.6\%`$ Foffi et al. (2003).

Finally, a feature of our simulation data should be emphasized which concerns the time regime where the correlators cross their plateaus. It deals with times larger than the ones describing the short-time transient but preceding the regime of validity of the alpha-relaxation scaling law. Within MCT, this regime is described for large densities by the $`\beta `$-relaxation scaling laws. In this respect, the MCT results for the hard-sphere mixtures behave as the ones for the simple hard-sphere system Götze and Th. Voigtmann (2003). Figure 20 shows, however, that the correlators for the $`\delta =0.83`$ system are close to the plateaus for time intervals exceeding the ones for corresponding MCT results by more than an order of magnitude. Hence, the $`\beta `$-relaxation theory cannot account for the simulation data dealing with the plateau crossing. In that respect our simulation data for the hard-sphere mixture are also qualitatively different from the ones measured for quasi-bidisperse hard-sphere colloids Williams and van Megen (2001) and from the simulation data for the binary Lennard-Jones mixture Kob (2003).

###### Acknowledgements.

W.G. and Th.V. thank their colleagues from the University of Rome for their kind hospitality during the time this work was performed. Our collaboration was supported in part by the European Community’s Human Potential Programme under contract HPRN-CT-2002-00307, DYGLAGEMEM, and the Deutsche Forschungsgemeinschaft through grant Go 154/12-2. G.F., F.S., and P.T. acknowledge support from MIUR Prin and Firb and INFM Pra-Genfdt. We thank S. Buldyrev for providing us the simulation code for the hard-sphere mixtures.


text:
# FROM UNINTEGRATED GLUON DISTRIBUTIONS TO PARTICLE PRODUCTION IN HADRONIC COLLISIONS AT HIGH ENERGIES

## 1 Introduction

The recent results from RHIC (see e.g. ) have attracted renewed interest in better understanding the dynamics of particle production, not only in nuclear collisions. Quite different approaches have been used to describe the particle spectra from the nuclear collisions . The model in Ref. with an educated guess for UGD describes surprisingly well the whole charged particle rapidity distribution by means of gluonic mechanisms only. Such a gluonic mechanism would lead to the identical production of positively and negatively charged hadrons. The recent results of the BRAHMS experiment  put into question the successful description of Ref.. In the light of this experiment, it becomes obvious that the large rapidity regions have more complicated flavour structure.

I discuss the relation between UGD’s in hadrons and the inclusive momentum distribution of particles produced in hadronic collisions. The results obtained with different UGD’s  are shown and compared.

## 2 Inclusive gluon production

At sufficiently high energy the cross section for inclusive gluon production in $`h_1+h_2g`$ can be written in terms of the UGD’s “in” both colliding hadrons

$$\frac{d\sigma }{dyd^2p_t}=\frac{16N_c}{N_c^21}\frac{1}{p_t^2}\alpha _s(\mathrm{\Omega }^2)_1(x_1,\kappa _1^2)_2(x_2,\kappa _2^2)\delta (\stackrel{}{\kappa }_1+\stackrel{}{\kappa }_2\stackrel{}{p}_t)d^2\kappa _1d^2\kappa _2.$$
(1)
Above $`_1`$ and $`_2`$ are UGD’s in hadron $`h_1`$ and $`h_2`$, respectively. The longitudinal momentum fractions are fixed by kinematics: $`x_{1/2}=\frac{p_t}{\sqrt{s}}\mathrm{exp}(\pm y)`$. The argument of the running coupling constant is taken as $`\mathrm{\Omega }^2=\mathrm{max}(\kappa _1^2,\kappa _2^2,p_t^2)`$.

Here I shall not discuss the distributions of “produced” gluons, which can be found in . Instead I shall discuss what are typical values of $`x_1`$ and $`x_2`$ in the jet (particle) production. Average value $`<x_1>`$ and $`<x_2>`$, shown in Fig.1, only weakly depend on the model of UGD. For $`y`$ 0 at the RHIC energy W = 200 GeV one tests UGD’s at $`x_g`$ = 10-3 \- 10-2. When $`|y|`$ grows one tests more and more asymmetric (in $`x_1`$ and $`x_2`$) configurations. For large $`|y|`$ either $`x_1`$ is extremely small ($`x_1<`$ 10-4) and $`x_2`$ 1 or $`x_1`$ 1 and $`x_2`$ is extremely small ($`x_2<`$ 10-4). These are regions of gluon momentum fraction where the UGD’s is rather poorly known. The approximation used in obtaining UGD’s are valid certainly only for $`x<`$ 0.1. In order to extrapolate the gluon distribution to $`x_g`$ 1 I multiply the gluon distributions from the previous section by a factor $`(1x_g)^n`$, where n = 5-7.

## 3 From gluon to particle distributions

In Ref. it was assumed, based on the concept of local parton-hadron duality, that the rapidity distribution of particles is identical to the rapidity distribution of gluons. In the present approach I follow a different approach which makes use of phenomenological fragmentation functions (FF’s). For our present exploratory study it seems sufficient to assume $`\theta _h=\theta _g`$. This is equivalent to $`\eta _h=\eta _g=y_g`$, where $`\eta _h`$ and $`\eta _g`$ are hadron and gluon pseudorapitity, respectively. Then

$$y_g=\mathrm{arsinh}\left(\frac{m_{t,h}}{p_{t,h}}\mathrm{sinh}y_h\right),$$
(2)
where the transverse mass $`m_{t,h}=\sqrt{m_h^2+p_{t,h}^2}`$. In order to introduce phenomenological FF’s one has to define a new kinematical variable. In accord with $`e^+e^{}`$ and $`ep`$ collisions I define a quantity $`z`$ by the equation $`E_h=zE_g`$. This leads to the relation

$$p_{t,g}=\frac{p_{t,h}}{z}J(m_{t,h},y_h),$$
(3)
where $`J(m_{t,h},y_h)`$ is given in Ref.. Now we can write the single particle distribution in terms of the gluon distribution as follows

$`{\displaystyle \frac{d\sigma (\eta _h,p_{t,h})}{d\eta _hd^2p_{t,h}}}={\displaystyle 𝑑y_gd^2p_{t,g}𝑑zD_{gh}(z,\mu _D^2)}`$ (4)
$`\delta (y_g\eta _h)\delta ^2\left(\stackrel{}{p}_{t,h}{\displaystyle \frac{z\stackrel{}{p}_{t,g}}{J}}\right){\displaystyle \frac{d\sigma (y_g,p_{t,g})}{dy_gd^2p_{t,g}}}.`$
In the present calculation I shall use only LO FF’s from .

Let us analyze now how the results for pseudorapidity distributions depend on the choice of the UGD. In Fig.2 I compare pseudorapidity distribution of charged pions for different models of UGD’s. In this calculation FF from  has been used.

In contrast to Ref., where the whole pseudorapidity distribution, including fragmentation regions, has been well described in an approach similar to the one presented here, in the present approach pions produced from the fragmentation of gluons in the $`ggg`$ mechanism populate only midrapidity region, leaving room for other mechanisms in the fragmentation regions. These mechanisms involve quark/antiquark degrees of freedom or leading protons among others. This strongly suggests that the agreement of the result of the $`ggg`$ approach with the PHOBOS distributions  in Ref. in the true fragmentation region is rather due to approximations made in  than due to correctness of the reaction mechanism. In principle, this can be verified experimentally at RHIC by measuring the $`\pi ^+/\pi ^{}`$ ratio in proton-proton scattering as a function of (pseudo)rapidity in possibly broad range. The BRAHMS experiment can do it even with the existing apparatus.

In Fig.3 I compare the theoretical transverse momentum distributions of charged pions obtained with different gluon distributions with the UA1 collaboration data . The best agreement is obtained with the Karzeev-Levin gluon distribution. The distribution with the GBW model is much too steep in comparison to experimental data. This is probably due to neglecting QCD evolution in .

## 4 Conclusions

I have calculated the inclusive distributions of gluons and associated charged $`\pi `$’s in the NN collisions through the $`ggg`$ mechanism in the $`k_t`$-factorization approach. The results for several UGD’s proposed recently have been compared. The results, especially $`p_{t,h}`$ distributions, obtained with different models of UGD’s differ considerably.

Contrary to a recent claim in Ref., we have found that the gluonic mechanism discussed does not describe the inclusive spectra of charged particles in the fragmentation region, i.e. in the region of large (pseudo)rapidities for any UGD from the literature. Clearly the gluonic mechanism is not the only one.

Since the mechanism considered is not complete, it is not possible at present to precisely verify different models of UGD’s. The existing UGD’s lead to the contributions which almost exhaust the strength at midrapidities and leave room for other mechanisms in the fragmentation regions. It seems that a measurement of $`p_t`$ distributions of particles at RHIC should be helpful to test better different UGD’s.