text:
# QUANTUM SUPERSTRINGS AND QUANTIZED FRACTAL SPACE TIME

## 1 Introduction

At the outset it must be pointed out that Quantum Mechanics (QM) and Quantum Field Theory (QFT) operate in the Minkowski space time - this is a differentiable manifold and space time points, infact point particles are legitimate. However this contradicts the very spirit of QM - arbitrarily small space time intervals imply arbitrarily large momenta and energies. Yet both QM and QFT have co-existed with this unhappy situation. From this point of view, two approaches are more satisfactory. One is string theory and the other is the discrete space time formulation. The former has evolved into the present theory of Quantum SuperStrings (QSS) while the latter leads to a formulation in terms of stochastic, fractal space-time. These two approaches appear very different, but as we will now show, they share several strikingly similar characteristics. This will throw light on some features of string theory which appear very puzzling, like the extra dimensions and the non-commutative geometry aspect.

## 2 String Theory

The origin of String Theory lies in two observed features. The first is the Regge trajectories, $`JM^2`$, which define the higher spin resonances, and which implies a finite size or spread in space for them.
The second feature is the dual resonance or $`st`$ channel scattering, which lead to the Venetziano model. It was realized that if the particles could be given an extension of the order of the Compton wavelength, that is treated as strings, then the puzzling $`st`$ feature could be explained.
All this leads to strings which are governed by the equation

$$\rho \ddot{y}Ty^{\prime \prime }=0$$
(1)
where the frequency $`\omega `$ is given by

$$\omega =\frac{\pi }{2}\sqrt{\frac{T}{\rho }}$$
(2)
$$T=\frac{mc^2}{l};\rho =\frac{m}{l}$$
(3)
$$\sqrt{T/\rho }=c$$
(4)
$`T`$ being the tension of the string, $`l`$ its length and $`\rho `$ the line density. The identification (3) gives (4) where $`c`$ is the velocity of light, and (1) then goes over to the usual d’Alembertian or massless Klein-Gordon equation. It is worth noting that as $`l0`$ the potential energy which is $`_0^lT\left(\frac{y}{x}\right)^2𝑑x`$ rapidly oscillates.
Further if the above string is quantized canonically, we get

$$\mathrm{\Delta }x^2l^2$$
(5)
The string effectively shows up as an infinite collection of Harmonic Oscillators . It must be mentioned that (5) and (3) both show that $`l`$ is of the order of the Compton wavelength. This has been called one of the miracles of String theory by Veniziano. Infact the minimum length $`l`$ turns out to be given by $`T/\mathrm{}^2=c/l^2`$ which from (3) and (4) is seen to give the Compton wavelength.
If the relativistic quantized string is given rotation, then we get back the equation for the Regge trajectories given above. Here we are dealing with objects of finite extension rotating with the velocity of light rather like spinning Black Holes. It must be pointed out that in Super String theory, there is an additional term $`a_0`$

$$J(2\pi T)^1M^2+a_0\mathrm{}\text{with}a_0=+1(+2)\text{for the open (closed) string}$$
(6)
In equation (6) $`a_0`$ comes from a zero point energy effect. When $`a_0=1`$ we have the usual guage bosons and when $`a_0=2`$ we have the gravitons.
It must be mentioned that in the relativistic case, even in Classical Theory it is known that there is an extension of the order of the Compton wavelength $`l`$ seen in (3) above such that within this extension there are negative energies, reminiscent of a Dirac particle with zitterbewegung, and indicative of the breakdown of the concept of space time points.
Finally, if we go over to the full theory of Quantum Superstrings, we deal with the Planck length, a non-commutative geometry and a total of ten space time dimensions, six of which are curled up in the Kaluza-Klein sense.

## 3 Quantized Fractal Space time

A starting point for QFST is the observation that the electron can be considered to be a Kerr-Newman Black Hole, which is otherwise classical, but which also describes the purely QM anomalous $`g=2`$ factor. However there is a naked singularity, in other words the radius of the horizon becomes complex and is given by

$$r_+=\frac{GM}{c^2}+ıb,b\left(\frac{G^2Q^2}{c^8}+a^2\frac{G^2M^2}{c^4}\right)^{1/2}$$
(7)
where $`a`$ is the angular momentum per unit mass.
The same circumstance is observed in the purely QM description of the electron by the Dirac equation. This time we have

$$x=(c^2p_1H^1t+a_1)+\frac{ı}{2}c\mathrm{}(\alpha _1cp_1H^1)H^1,$$
(8)
The bridge between equations (7) and (8) is the fact that the imaginary term in (8) is the zitterbewegung term which again is symptomatic of the fact that space time points are not meaningful. Infact an average over Compton scales removes the naked singularity.
Another way of looking at this is that for the Dirac or Klein-Gordan operator we have a non Hermitian position operator given by

$$\stackrel{}{X}_{op}=\stackrel{}{x}_{op}\frac{ı\mathrm{}c^2}{2}\frac{\stackrel{}{p}}{E^2}$$
(9)
The Hermitianization of (9) leads to a quantization of space time at the Compton scale. As with strings, we have in this description a background ZPF  corresponding to a collection of oscillation bounded by the Compton wavelength. Interestingly an equation like (2) holds.
Thus we have to deal with space time intervals, which are a sort of a minimum cut off. This leads to the non commutative geometry

$$[x,y]=0(l^2),[x,p_x]=ı\mathrm{}[1+l^2],[t,E]=ı\mathrm{}[1+\tau ^2]$$
(10)
. Interestingly using (10) we can get back the Dirac representation. Next, from the Dirac equation we can get the Klein-Gordon equation  or alternatively the massless string equation (1) subject to (4). (Physically this means that bosons are bound states of Fermions.) Further with bound Fermions we recover the equation of Regge trajectories . Interestingly, in Superstring theory also Dirac spinors are introduced ad hoc, averages over an internal time $`\tau `$ are taken, and then the Klein-Gordon equation of the bosonic theory is deduced. At the same time this space time fudge is exactly of the type required to explain the $`st`$ channel feature. It must be observed that in (10) if $`l^2=0`$, then we are back with usual space time and Quantum Mechanics. It is only when $`l^20`$ that we have a non commutativity, leading to Fermions. Indeed, as also noted by Witten, our usual Lorentzian space time differentiable manifold is bosonic and leads to the bosonic string or Klein-Gordon equation (1), \- but there is a substructure given by (10) which leads to the Fermionic description and quantized space time. The former is an approximation when $`l^2=0`$ and we have equation (1) and what Witten calls bosonic space time (Cf. ). The latter gives the spinorial representation for the Lorentz group which as noted by Einstein(Cf. ref.) is the more fundamental representation and also the Dirac equation, rather than (1).
It is quite significant that the above discrete space time features can be understood in terms of the Nelsonian Theory in which the diffusion constant, in the light of the above remarks has been meaningfully identified with $`\mathrm{}/m`$ . Here we are dealing with a double Weiner process, and the minimum space time cut offs are stochastic in nature and may be called QFST.
Infact, a quick way to Quantum Mechanics from the Nelsonian theory is by starting with the diffusion equation,

$$\mathrm{\Delta }x\mathrm{\Delta }x=\frac{h}{m}\mathrm{\Delta }t$$

This can be rewritten as the usual Quantum Mechanical relation,

$$m\frac{\mathrm{\Delta }x}{\mathrm{\Delta }t}\mathrm{\Delta }x=h=\mathrm{\Delta }p\mathrm{\Delta }x$$

To throw further light on this, we observe that it is well known that from a two state formulation we can recover the Schrodinger equation by considering the diffusion of a particle to neighbouring points on either side, separated by a distance of the order of the Compton wavelength. Similarly if we could consider a diffusion forwards and backwards in time through a Compton time interval reminiscent of the double Weiner process referred to, we can get back the d’Alembertian equation, or the string equation (1). To see this more explicitly, let us consider a point $`x`$ and its neighbouring points $`x\pm l`$ and a time $`t\pm \tau `$. Then a typical simple diffusion equation would be given by , following Feynman,

$$C_ı(t\tau )+C_ı(t+\tau )=\underset{j}{}[\delta _{ıj}\frac{ı}{\mathrm{}}H_{ıj}(t)\tau ]C_j(t)$$
(11)
where we denote by, $`C_ı(t)<ı|\psi (t)>,`$ the amplitude for the state $`|\psi (t)>`$ to be in the state $`|ı>,`$ and

$$<ı|U|j>U_{ıj}(t+\tau t,t)\delta _{ıj}\frac{ı}{\mathrm{}}H_{ıj}(t)\tau t.$$

and $`ı`$ denotes the point $`x`$ and the points $`x\pm l`$. In the limit $`l=c\tau 0`$ (11) gives the Klein-Gordon equation or alternatively the equation (1) subject to the condition (4). This also shows the close relation of the diffusion process and special relativity, as noted in .
Finally in the above formulation using the spirit of fluctuations and the non commutative relation (10), we get, $`\mathrm{\Delta }xl^3`$, where $`l`$ is the electron Compton wavelength. We are thus lead firstly to the Planck length and secondly to the fact that the single $`x`$ dimension has now become three dimensional. The two extra dimensions are fluctuations and curled up in the Kaluza-Klein sense as shown elsewhere. It must be stressed that in this formulation we use as seen in the Super String case, the zero point fluctuations.

## 4 Conclusion

It is not surprising that QSS and QFST reach the same extra dimensions, non commutative geometry and Planck scale, because the former starts with extension and then the QM input, while the latter starts with extension and the Nelsonian QM input. From this point of view we get an insight into the otherwise inexplicable features of QSS. As Veneziano put it , ”we thus face a kind of paradoxical situation. On the other hand quantum mechanics is essential to the success of the KK idea. At the same time, QFT gives meaningless infinities and spoils the nice semiclassical results. If the beautiful KK idea is to be saved we need a better quantum theory than QFT. I will argue below that such a theory already exists: it is called (super) string theory.” The key lies in space time intervals which underlie a Fermionic sub structure of space time given by (10), which has a perfectly natural origin in stochastic process or in the spirit of which lies ”Law without Law”, and is exemplified by QFST.


format: text

text:
# Synchronization learning of coupled chaotic maps

## I Introduction

Synchronization is a form of macroscopic evolution observed in a wide class of complex systems. Typically, it appears when the range of the interactions inside the system is of same order as the system size. Mechanical and electronic devices, as well as certain chemical reactions  are known to exhibit synchronized dynamics. In the realm of biology, among many other instances , synchronization appears at the cellular level in neural networks  and in heart tissues . Animal populations show also complex forms of synchronous behavior, the most spectacular example being probably the synchronous flashing of certain fireflies . The spontaneous occurrence of synchronization in biological systems suggests that this form of collective behavior develops as a consequence of evolutionary selection or of some kind of adaptive learning. In this paper we explore a model where a set of coupled chaotic elements is added with a feedback learning process targeting a collective synchronized state. We find that, under suitable conditions, the system evolves from completely incoherent behavior to a state of full synchronization.

Collective behavior under the effect of long-ranged interactions can be modelled by means of ensembles of globally coupled dynamical systems. Introduced by Kaneko a decade ago , globally coupled logistic maps have proven to be an appropriate paradigm for such kind of emerging evolution. The system consists of $`N`$ identical mappings whose individual dynamics, in the absence of interactions, is given by $`x(t+1)=F[x(t)]`$ with $`F(x)=rx(1x)`$. The individual dynamics are coupled according to

$$x_i(t+1)=(1ϵ)F[x_i(t)]+\frac{ϵ}{N}\underset{j=1}{\overset{N}{}}F[x_j(t)],$$
(1)
($`i=1,\mathrm{},N`$), where $`ϵ(0,1)`$ is the coupling intensity.

The system defined by Eqs. (1) features a variety of collective behaviors in the space spanned by the nonlinear parameter $`r`$ and the coupling constant $`ϵ`$ . At large values of $`ϵ`$ the ensemble is in a coherence phase, where all the elements tend asymptotically to exactly the same trajectory $`x(t)`$. For sufficiently long times, thus, the system is fully synchronized. Note that the trajectory of the synchronized ensemble is governed by the dynamics of a single element. At low $`ϵ`$, when coupling is weak, a turbulence phase is observed, where the evolution is completely unsynchronized. At intermediate values of the coupling intensity a clustering phase with different groups of mutually synchronized elements appears.

Bearing in mind the role of synchronization as a collective acquired behavior of biological systems, we incorporate to model (1) an additional evolutionary mechanism. Concretely, each element is allowed to vary its coupling constant in time, so that the effect of the collective evolution on its dynamics—defined by the last term in the right-hand side of (1)—is modified according to a given criterion. Thus, the system may be able to learn to perform a specific collective task, in particular to evolve towards a coherent synchronized state. We focus the attention in the properties of the learning process, and find that a regime of successful learning exists for a sufficiently stiff algorithm. Within this regime, the time necessary to achieve synchronization is determined by the stiffness.

## II Modeling

We consider a variation of model (1) where the coupling constants depend on time and may be different for each element, namely,

$$x_i(t+1)=[1ϵ_i(t)]F[x_i(t)]+\frac{ϵ_i(t)}{N}\underset{j=1}{\overset{N}{}}F[x_j(t)].$$
(2)
Inhomogeneous global coupling in ensembles of logistic maps with time-independent coupling constants has been considered in previous work . It can be shown that the system exhibits full synchronization if and only if all the coupling constants $`ϵ_i`$ are larger than a certain value $`ϵ_c`$. This value turns out to coincide with the critical point for the onset of full synchronization in homogeneous ensembles, and can be given in terms of the Lyapunov exponent $`\lambda `$ of a single map as $`ϵ_c=1\mathrm{exp}(\lambda )`$ . In particular, for a nonlinear parameter $`r=4`$ we have $`\lambda =\mathrm{ln}2`$ and therefore $`ϵ_c=1/2`$.

Time variation of coupling intensities as a form of adaptive behavior has been considered in globally coupled maps by analogy with synaptic evolution in neural networks , and in models of asymmetric imitative dynamics for two-element systems . Here, we consider a learning algorithm based on a comparison of the instantaneous state $`x_i(t)`$ of each element with a global property of the ensemble, namely the instantaneous average state $`x(t)=N^1_jx_j(t)`$. If the distance from the individual state to the average is larger than a certain threshold $`u`$ the learning process acts, and the coupling constant of the element is changed to a new value, chosen at random from a uniform distribution in $`(0,1)`$. Otherwise, $`ϵ_i`$ remains unchanged. Explicitly,

$$ϵ_i(t+1)=\{\begin{array}{cc}\xi _i(t+1)\hfill & \text{ if }|x_i(t+1)x(t+1)|>u\hfill \\ ϵ_i(t)\hfill & \text{ otherwise},\hfill \end{array}$$
(3)
where $`\xi _i`$ is a random number with uniform distribution in $`(0,1)`$. Note that the threshold $`u`$ can be interpreted as an inverse measure of the stiffness of learning. The evolution proceeds according to the following dynamical rules. First, the state $`x_i(t)`$ of every element is updated to $`x_i(t+1)`$ applying map (2). The average state $`x(t+1)`$ is then calculated. Finally, the learning algorithm (3) is applied to every element. This procedure is successively iterated, so that each evolution step consists of two substeps where coupled dynamics and learning act sequentially. Both processes are applied synchronously to the whole system.

This is a form of stochastic unsupervised learning  where the whole ensemble is expected to selforganize into a coherent state where the orbit of every element coincides with the average trajectory. The algorithm can be interpreted as an adaptive control mechanism . Recently, control techniques have been proposed to drive both low-dimensional and extended dynamical systems towards a prescribed state, such as a particular spatiotemporal pattern . The present variant is inspired in arguments of biological plausibility. First, the target of learning is not a specific dynamical state, but a wide class defined by a collective property, namely synchronization. This class includes not only infinitely many synchronized orbits but also a variety of different configurations of the set of couplings $`ϵ_i`$. Second, the learning algorithm acts at the individual level. That is, at each time step the criterion (3) is applied to each element. According to its individual state a modification is introduced to its coupling. The collective state achieved through learning emerges thus as a consequence of individual evolution. Finally, the modification applied to the coupling intensities is random and unbiased. No hints are given on the desired values of the coupling intensities, which must be adaptively found by the system through iterations of trials and errors. Successive values of $`ϵ_i`$ are completely uncorrelated.

In the following section, we report results of extensive numerical realizations of the above model. They correspond to an ensemble of logistic maps with $`r=4`$, i.e. at the fully developed chaotic regime. The system is investigated as a function of the threshold $`u`$. We detect a sharp transition at $`u=0.5`$, between a regime of successful learning and a regime where learning fails. The origin of this transition is identified by studying the intermittent dynamics just before the state of full synchronization is reached. We define a parameter that measures the performance of learning and find an optimal value of $`u`$ for which learning is fastest.

## III Results

The numerical results reported in this section correspond to systems with $`N=10^3`$ elements. The initial states $`x_i(0)`$ are distributed at random, with uniform density, in the interval $`(0,1)`$. All the coupling constants have initially the same value, $`ϵ_i(0)=ϵ_0=0.1`$. For this coupling intensity and a nonlinear parameter $`r=4`$, the initial state of the system is well within the turbulence phase . During a first stage, the ensemble is left to evolve $`10^3`$ steps without applying the learning dynamics, such that the individual states $`x_i`$ adopt the characteristic distribution of an incoherent state. After this, time is reset to zero, learning is switched on, and the states and coupling constants of all elements are recorded during the following $`10^4`$ to $`10^6`$ steps.

A suitable way of representing the instantaneous state of the system is a plot where each element is shown as a dot in the plane spanned by the individual state $`x_i`$ and the coupling constant $`ϵ_i`$ . Figure 1 shows four such snapshots for a single system at different times, corresponding to a threshold $`u=0.2`$. At the initial time $`t=0`$ all elements have the same coupling constant $`ϵ_0`$, and form an incoherent, extended cloud in $`x_i`$. At subsequent times, $`t=20`$ and $`t=80`$, the effects of learning are clearly visible. The elements that have migrated to larger values of $`ϵ_i`$ form now a rather compact cluster though, for $`t=20`$, many elements with relatively large coupling constants are still far from the main cluster. For $`t=80`$, almost all the elements have got coupling constants above the critical value $`ϵ_c=0.5`$. The remaining elements form a small cloud just below $`ϵ_c`$, and approximately follow the motion of the main cluster. Finally, for $`t=10^5`$ all the elements have $`ϵ_i>0.5`$ and the same value of $`x_i`$. For this value of $`u`$, learning has been successful and the ensemble has become fully synchronized.

We show in Fig. 2 a situation where, on the other hand, learning is not able to lead the system to the synchronization phase. Here $`u=0.6`$ and the system has been left to evolve for $`t=10^4`$ steps. After an initial redistribution of the coupling constants, the ensemble appears to have reached a stationary state with a complex organization in the $`(x_i,ϵ_i)`$-plane (compare with analogous distributions reported in ) but with no traces of synchronization. According to our simulations, this state is preserved at longer times, of the order of $`10^6`$ steps.

These preliminary results suggest that some kind of transition occurs at an intermediate value $`u_c`$ of the threshold, separating two regions where learning is respectively successful and unsuccessful. This transition is characterized in the following.

### A Learning transition

As a global measure of the collective state of the ensemble during its evolution we have chosen the mean dispersion of the individual states $`x_i`$, namely, $`\sigma _x=\sqrt{x^2x^2}`$. For the fully synchronized state, $`\sigma _x=0`$. We have studied the evolution of $`\sigma _x(t)`$ for several values of the threshold $`u`$, and found that for $`u>0.5`$ the dispersion asymptotically approaches a finite value. The average of this asymptotic value over different realizations for a fixed threshold is practically independent of $`u`$, $`\overline{\sigma }_x0.25`$. This is the dispersion that corresponds to the state shown in Fig. 2.

For $`u<0.5`$, on the other hand, we have always found that $`\sigma _x0`$ for sufficiently long times. We stress, however, that the typical times associated with this evolution depend strongly on the threshold, as shown in detail later. In any case, for such values of $`u`$, the ensemble approaches asymptotically the state of full synchronization. Thus, the critical threshold $`u_c0.5`$ is the boundary between the zone of successful learning ($`u<u_c`$) and the zone where the system fails to learn how to synchronize ($`u>u_c`$). As discussed, the average asymptotic dispersion $`\overline{\sigma }_x`$ is sharply discontinuous at $`u_c`$.

Actually, this abrupt transition in the performance of learning can be explained by studying the dynamics of our system in the intermediate stages of evolution, when a considerable fraction of the ensemble already defines a compact cluster whereas the remaining elements, whose coupling constants are just below $`ϵ_c=0.5`$, form the small cloud depicted in Fig. 1 for $`t=80`$. In this situation, the system is in the threshold of a bifurcation where the stable state of full synchronization appears. Indeed, it would suffice to slightly change the coupling constants of the elements in the cloud to values above $`ϵ_c`$ in order to create an attractor corresponding to the synchronized state. In the threshold of such bifurcation the system is expected to display intermittent evolution . We have in fact verified in the numerical simulations that the dynamics of each element in the cloud exhibits two distinct regimes. Most of the time, the element is found in a “laminar” regime, where its state is practically equal to that of the main cluster. Occasionally, however, the evolution exhibits “turbulent” bursts during which the element performs short excursions far away from the cluster. In order to illustrate this behavior, we have recorded the evolution of a single element in the ensemble, whose coupling constant is initially fixed at $`ϵ_i=0.4997`$. Learning is subsequently prevented for this special element, in such a way that its coupling constant remains fixed as time elapses. For long times, $`t10^6`$, the coupling constants of all the other elements are found above $`ϵ_c`$. Intermittency is apparent in Fig. 3, where the difference $`\mathrm{\Delta }x=|x_ix|`$ has been plotted for the special element as a function of time. Here, the average state $`x`$ is essentially determined by the position of the main cluster.

As expected, the intervals of “laminar” behavior are found to grow in length as the coupling constant of the element under study approaches $`ϵ_c`$. Conversely, the frequency of “turbulent” bursts decreases. Note that, for the elements in the cloud, learning is possible during these bursts only. In fact, according to (3), learning acts when the difference between the individual state and the average state is large enough. Thus, the closer the coupling constant is to $`ϵ_c`$, the later an element undergoes a learning step.

The dynamics of a single element in the cloud below $`ϵ_c`$ can be well approximated as follows. We disregard the effect of the remaining elements in the cloud and suppose that the interaction of the element under study with the ensemble occurs only through the main cluster. Conversely, we suppose that the main cluster contains essentially all the elements of the system in a synchronized state $`x_0(t)`$, and that is not affected by the dynamics of the cloud. Within these assumptions, the state $`x(t)`$ of the element under study, whose coupling constant is $`ϵϵ_c`$, evolves according to

$$x(t+1)=(1ϵ)F[x(t)]+ϵF[x_0(t)].$$
(4)
Meanwhile, the state of the main cluster obeys the dynamics of a single independent map, $`x_0(t+1)=F[x_0(t)]`$. The analytical study of the intermittent evolution of $`x(t)`$ from Eq. (4) may be difficult. However, we can easily find bounds for the excursions of $`x(t)`$ during the bursts. Note in fact that we can write

$`\mathrm{\Delta }x(t+1)`$ $`=`$ $`|x(t+1)x_0(t+1)|`$ (5)
$`=`$ $`(1ϵ)|F[x(t)]x_0(t+1)|.`$ (6)
Putting $`F(x)=4x(1x)`$ and taking into account that at any time both $`x(t)`$ and $`x_0(t)`$ are in the interval $`(0,1)`$, we find $`\mathrm{\Delta }x<1ϵ`$ which, for $`ϵϵ_c`$ reduces to

$$\mathrm{\Delta }x<1/2.$$
(7)
This bound is clearly seen in Fig. 3. We conclude that the distance from the main cluster to an element in the cloud during a burst cannot be larger than $`1/2`$. For such an element, consequently, learning will occur during a sufficiently ample burst only if $`u<0.5`$. If $`u>0.5`$, the learning algorithm will never be applied and the system will not reach the fully synchronized state. This fixes the learning transition at $`u_c=0.5`$.

### B Learning times

In order to give a more detailed description of the evolution of our system under the action of learning, we focus now the attention on the time needed to approach the synchronization state. To define such a time, we study first the fraction $`n(t)`$ of elements that, at a certain moment, have their coupling constants below $`ϵ_c`$. We have $`n(0)=1`$ and, for $`u<u_c`$, $`n(t)0`$ as $`t\mathrm{}`$.

Figure 4 shows the decay of $`n(t)`$ as a function of time, for several values of the threshold $`u<u_c`$. Each curve is the average of several hundred realizations. We see that the decay is monotonous, though for some values of $`u`$ the evolution is extremely slow. This happens near $`u=0`$ and $`u=0.5`$, whereas at intermediate values of the threshold the decay is faster.

It is apparent from Fig. 4 that, for different values of the threshold, the functional form of $`n(t)`$ is not uniform. In order to define a characteristic time associated with learning, then, we fix a reference level $`n_0`$ and measure the time $`T(u)`$ needed for $`n(t)`$ to reach that level for each value of $`u`$. A plausible value for the reference level is $`n_0N^1`$, which indicates that for $`t>T(u)`$ only a few elements remain in the unsynchronized cloud. In our simulations, the learning time $`T(u)`$ has been determined taking $`n_0=N^1=10^3`$. It thus corresponds to the time taken by all elements but one to migrate to the main cluster. The learning time as a function of the threshold $`u`$ is shown in Fig. 5. Each dot stands for the average of several hundred realizations. There is a clear minimum at $`u0.01`$ showing that, as for the performance of learning, there is an optimal choice for the threshold $`u`$. Note that, as far as $`n_0N^1`$, the position of the minimum does not depend on the reference level.

The presence of an optimal threshold for our learning algorithm can be explained as follows. For large values of the threshold, $`u0.5`$, the formation of the small cloud at $`ϵ0.5`$ and the ensuing appearance of intermittent evolution occur relatively fast. Once this situation has been established, however, the elements in the cloud spend very long times in the “laminar” regime and in unsuccessful bursts, whose amplitude is not enough to drive the elements beyond the threshold. Successful bursts, where $`\mathrm{\Delta }x>u`$ \[see Eq. (6)\], become in fact increasingly rare as $`uu_c`$ and, consequently, the learning time is expected to diverge in such limit. At the other end, $`u0`$, the threshold is very narrow and the elements keep changing their coupling constants for long times. They need many attempts to approach the average behavior. Therefore, even the initial stage during which the main cluster is formed lasts asymptotically large times as $`u0`$. Since the learning time $`T(u)`$ should diverge both at $`u=0`$ and at $`u=u_c`$, it must reach (at least) a minimum for an intermediate value of the threshold, as fully confirmed by our numerical results.

## IV Summary and Conclusion

We have studied an ensemble of globally coupled chaotic maps able to change their individual couplings in order to evolve from an incoherent collective dynamics to a completely synchronized state. The learning procedure is implemented by means of a stochastic unsupervised algorithm characterized by a single parameter $`u`$ than measures the stiffness of learning. Our numerical results show that the emergence of synchronization is only possible for a specific range of the parameter $`u`$. In fact, a sharp transition at $`u_c=0.5`$ has been found, separating a regime of successful learning ($`u<u_c`$) from a regime where the algorithm fails to drive the system to full synchronization ($`u>u_c`$). In the zone of successful learning, in turn, the time needed to reach a certain level of the learning process has been shown to strongly depend on $`u`$. In particular, we have found that there is an optimal value, $`u0.01`$, for which the learning time is minimum, i.e. learning is fastest.

The sharp transition between the regimes of successful and unsuccessful learning can be explained taking into account the intermittent evolution observed just before the state of full synchronization has been reached. This kind of evolution is in fact typical at the threshold of a synchronization transition . In this regime it is possible to formulate an approximate dynamical description that yields the bounds for the amplitude of intermittent bursts, during which the elements are subject to the learning process. These bounds define in turn the maximal value $`u_c`$ for which learning is possible. As for the optimal value of $`u`$, we have argued that both for $`uu_c`$ and for $`u0`$ the learning time is expected to diverge, so that at least one minimum should be found for intermediate values. The existence of an optimal value for the learning stiffness should be a generic property of a large class of learning algorithms. This point has already been discussed to some extent in connection with several training algorithms for neural networks . Indeed, every teacher should know that there is an optimal “pressure” to be applied on the average student to obtain the best and fastest results in learning.

Several generalizations to the present model can be foreseen, attempting to describe other situations found in real systems. Our learning algorithm, in fact, is based on the comparison of the individual state of each element with a global quantity, namely, the average state over the ensemble. This could be replaced by a sort of “local” criterion, where the comparison takes place between pairs  or small groups of elements. Moreover, coupling constants could be subject to smoother changes, representing a smarter learning process, instead of the trial-and-error method used here. In the line of some training algorithms for neural networks and of optimization schemes, a variation of the present model would consist in allowing the learning stiffness to change with time. A suitably controlled temporal variation for $`u`$ could in fact result in a substantial decrease of the learning times. Learning itself could be fully replaced by an evolutionary mechanism, in the spirit of genetic algorithms. In this case, unsuccessful elements should be eliminated and replaced by slightly modified copies of successful elements, which are the effects expected from natural selection and mutation, respectively. Finally, the individual dynamics of the coupled elements admits to be varied within an ample class of behaviors, including discrete and continuous evolution. The present work, in summary, is a first step in the study of a wide spectrum of problems, of interest from the viewpoint of biology, optimization techniques, and artificial intelligence.


text:
# Nonlinear localized waves in a periodic medium

## Abstract

We analyze the existence and stability of nonlinear localized waves in a periodic medium described by the Kronig-Penney model with a nonlinear defect. We demonstrate the existence of a novel type of stable nonlinear band-gap localized states, and also reveal an important physical mechanism of the oscillatory wave instabilities associated with the band-gap resonances.

Wave propagation in periodic media is associated with many interesting physical phenomena . Modern technology allows to create different kinds of macro- and mesoscopic periodic and layered structures such as semiconductor superlattices and heterostructures, magnetic multilayers possessing the giant magnetoresistance, multiple-quantum-well structures, optical waveguide arrays, photonic band-gap materials and photonic crystal fibers, etc. The main feature of different periodic structures (which follows from the classical Floquet-Bloch theory) is the existence of forbidden frequency band gaps (or stop bands) where linear waves undergo Bragg reflection from the periodic structure . However, many of the recently fabricated periodic structures exhibit pronounced nonlinear properties that give rise to qualitatively new physical effects such as multistability of a finite nonlinear periodic medium  and energy localization in the form of gap solitons . Such effects are usually analyzed in the framework of the coupled-mode theory , and they are associated with the nonlinearity-dependent tuning of the stop band as the wave intensity is increased.

Another fundamental reason for the wave localization is the presence of defects. As a matter of fact, artificially introduced inhomogeneities are often used as a powerful means of controlling wave scattering. However, the defect-induced localization in nonlinear and periodic media is largely an open area of research. We stress that the defect response is determined by the local field amplitude, and therefore the standard averaging procedure can not be directly applied to this kind of problems.

Recent experiments, e.g. the observation of optical gap solitons  and the control of coherent matter waves in optical lattices , as well as theoretical results such as the discovery of the oscillatory instability of gap solitons , call for a systematic analysis of nonlinear effects in periodic structures and band-gap localized states beyond the approximation provided by the coupled-mode theory. Such an analysis is crucially important for determining stability of nonlinear waves in periodic structures because the wave instabilities can appear due to the mode coupling to other bands.

In this Letter, we present the first analysis of the existence and stability of nonlinear localized waves in a periodic medium with multiple gaps in the transmission spectrum that is valid beyond the coupled-mode theory. We consider a simple model where waves are localized in a layered medium by an intensity-dependent defect (or, in other words, they are guided by a thin-layer nonlinear waveguide). Assuming the applicability of our results to a variety of different physical systems (see below), we describe four qualitatively different situations and characterize, in the framework of a unified and systematic approach, the properties of two types of nonlinear localized waves . For the first time to our knowledge, we analyze stability of nonlinear localized waves in a periodic medium and reveal an important physical mechanism of wave instability associated with the band-gap resonances. We demonstrate also that several types of stable band-gap localized states can exist in the presence of nonlinearity.

Model. We describe localized waves in superlattices using the nonlinear Schrödinger equation for the wave envelope $`\psi (x,t)`$,

$$i\frac{\psi }{t}+\frac{^2\psi }{x^2}+(I;x)\psi =0,$$
(1)
where $`I|\psi |^2`$ is the wave intensity, $`t`$ is time (or propagation variable), $`x`$ is the spatial coordinate, and the real function $`(I;x)`$ describes both nonlinear and periodic properties of the medium. We note that the system (1) is Hamiltonian, and for localized solutions the power $`P=_{\mathrm{}}^+\mathrm{}I(x)𝑑x`$ is conserved.

We seek stationary localized solutions of Eq. (1) in the form $`\psi (x,t)=u(x)e^{i\omega t}`$, where $`\omega `$ is the normalized frequency (or the propagation constant, in optics), and the real function $`u(x)`$ satisfies the equation:

$$\omega u+\frac{d^2u}{dx^2}+(I;x)u=0.$$
(2)
We assume that the superlattice is linear, and nonlinearity appears only through the properties of an embedded localized defect. Then, if the corresponding width of the wave envelope is much larger than that of the defect, the inhomogeneity can be modeled by a delta-function and, in the simplest case, we can write $`(I;x)=\nu (x)+\delta (x)G(I)`$, where the function $`G(I)`$ characterises the properties of the defect, and $`\nu (x)\nu (x+h)`$ describes an effective potential of the superlattice with the spatial period $`h`$.

For such a local nonlinearity, the localized waves can be constructed with the help of a matching condition, by using the solutions of Eq. (2) with $`(I;x)=\nu (x)`$, presented in the form of the Bloch-type wave functions .

Band-gap structure and localized waves. If the effective periodic potential $`\nu (x)`$ is approximated by a piecewise-constant function (the so-called Kronig-Penney model), the solution on both sides of the defect can be decomposed into a pair of counter-propagating waves with the amplitudes $`a(x)`$ and $`b(x)`$,

$`u_b(x)=a(x)e^{\mu (x)x}+b(x)e^{\mu (x)x},`$

where $`\mu (x)=\sqrt{\omega \nu (x)}`$ is the local wavenumber. As follows from the Floquet-Bloch theory, for a Bloch-wave solution the reflection coefficient $`r(x)=b(x)/a(x)`$ is a periodic function, i.e. $`r(x)=r(x+h)`$, and it is found from the following eigenvalue problem:

$$T(x)\left(\begin{array}{c}1\hfill \\ r(x)\hfill \end{array}\right)=\tau (x)\left(\begin{array}{c}1\hfill \\ r(x)\hfill \end{array}\right),$$
(3)
where $`T(x)`$ is a transfer matrix that describes a change of the wave amplitudes $`\{a,b\}`$ after one period $`(x,x+h)`$. An explicit expression for the transfer matrix can be found in Ref. , where it was also proven that $`\mathrm{det}T1`$, and two linearly independent solutions of Eq. (3) correspond to a pair of the eigenvalues $`\tau `$ and $`\tau ^1`$. Relation $`\tau (\omega )`$ determines a band-gap structure of the superlattice spectrum: the waves are propagating, if $`|\tau |=1`$, and they are localized, if $`|\tau |1`$. In the latter case, a nonlinear defect can support nonlinear localized waves, and the wave amplitude at the defect is determined from the continuity condition  at $`x=0`$ \[$`I_0I(0)`$\]:

$$G_0G(I_0)=\zeta (\omega ).$$
(4)
Here $`\zeta =(\zeta ^++\zeta ^{})_{x=0}`$, $`\zeta ^\pm =\mu ^\pm \left(1r^\pm \right)\left(1+r^\pm \right)^1`$, and we denote with “$`+`$” and “$``$” the lattice characteristics on the right and the left side of the nonlinear defect, respectively, i.e. $`\nu (x)=\nu ^+(|x|)`$, for $`x>0`$, and $`\nu (x)=\nu ^{}(|x|)`$, for $`x<0`$. Note that the reflection coefficients should correspond to the Bloch-wave solutions with the asymptotics $`|u(x\pm \mathrm{})|0`$.

Relation (4) allows us to identify different types of nonlinear localized states. We notice that such localized states, supported by an attractive nonlinear defect ($`G_0>0`$), exist in the so-called waveguiding regime, when $`\zeta (\omega )>0`$. Additionally, localization can occur at a repulsive nonlinear defect ($`G_0<0`$) in the anti-waveguiding regime, provided $`\zeta (\omega )<0`$. A characteristic band-gap structure for a two-component superlattice is presented in Fig. 1 (top). Using the guided-wave terminology, we notice that the first (semi-infinite) band-gap corresponds to the conditions of the total internal reflection (IR) and, therefore, its dispersion properties should be similar to those of the conventional waveguiding regime. Localized waves in this band resemble conventional solitary waves modulated by a periodic structure \[see Fig. 1(b)\]. In contrast, both waveguiding and anti-waveguiding regimes can occur at smaller $`\omega `$, where band-gaps appear due to the resonant Bragg-type reflection (BR) from the periodic structure, so that the localized waves are somewhat similar to the gap solitons composed of mutually coupled backward and forward propagating waves \[see Fig. 1(a)\].

To study the linear stability, we consider the evolution of small-amplitude perturbations of the localized state presenting the solution in the form

$`\psi (x,t)=\left\{u(x)+v(x)e^{i\gamma t}+w^{}(x)e^{i\gamma ^{}t}\right\}e^{i\omega t},`$

and obtain the linear eigenvalue problem for $`v(x)`$ and $`w(x)`$. Then, from the condition of mode localization we define the so-called Evans function, $`Y(\gamma )=\left[G_1\zeta (\omega +\gamma )\right]\left[G_1\zeta (\omega \gamma )\right]\left(G_1G_0\right)^2`$, where $`G_1G_0+I_0G^{}(I_0)`$, the intensity $`I_0`$ is calculated for an unperturbed solution, and the prime stands for the derivative. Eigenvalues are found as zeros of the Evans function, $`Y(\gamma )=0`$, and the corresponding eigenmode solutions fall into one of the following categories: (i) internal modes with real eigenvalues describe periodic oscillations (“breathing”) of the localized state, (ii) instability modes correspond to purely imaginary eigenvalues, and (iii) oscillatory instabilities can occur when the eigenvalues are complex.

To demonstrate the basic stability results, we consider a localized defect possessing a cubic nonlinear response, $`G(I)=\alpha +\beta I`$. Under proper scaling the absolute value of the nonlinear coefficient $`\beta `$ can be normalized to unity, so that $`\beta =+1`$ corresponds to self-focusing and $`\beta =1`$ to self-defocusing nonlinearity. Localization depends also on the sign of the linear coefficient $`\alpha `$, which defines the defect response at small intensities: attractive if $`\alpha >0`$, and repulsive, otherwise. Therefore, below we consider four qualitatively different examples that correspond to the different signs of $`\alpha `$ and $`\beta `$.

Self-focusing nonlinearity. First we consider the properties of the localized waves supported by a defect with a self-focusing nonlinearity ($`\beta =+1`$), and attractive linear response ($`\alpha >0`$). Such waves already exist in the linear limit at the frequencies $`\omega _b`$ defined by the equality $`\zeta (\omega _b)=\alpha `$, and they correspond to the waveguiding regime only (white regions in Fig. 1, top). Since an IR wave is a fundamental eigenstate of the self-induced waveguide, it can be demonstrated that the conditions of the Vakhitov-Kolokolov (VK) stability theorem  are satisfied, and the IR states are unstable if and only if $`dP/d\omega <0`$. The critical point $`dP/d\omega =0`$ corresponds to a “collision” between two internal modes of the localized wave at the origin, as illustrated in Fig. 1 (middle). For the BR states, the VK criterion provides only a necessary condition for stability, since the higher-order localized states can also exhibit oscillatory instabilities. Indeed, we notice that in the linear limit there always exists an internal mode corresponding to a resonant coupling between the BR and IR band-gaps, since $`Y\left(\omega _b^{(\mathrm{IR})}\omega _b^{(\mathrm{BR})}\right)=0`$. We perform extensive numerical calculations and find that this mode leads to an oscillatory instability of BR waves when the value $`(\omega \mathrm{Re}\gamma )`$ moves outside the band gap; it happens when the intensity exceeds a threshold value (see Figs. 1 and 2).

Next, we study the linearly repulsive defect ($`\alpha <0`$) with self-focusing nonlinearity. In this case, the waveguiding-type states bifurcate from the band-gap edges \[with $`G_0=\zeta (\omega _b)=0`$\] and the stability criteria are the same as in the case $`\alpha >0`$ discussed above, see Fig. 3. The major difference from the previous case is the existence of the anti-waveguiding states at small intensities, $`I_0<|\alpha /\beta |`$. We find that such waves are stable in the lowest photonic band gap \[such as one shown in Fig. 3(b)\], while in higher band gaps they exhibit oscillatory instabilities due to a resonance with the higher-frequency states. Different types of the localized waves and their stability are summarized in Fig. 3.

Self-defocusing nonlinearity. Let us now consider the case $`\alpha <0`$ and $`\beta =1`$, when the total response of the defect is negative for any intensity, i.e. $`G(I)<0`$. In this case, the localized waves can only exist in the anti-waveguiding regime (see Fig. 4, top), the nonlinear waves continue the linear impurity states, and their frequency decreases at higher powers. We find that the localized states corresponding to the lowest BR gap, such as that shown in Fig. 4(b), are always stable. However, in other band gaps the localized states \[e.g., the one shown in Fig. 4(a)\] can exhibit oscillatory instability due to the physical mechanism identified in other cases — a resonance with the lower band gaps.

Finally, we study the case of a linearly attractive defect ($`\alpha >0`$) possessing a self-defocusing nonlinearity ($`\beta =1`$). At small intensities, the defect is attractive, and it can support localized waves in a waveguiding regime, in both IR and BR gaps, as shown in Fig. 5. Performing the stability analysis, we find that IR waves are always stable, while oscillatory instabilities appear for higher band-gap states. Since at higher intensities, i.e. for $`I_0>|\alpha /\beta |`$, the defect response changes its sign, a new type of localized waves can exist in the anti-waveguiding dispersion region, bifurcating from the band-gap edge. It is possible to demonstrate that for such waves the VK stability criterion becomes inverted, i.e. the waves are unstable if $`dP/d\omega >0`$, as shown in Fig. 5 (top). This happens because the signs of both the nonlinear response and effective dispersion are altered compared to the IR waves supported by a self-focusing defect. All higher band-gap states can exhibit oscillatory instabilities.

Concluding remarks. We have analyzed the existence and stability of nonlinear localized waves in one-dimensional periodic structures considering the simplest nonlinear generalization of the Kronig-Penney model with a nonlinear impurity for which both the nonlinear wave classification and stability analysis can be carried out in a complete and systematic way. Taking into account many common features of nonlinear guided waves and impurity modes in stratified and disordered media, on one side, and the self-trapped states and solitary waves in homogeneous nonlinear media, on the other side, we expect that many of the results described above will be found in more realistic physical models of nonlinear periodic media. Such cases include the electron self-trapping and locking states in cuprates and semiconductor superlattices ; nonlinear guided waves in optical superlattices ; impurity modes in photonic band-gap materials , magneto-optical periodic structures , and photonic crystal fibers ; coherent matter waves in optical lattices , etc. In particular, the stability of nonlinear localized impurity modes is a crucial issue for creating tunable band-gap materials where gaps could be controlled by changing the input light intensity.

We are indebted to O. Bang and C. M. Soukoulis for useful collaboration. The work has been partially supported by the Planning and Performance Fund of the Institute of Advanced Studies.


text:
# 1 Introduction

## 1 Introduction

The idea of parton saturation in large nuclei has recently obtained much attention \[1-5\]. In particular A.Mueller has discussed in  the behaviour of the quark and gluon densities for different values of momenta in a large nucleus. In his treatment the interaction of the probing current and the nucleons of the target was taken in the Glauber form, with the parton-nucleon amplitude taken in the two-gluon exchange approximation. The evolution of the gluonic distribution in $`x`$ has not been studied. The result of these study is that the behaviour of both the quark and gluonic distribution in the momentum $`l`$ is determined by a cerain critical ”saturation momentum” defined as

$$Q_s^2=c\frac{8\pi ^2\alpha _sN_c}{N_c^21}AT(b),$$
(1)
with $`T(b)`$ the standard nuclear profile function and $`c`$ essentially a constant. At high momenta $`l^2>>Q_s^2`$ the distributions retain their perturbative form

$$\frac{d\left(xq(x,Q^2)\right)}{d^2bd^2\mathrm{}}=\frac{N_c}{6\pi ^4}\frac{Q_s^2}{\mathrm{}^2}$$
(2)
and

$$\frac{d\left(xG(x,Q^2)\right)}{d^2bd^2\mathrm{}}=\frac{N_c^21}{4\pi ^4}\frac{Q_s^2}{\mathrm{}^2}\mathrm{ln}\frac{1}{x}\mathrm{for}\mathrm{}^2>>Q_s^2.$$
(3)
for the quark and gluon distributions respectively. However as the momentum $`l`$ diminishes both distributions do not infinitely grow (as would follow from (2) and (3)) but saturate at a certain finite value. In fact at $`l^2<<Q_s^2`$ A.Mueller found

$$\frac{d\left(xq(x,Q^2)\right)}{d^2bd^2\mathrm{}}=\frac{N_c}{2\pi ^4}$$
(4)
and

$$\frac{d\left(xG(x,Q^2)\right)}{d^2bd^2\mathrm{}}\frac{N_c^21}{4\pi ^4}\mathrm{ln}\frac{1}{x}\mathrm{for}\mathrm{}^2<<Q_s^2.$$
(5)
Independently of the saturation problem, in the last year the interaction of a probe with a large nucleus at small $`x`$ was studied in the perturbative QCD with reggeized gluons. It reduces to summing fan diagrams constructed of BFKL pomerons with a splitting triple Pomeron vertex \[6-9\]. In  the resulting equation was solved numerically, which enabled us to calculate the corresponding gluon distribution. Its features turned out to be somewhat different from what has been found by A.Mueller. Although at large $`l`$ the gluonic distribution indeed falls down similarly to the perturbative result and to what has been found in , it also falls down to zero at small $`l^2`$. The form of the found gluonic distribution as a function of $`\mathrm{log}l`$ is close to a Gaussian with a constant height and a center moving with a constant velocity to higher $`l`$ as the rapidity $`y`$ increases . Thus in the $`(\mathrm{log}l,y)`$ space the gluonic density turns out to be a soliton wave.

This unexpected behaviour of the gluonic density caused some mistrust due to the seemingly sound base of the analytic derivation of A.Mueller and a purely numerical character of the calculations in . Here we intend to study the reasons for the difference of these results.

One might think that this difference is due to a complicated gluonic interaction introduced in the approach pursued in our numerical calculation and absent in the picture of A.Mueller. However the answer seems to be much simpler: the gluonic densities studied in both papers are in fact different quantities. While the gluonic density numerically calculated in  is more or less standard, the one introduced by A.Mueller is a new quantity. We show that for this latter density our numerical results indeed lead to the same saturation as in . Moreover, A.Mueller’s quark density found from our solitonic gluon density also exhibits the same saturation properties as found in .

Besides showing that there is no contradiction between  and  at all, this confirms that the gluonic interaction at small $`x`$ and the resulting very sophisticated dynamical picture do not actually influence the saturation properties of A.Mueller’s densities, as conjectured in .

Of course there remains a question, which is the “correct” parton density. As we argue in the last section, we beleive that the densities introduced by A.Mueller do not in general describe partonic distributions in the target.

## 2 Gluonic distributions: standard and of A.Mueller

We begin with the definition of the gluonic density adopted in our numerical calculations, which we call ”standard”. It can be taken from the standard expression for the structure function of the nucleus in its terms (see e.g \[10 \])

$$F_2(x,Q^2)=\frac{4\alpha _sQ^2}{\pi ^2N_ce^2}d^2bd^2r[\rho _T(r)+\rho _L(r)]f(x,r,b),$$
(6)
where $`\rho _{T,L}`$ are the transverse and longitudinal colour densities created by the incoming photon and $`f`$ is expressed via the double density of gluons in momentum and impact parameter:

$$f(x,r,b)=\frac{d^2k}{(2\pi )^2k^2}\frac{d\left(xG(x,k^2,b)\right)}{d^2bdk^2}\left(1e^{ikr}\right)\left(1e^{ikr}\right).$$
(7)
In  it was found that, up to a factor, fuction $`f`$ coincides with a sum of BFKL fan diagrams $`\mathrm{\Phi }`$

$$f(x,r,b)=\frac{N_c}{2\pi ^3\alpha _s}\mathrm{\Phi }(x,r,b).$$
(8)
Taking a Fourier transform of (7) and neglecting the term proportional to $`\delta ^2(k)`$ we obtained

$$\frac{d\left(xG(x,k^2,b)\right)}{d^2bdk^2}=\frac{N_c}{2\pi ^2\alpha _s}k^2_k^2\varphi (x,k,b)\frac{N_c}{2\pi ^2\alpha _s}h(x,k,b),$$
(9)
where in the coordinate space

$$\varphi (x,r,b)=\frac{\mathrm{\Phi }(x,r,b)}{2\pi r^2}.$$
(10)
Function $`\varphi `$ was found in  by numerical solution of the BFKL fan diagram equation. Applying $`k^2_k^2`$ to it we found the gluon density. As mentioned in the Introduction, at any given rapidity the density turned out to have the same, roughly Gaussian, shape in variable $`\xi =\mathrm{ln}k`$, centered at the point $`\xi =\xi _0(y)`$, which moves to the right with a nearly constant velocity. Approximately the distribution can be described by

$$h(k)=h_0e^{a(\xi \xi _0(y))^2},$$
(11)
where $`h_0`$ and $`a`$ are practically independent of $`y`$ and $`\xi _0(y)`$ linearly grows with it:

$$h_00.3,\xi _0(y)=c+0.222\mathrm{ln}A+2.23\frac{\alpha _sN_c}{\pi }y.$$
(12)
It can be shown that it is normalized according to

$$_0^{\mathrm{}}\frac{dk}{k}h(y,k,b)=1,$$
(13)
which condition relates $`a`$ and $`h_0`$. The starting position of the center $`c`$ depends on the initial distribution at small rapidities.

Evidently at a given value of $`k`$ the density always remains limited, irrespective of the form of the initial distribution (and on the atomic number $`A`$, in particular). In this sense we have saturation as discussed in \[1-5\]. However with the growth of $`y`$ the strongly peaked density moves away toward higher values of $`k`$ so that the density at a fixed point tends to zero at high values of rapidity. We thus find “supersaturation”: with $`y\mathrm{}`$ the gluon density at an arbitrary finite momentum tends to zero.

Now we pass to the gluon density introduced by A.Mueller. He related it to a particular process, analogous to the $`\gamma ^{}`$-target scattering, in which the electromagnetic current is changed to a “gluonic” current

$$j(x)=\frac{1}{4}F_{\mu \nu }^iF_{\mu \nu }^i.$$
(14)
This current generates a colourless pair of gluons, which according to  interact with the target very much in the same way as a $`q\overline{q}`$ pair generated by the electromagnetic current (see Fig. 1, diagrams $`b`$ and $`c`$). The gluon density is then defined as a contribution to the cross-section from a given value of the transverse momentum $`l`$ of the struck gluon on the left part of the diagrams.

For the case when the target is just a nucleon, in the limit $`Q^2\mathrm{}`$ and in the single gluon exchange exchange approximation, A.Mueller found for his density, which we denote $`xG_M(x,Q^2,l^2)`$ (it generally depends on two different momenta):

$$\frac{d\left(xG_M(x,Q^2,\mathrm{}^2)\right)}{d^2\mathrm{}}=\frac{\alpha N_c}{2\pi ^5}_{x\mathrm{}^2/Q^2}^{\mathrm{}^2/Q^2}\frac{dz}{z}\frac{d^2r_1d^2r_2}{r_1^2r_2^2}e^{i\mathrm{}(r_1r_2)}(\frac{2(r_1r_2)^2}{r_1^2r_2^2}1)$$

$$(1+e^{ik(r_1r_2)}e^{ikr_1}e^{ikr_2})\frac{d^2k}{k^2}\frac{d\left(xG(x,k^2)\right)}{dk^2}.$$
(15)
On the right-hand side there appeared $`xG(x,k^2)`$, which is just the standard gluonic density. In the approximations made in deriving (15) it does not depend on $`x`$, nor the whole expression depends on $`Q^2`$. As a result, the gluon density defined by (15) also depends on only one momentum squared $`\mathrm{}^2`$ and superficially looks similar to the standard gluon density. But in fact it is a different quantity (factor $`\alpha _s`$ in front of (15) is already sufficient to see it). Its physical meaning has been explained before: it is a probability to find one of the struck gluons (the left one in the Fig. 1) with momentum $`\mathrm{}`$. Later we shall discuss to what extent this quantity probes the gluon distribution in the nucleon.

We see that once the standard $`xG`$ is known, the density $`xG_M`$ can be found directly from Eq. (15). This exersize is done in the next section.

## 3 The soliton density and saturation

The first thing to do is to generalize $`xG_M`$ to the nucleus target. Having in mind the fan diagram picture for the interaction of a probe with the nucleus, on the one hand, and the scattering process considered by A.Mueller to define his density, on the other, all we have to do is to change the standard density $`xG`$ on the right- hand side of Eq. (15) to the same quantity for the nucleus:

$$\frac{d\left(xG_M(x,Q^2,\mathrm{}^2,b)\right)}{d^2bd^2\mathrm{}}=\frac{\alpha N_c}{2\pi ^5}_{x\mathrm{}^2/Q^2}^{\mathrm{}^2/Q^2}\frac{dz}{z}\frac{d^2r_1d^2r_2}{r_1^2r_2^2}e^{i\mathrm{}(r_1r_2)}(\frac{2(r_1r_2)^2}{r_1^2r_2^2}1)$$

$$(1+e^{ik(r_1r_2)}e^{ikr_1}e^{ikr_2})\frac{d^2k}{k^2}\frac{d\left(xG\right)}{d^2bdk^2}.$$
(16)
Here we suppresed all evident arguments in $`G`$ on the right-hand side. This will be also often done in the following when the arguments are obvious.

We want to calculate this expression taking for $`d(xG)/d^2bdk^2`$ the solitonic gluon distribution in a nucleus obtained by our numerical studies. In principle this can be done rigorously, by numerical integration. However, since we are interested in the qualitative features of the saturation behaviour, we shall adopt here a simpler approach. We shall use the fact that at large rapidities the solitonic distribution is strongly peaked at the momena determined by

$$k_0=e^{\xi _0}A^{2/9}\mathrm{exp}\left(2.23\frac{N_c\alpha _s}{\pi }y\right)$$
(17)
(which are very high already at reasonable $`y`$ for not very small $`\alpha _s`$). This momentum will play the role of Mueller’s saturation momentum $`Q_s`$. Assuming that $`k_0^2<<Q^2`$ one can take all the integrand out of the integral over $`k`$ at point $`k=k_0`$ except for the density factor. The latter, integrated over $`k^2`$ gives $`N_c/(2\alpha _s\pi ^2)`$ according to (9) and the normalization condition (13). So we are left with the expression

$$\frac{d\left(xG_M\right)}{d^2bd^2\mathrm{}}=\frac{N_c^2}{4\pi ^7}\mathrm{ln}\frac{1}{x}𝑑\varphi \frac{d^2r_1d^2r_2}{r_1^2r_2^2}e^{i\mathrm{}(r_1r_2)}(\frac{2(r_1r_2)^2}{r_1^2r_2^2}1)$$

$$(1+e^{ik(r_1r_2)}e^{ikr_1}e^{ik,r_2})$$
(18)
where one should take $`k=k_0n`$ with a two dimensional unit vector $`n`$.

Calculation of the integrals over $`r_1`$ and $`r_2`$ is simplified by using the identity

$$\frac{1}{r_1^2r_2^2}\left(\frac{2(r_1r_2)^2}{r_1^2r_2^2}1\right)=\frac{1}{2}(_1_2)^2\mathrm{ln}r_1\mathrm{ln}r_2.$$
(19)
One obtains

$$\frac{d\left(xG_M\right)}{d^2bd^2\mathrm{}}=\frac{N_c^2}{8\pi ^7}\mathrm{ln}\frac{1}{x}𝑑\varphi \left[X(l,l)+X(l+k,l+k)2X(l+k,l)\right],$$
(20)
where

$$X(l_1,l_2)=d^2r_1d^2r_2e^{il_1r_1+il_2r_2}(_1_2)^2\mathrm{ln}r_1\mathrm{ln}r_2=4\pi ^2\frac{(l_1l_2)^2}{l_1^2l_2^2}.$$
(21)
This leads to

$$\frac{d\left(xG_M\right)}{d^2bd^2\mathrm{}}=\frac{N_c^2}{\pi ^5}\mathrm{ln}\frac{1}{x}𝑑\varphi \frac{l^2k_0^2(lk)^2}{l^2(l+k)^2}.$$
(22)
Let us consider the two limiting cases studied in . If $`l>>k_0`$ then one finds from (22)

$$\frac{d(xG_M)}{d^2bd^2\mathrm{}}=\frac{N_c^2}{\pi ^4}\frac{k_0^2}{l^2}\mathrm{ln}\frac{1}{x}.$$
(23)
In the opposite case $`l^2<<k_0^2`$ we find

$$\frac{d(xG_M)}{d^2bd^2\mathrm{}}=\frac{N_c^2}{\pi ^4}\mathrm{ln}\frac{1}{x}$$
(24)
This behaviour is identical with (3) and (5) (for $`N_c>>1`$ and except for factor 4 in (24)).

So Mueller’s density indeed saturates to a finite value at small $`l^2`$ irrespective of the dynamical picture. The standard density however behaves very differently and for a heavy nucleus goes to zero at any fixed value of $`l`$ as the rapidity gets large enough.

## 4 The quark density from the soliton gluonic density

Following  and similarly to the gluon case, one may introduce a quark distribution depending on two momenta squared $`xq_M(x,Q^2,\mathrm{}^2)`$, as a contribution from the left struck quark to the cross-section corresponding to diagrams $`b`$ and $`c`$ in Fig. 1. In the two-gluon exchange approximation for the interaction, this quantity again results independent of $`Q^2`$ at high virtualities. We postpone discussion of the relation of this limiting distribution in $`\mathrm{}`$, $`xq_M(x,\mathrm{}^2)`$ to the standard quark density as a function of $`Q^2`$, $`xq(x,Q^2)`$, until the last section. Here we shall consider the properties of $`xq_M(x,\mathrm{}^2)`$ which follow from the found standard gluon distribution in a large nucleus.

The definition of $`xq_M(x,Q^2,\mathrm{}^2)`$ is very similar to that of (15). For a single nucleon as a target A.Mueller found:

$$\frac{d\left(xq(x,Q^2,\mathrm{}^2)\right)}{d^2\mathrm{}}=\frac{\alpha Q^2}{4\pi ^3}_0^1𝑑z[z^2+(1z)^2]\left(\frac{\mathrm{}^2}{\mathrm{}^2+ϵ^2}\frac{(\mathrm{},\mathrm{}+k)}{(\mathrm{}^2+ϵ^2)[(\mathrm{}+k)^2+ϵ^2]}\right)\frac{d^2k}{k^2}\frac{dxG}{dk^2},$$
(25)
where

$$ϵ^2=Q^2z(1z).$$
(26)
As one observes, the quark density $`xq_M`$ is again defined via the standard gluon density.

For the nucleus target, as explained before, all we have to do is to take the corresponing nuclear density on the right-hand side of (25):

$$\frac{d\left(xq_M\right)}{d^2bd^2\mathrm{}}=\frac{\alpha Q^2}{4\pi ^3}_0^1𝑑z[z^2+(1z)^2]\left(\frac{\mathrm{}^2}{\mathrm{}^2+ϵ^2}\frac{(\mathrm{},\mathrm{}+k)}{(\mathrm{}^2+ϵ^2)[(\mathrm{}+k)^2+ϵ^2]}\right)\frac{d^2k}{k^2}\frac{dxG}{d^2bdk^2}.$$
(27)
At high $`Q^2`$ evidently only small $`z`$ or $`1z`$ contribute. Then integrating over $`ϵ^2Q^2z`$ we find

$$\frac{d\left(xq_M\right)}{d^2bd^2\mathrm{}}=\frac{\alpha }{2\pi ^3}\frac{d^2k}{k^2}\frac{dxG}{d^2bdk^2}\left(1\frac{(l,l+k)}{(l+k)^2l^2}\mathrm{ln}\frac{(l+k)^2}{l^2}\right).$$
(28)
Let us take for $`xG`$ the soliton-like gluon density found in . As before we take all smooth functions out of the integral over $`k`$ and afterwards integrate over $`k`$ using (13). We obtain

$$\frac{d\left(xq_M\right)}{d^2bd^2\mathrm{}}=\frac{N_c}{4\pi ^5}𝑑\varphi \left[1\frac{(l,l+k)}{(l+k)^2l^2}\mathrm{ln}\frac{(l+k)^2}{l^2}\right]$$
(29)
where $`k=k_0n`$ as in the preceding section.

Again we shall concentrate on the two limiting cases. If $`l>>k_0`$ using

$$\frac{1}{(l+k)^2l^2}\mathrm{ln}\frac{(l+k)^2}{l^2}\frac{1}{l^2}\left(1\frac{1}{2}\frac{2(lk)+k^2}{l^2}+\frac{4}{3}\frac{(lk)^2}{l^4}\right)$$

and integrating over the angle $`\varphi `$ we get

$$\frac{d\left(xq_M\right)}{d^2bd^2\mathrm{}}=\frac{N_c}{6\pi ^4}\frac{k_0^2}{l^2}.$$
(30)
In the opposite case $`l^2<<k_0`$ the second term in the bracket of (29) does not contribute and we get

$$\frac{d\left(xq_M\right)}{d^2bd^2\mathrm{}}=\frac{N_c}{2\pi ^4}.$$
(31)
This behaviour exactly coincides with A.Mueller results (2) and (4).

## 5 Discussion

As we have seen, there are in fact two different gluon densities, which exhibit different behaviour as functions of momenta. The standard one defined via Eq (7) has a soliton-like behaviour and vanishes at any given $`\mathrm{}`$ as $`y\mathrm{}`$. The other, defined in , saturates at small momenta to a finite value, which turns out to be universal in the sense that it does not depend on the complexity of the interaction mechanism. Which of the two densities is the “correct” one?

The definition of the standard gluon density via Eq. (7) seems to us perfectly adequate. Its validity is clear from the inspection of the diagrams $`b`$ and $`c`$ of the Fig. 1 corresponding to (7). The lower gluonic blob, to which the quark loop is coupled, has to be integrated over the “-” component of the 4-momentum k (in a system, in which the target has momentum with a large “+” component). Evidently ihe blob then represents just the gluonic density at a given “time” $`x_+`$.

Interpretation of Mueller’s density, in contrast, is not so obvious. First, in the presence of the interference diagram $`c`$ of Fig. 1, one may hardly say that the contribution from the left struck gluon to the cross-section probes the gluon distribution in the target. In the interference diagram the right struck gluon has a different transverse momentum $`\mathrm{}+k`$, so that the interpretation of the whole contribution as a gluon density is doubtful. Then it is not clear why the lowest order diagram Fig. 1$`a`$ has not been included. Finally, in the BFKL kinematics the virtual gluons tend to reggeize. One then expects that the interference digram can be neglected and diagram Fig. 1$`b`$ will be absorbed into $`a`$ to form a BFKL pomeron coupled to the target. So in the end one will find only the diagram Fig. 1$`a`$, which leads to the standard gluon density as defined through Eq. (7).

The same doubts refer to the quark density in the definition of A.Mueller. Again, so long as one cannot neglect the interference diagram, the interpretation of the contribution to the cross-section from the left struck quark with a given momentum as a quark density of the target does not seem possible. True, in the DGLAP kinematics, with finite $`x`$ , the interference diagram can indeed be neglected, substituted by the requirement that $`\mathrm{}^2<<Q^2`$. However then one has the transverse momentum ordering and the density found from Eq. (25) would be trivial (see , Eq. (14)). In the BFKL regime, at small $`x`$, the virtual quarks also tend to reggeize and the contribution from digram Fig. 1$`c`$ hopefully becomes neglegible as compared to $`a`$ and $`b`$, which actually become identical. Then the interpretation of (25) as the quark density may aquire sense. However to calculate it one should couple the found gluonic distribution to a reggeized quark. We leave this problem for future studies.

## 6 References

1. L.V. Gribov, E.M. Levin and M.G. Ryskin, Phys.Rep.100 (1983)1.

2. A.H. Mueller, Nucl.Phys.B335 (1990) 115.

3. J.Jalilian-Marian, A. Kovner, L.McLerran and H. Weigert, Phys.Rev.D55 (1997) 5414.

4. Yu. V. Kovchegov and A.H. Mueller, Nucl.Phys.B529 (1998) 451.

5. A.H.Mueller, hep-ph/9904404

6. I.Balitsky, hep-ph/9706411; Nucl. Phys. B463 (1996) 99.

7. Yu. Kovchegov, Phys. Rev D60 (1999) 034008;
preprint CERN-TH/99-166 (hep-ph/9905214).

8.E.Levin and K.Tuchin, preprint DESY 99-108, TAUP 2592-99 (hep-ph/9908317).

9. M.A.Braun, Eur.Phys. J. C16 (2000) 337-348.

10. M.A.Braun, Eur. Phys. J C6 (1999) 343.