format: text

text:
# Combining Hebbian and reinforcement learning in a minibrain model

## 1 Introduction

The central question which we address in this article is in what way a biological neural network, i.e., the brain, or, more generally, a part of the nervous system of an animal or human being, may learn to realize input–output relations. By biological we here mean: realizable with the help of elements occurring in nature, e.g., neurons or chemical substances that may influence other neurons or the synaptic efficacy.

An example of an input–output relation is a motor task, like catching a prey, in reaction to visual, auditive, or other input. Many attempts to explain the way input–output relations of this kind might be realized by (artificial) neural nets are encountered in the literature, most of which are not satisfactory from a biological point of view as we will illustrate in subsection 1.1.

It is the purpose of this article to combine ideas which do satisfy certain biological constraints and study a toy model, in particular with respect to its ability to learn and realize input–output relations.

The widely accepted idea of Hebbian learning  at the one hand will be combined with some rule that implements a feedback signal at the other hand, in a way that, in principle, might be biologically realizable. Without the addition of any feedback-signal, learning of prescribed input–output relations —whether in reality or in a model— is, of course, impossible.

### 1.1 Artificial learning rules

If one wants a network to learn to realize input–output relations, there are various well-known prescriptions, associated with names like perceptron learning rule, back-propagation or Boltzmann machines . None of these, however, model the functioning of real brains, since the learning rules in question violate the existing biological limitations. In order to illustrate this statement, let us give an example.

Consider a single layered feed-forward network, i.e., a network consisting of an input and an output layer only, in which signals are sent by neurons of the input layer to neurons of the output layer, and not the other way around. Let $`w_{ij}`$ be the strengths or ‘weights’ of the connections in this simple net. In 1962, Rosenblatt  proved that such a network will realize desired input–output relations if, a finite number of times, the weights are adapted according to the rule

$$w_{ij}w_{ij}+\mathrm{\Delta }w_{ij}$$
(1)
with

$$\mathrm{\Delta }w_{ij}=\epsilon (x_{\mathrm{T}i}x_{\mathrm{O}i})x_j$$
(2)
where $`x_{\mathrm{T}i}`$ is the desired or target output of neuron $`i`$, and $`x_{\mathrm{O}i}`$ is its actual output. Furthermore, $`x_j`$ is the state of the pre-synaptic input neuron $`j`$ and $`\epsilon `$ is some function of the neuron states and properties of neurons $`i`$ and $`j`$. This learning rule can not be realized by a biological neural net since neuron $`i`$, producing $`x_{\mathrm{O}i}`$, cannot know that it should produce $`x_{\mathrm{T}i}`$. If, e.g., an animal does not succeed in catching a prey, its neurons get no specific feedback individually, on what the right output $`x_{\mathrm{T}i}`$ should have been. Hence, $`x_{\mathrm{T}i}x_{\mathrm{O}i}`$ cannot be determined by the biological system, and, therefore, neither can it adapt the weights according to (2). Consequently, the perceptron learning rule (2) is unsuitable for a realistic modeling of the way in which a biological neural net can learn and realize input–output relations. Similar observations can be made for back-propagation or Boltzmann machines.

### 1.2 Biological learning rules; Hebbian learning and reinforcement learning

Already in 1949, Hebb suggested  that, in biological systems, learning takes place through the adaptation of the strengths of the synaptic interactions between neurons, depending on the activities of the neurons involved. In a model using binary neurons, i.e., $`x_i=0`$ or $`x_i=1`$, the most general form of a learning rule based on this principle is a linear function in $`x_i`$ and $`x_j`$ since $`x_i^2=x_i`$ and $`x_j^2=x_j`$. It therefore reads

$$\mathrm{\Delta }w_{ij}=a_{ij}+b_{ij}x_i+c_{ij}x_j+d_{ij}x_ix_j$$
(3)
In a biological setting, the coefficients $`a_{ij},b_{ij},c_{ij}`$ and $`d_{ij}`$ in this learning rule can only depend on locally available information, such as the values of the membrane potential

$$h_i=\underset{j}{}w_{ij}x_j$$
(4)
and the threshold potential $`\theta _i`$ of neuron $`i`$. In this way, the system adapts its weights without making use of neuron-specific information, like, e.g., $`x_{\mathrm{T}i}`$, of which there can be no knowledge, locally, at the position of the synapse.

In a recurrent neural net, a Hebbian learning rule suffices to store patterns  if all neurons are clamped to the patterns which are to be learned during a certain period, the ‘learning stage’. In feed-forward networks, however, only the neurons of the input layers are clamped, and some kind of feed-back signal, governing the direction of adaptation of the synapses during the learning procedure, is indispensable. Probably the simplest form of such a signal one can think of is a ‘success’ or ‘failure’ signal that is returned to the network after each attempt to associate the correct output to given input. On the basis of trial and error, a neural network can then indeed learn certain tasks, the principle on which ‘reinforcement learning’ is based . This principle of reinforcement learning has a rather natural interpretation: satisfactory behavior is rewarded, or *reinforced*, causing this behavior to occur more frequently. The reinforcement signal is supplied by the subject’s environment, or by its own judgment of the effect of its behavior. In a biological perspective, one could think of the synaptic change being influenced by some chemical substance, which is released depending on whether the evaluation by the subject of the effect of the output is positive or negative, i.e., whether it is happy or unhappy with the attempt it made.

Note that, in learning by reinforcement, the search for the correct output is more difficult, and, hence, slower, than for non-biologically realizable algorithms like the perceptron learning rule or back-propagation. This is not surprising, since the latter give the system *locally* specific information on how to adjust individual weights, while reinforcement rules only depend upon a *global* ‘measure of correctness’.

The most general form of a reinforcement learning algorithm is given by the prescription

$$\mathrm{\Delta }w_{ij}=a_{ij}(r)+b_{ij}(r)x_i+c_{ij}(r)x_j+d_{ij}(r)x_ix_j$$
(5)
Here the coefficients $`a_{ij},b_{ij},c_{ij},d_{ij}`$, besides their dependence on the local variables such as the membrane potential, will in principle also depend on a reinforcement signal, denoted by $`r`$. The value of $`r`$ is usually a real number between $`0`$ and $`1`$, denoting the degree of success ($`r=1`$ means success, $`r=0`$ means failure).

An important issue in the literature on reinforcement learning is the so called ‘credit assignment problem’ . It refers to the question how a neural net knows which connections $`w_{ij}`$ were responsible for a successful or unsuccessful trial, and, as a consequence, which connections should be ‘rewarded’, and which should be ‘punished’, respectively.

In their article ‘Learning from mistakes’ (1999), Chialvo and Bak , proposed a class of networks, in which learning takes place on the basis of a ‘deinforcement signal’ only, i.e., the weights of active synapses are decreased if the output is wrong, they are ‘punished’, so to say, in case of wrong performance of the network. If the output is right nothing happens. This procedure works as long as the average activity in the net is kept very low: when only a few neurons are active at an unsuccessful attempt, one can be sure that the connections between these active neurons were the ones which were responsible, and thus should be ‘punished’. In this way Chialvo and Bak obtained an elegant solution to the credit assignment problem.

The absence of a reinforcement signal (nothing happens if $`r=1`$) makes their learning rule relatively simple. It is a version of the general rule (5) with $`b_{ij}=0`$ and $`c_{ij}=0`$:

$$\mathrm{\Delta }w_{ij}=(1r)(\rho x_ix_j\phi )$$
(6)
where $`\rho `$ and $`\phi `$ are positive constants; in this article we will suppose $`\phi <<\rho `$. A biological mechanism that could implement the learning rule (6) is the following: if the output is correct, nothing happens, since the network obviously performs satisfactory. If not, a chemical substance is released, which has the effect that synapses between neurons that have just been active, and thereby are ‘tagged’ in some electro-chemical way, are depressed.

### 1.3 Purpose

The success of the ‘minibrain model’ of Chialvo and Bak  (as Wakeling and Bak referred to it in ), is limited to feed-forward neural nets in which the number of input and output neurons (or, equivalently in this model, the number of patterns) is small compared to the number of neurons in the hidden layer. As the number of neurons in the hidden layer decreases, learning, at a certain moment, becomes impossible: ‘path interference’ is the phenomenon which causes this effect . Essentially, it amounts to the following. If, in each layer of the feed-forward neural net, only one neuron is active at each time step, an input–output relation corresponds to a path of activity along the strongest connections between the neurons. Basically, path interference comes down to the erasure of an existing path of activity, which was correct at a previous learning step, by a change due to a punishment of a connection while trying to learn a different input–output relation. If the probability for this path interference to occur becomes too large, learning times tend to infinity.

In this article we attempt to improve the performance of the minibrain model of Chialvo and Bak —in the sense of decreasing the learning time— by making sure that, at the occurrence of path interference, the punishment of a correct activity path is no longer such that the memory is erased. We achieve this by adding to the deinforcement term in the learning rule (6), which is proportional to $`\rho `$, a Hebbian term proportional to $`\eta `$. The latter term has the effect that an active path is strengthened, mitigating in this way the punishment. By choosing the ratio between the coefficients $`\eta `$ and $`\rho `$ of both terms in the learning rule appropriately, we are able to reduce the number of learning steps significantly, without making the model less realistic from a biological point of view. In fact, in the class of models we study, Hebbian learning is a most appropriate way to account for biological observations like LTP and LTD . In section 4 we explain that if the quotient of the Hebbian learning rate and the coefficient of the deinforcement term is in the range

$$\frac{1}{4}<\frac{\eta }{\rho }<\frac{1}{2}$$
(7)
learning times are reduced considerably.

In their article , Chialvo and Bak proposed a different way to solve the problem of path interference. They reduced the amount of punishment of the connections that once had been active in a correct response. In this model a neuron needs to memorize whether it previously was active in a successful trial. In our model such a neuron memory is not needed.

Let us denote the deinforcement contribution to learning by $`\mathrm{\Delta }w_{ij}^{}`$ and denote the Hebbian part by $`\mathrm{\Delta }w_{ij}^{\prime \prime }`$. We will study a learning rule of the form

$$\mathrm{\Delta }w_{ij}=\mathrm{\Delta }w_{ij}^{}+\mathrm{\Delta }w_{ij}^{\prime \prime }$$
(8)
From all possibilities for Hebbian learning summarized by equation (3), we choose for $`\mathrm{\Delta }w_{ij}^{\prime \prime }`$ a rule in which the coefficients $`a_{ij}`$ and $`b_{ij}`$ both are zero:

$$\mathrm{\Delta }w_{ij}^{\prime \prime }=\epsilon (x_i,x_j)(2x_i1)x_j$$
(9)
We choose this particular form since it can be argued that this form is a most plausible candidate from a biological point of view .

Our paper has been set up as follows. In section 2, we describe a feed-forward network with one hidden layer, which we will use to study the learning rule (8), with (6) and (9). In section 3, numerical studies of various situations are given and explained. It turns out, in general, that taking into account Hebbian learning, and viewing it as a process which is permanently active, irrespective of the occurrence of reinforcement learning, has a positive influence on the learning time of the neural net. This is a new result, which, to the best of our knowledge, has not been noticed before.

## 2 Description of the model: updating rules for neuron activities and connection weights

In order to explore a simplified model of the brain we consider a fully connected, feed-forward neural network with an input layer of $`N_\mathrm{I}`$ neurons, one hidden layer of $`N_\mathrm{H}`$ neurons, and an output layer of $`N_\mathrm{O}`$ neurons, see figure 1.

The state $`x_i`$ of neuron $`i`$ is $`1`$ if neuron $`i`$ fires and $`0`$ if it does not. In general, a neuron fires if its potential $`h_i`$ is sufficiently high, where $`h_i`$ stands for the membrane potential $`V_{\mathrm{ex}}V_{\mathrm{in}}`$, the difference between the intra- and extra cellular potentials $`V_{\mathrm{in}}`$ and $`V_{\mathrm{ex}}`$. Following Chialvo and Bak , we model the dynamics of the neural net by simply stating that in the hidden and output layers a *given* number of neurons having the highest potentials $`h_i`$ —in their respective layers— are those that will fire. In their terminology: we use *extremal dynamics*.

For McCulloch and Pitts neurons a way to control the average activity might be realized via suitably chosen threshold potentials $`\vartheta _i`$ (see e.g. , ). In nature, the average activity will depend on the threshold potentials and may, moreover, be influenced by chemical substances or the network structure . In our model we put the threshold potentials $`\theta _i`$ equal to zero:

$$\theta _i=0$$
(10)
The number of neurons in the input, hidden and output layers that we choose to be *active*, will be denoted by $`N_\mathrm{I}^{(\mathrm{a})}`$, $`N_\mathrm{H}^{(\mathrm{a})}`$ and $`N_\mathrm{O}^{(\mathrm{a})}`$, respectively.

The input pattern, a specific set of states of neurons in the input layer, will be denoted by $`𝝃_\mathrm{I}=(\xi _{\mathrm{I1}},\mathrm{},\xi _{\mathrm{I}N_\mathrm{I}})`$. The network is to associate every input pattern with a desired, or target, output pattern, $`𝝃_\mathrm{T}=(\xi _{\mathrm{T1}},\mathrm{},\xi _{\mathrm{T}N_\mathrm{O}})`$. The $`𝝃_\mathrm{I}`$ and $`𝝃_\mathrm{T}`$ are vectors with components equal to $`0`$ or $`1`$. Consequently, the number of active neurons of the input and output layer are given by

$`N_\mathrm{I}^{(\mathrm{a})}`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{N_\mathrm{I}}{}}}\xi _{\mathrm{I}j}`$ (11)
$`N_\mathrm{O}^{(\mathrm{a})}`$ $`=`$ $`{\displaystyle \underset{j=1}{\overset{N_\mathrm{O}}{}}}\xi _{\mathrm{T}j}`$ (12)
In our numerical experiments, these numbers will be taken to be equal. Moreover, the number of active neurons in the hidden layer,

$$N_\mathrm{H}^{(\mathrm{a})}=\underset{j=1}{\overset{N_\mathrm{H}}{}}x_{\mathrm{H}j}$$
(13)
where $`x_{\mathrm{H}j}`$ is the neuron state of neuron $`j`$ in the hidden layer, will also be equal to the number of active neurons in the other layers. Hence, we *choose* $`N_\mathrm{I}^{(\mathrm{a})}=N_\mathrm{O}^{(\mathrm{a})}=N_\mathrm{H}^{(\mathrm{a})}`$.

We thus have explicitated the network dynamics. We now come to the update rules for the synaptic weights $`w_{ij}`$. Updating of the network state will happen at discrete time steps. At every time step $`t_n`$, all neuron states are updated in the order: input layer – hidden layer – output layer. This being done, the values of the weights are updated, according to

$$w_{ij}(t_{n+1})=w_{ij}(t_n)+\mathrm{\Delta }w_{ij}(t_n)$$
(14)
Substituting (9) and (6) into (8), we find

$$\mathrm{\Delta }w_{ij}=\epsilon (x_i,x_j)(2x_i1)x_j(1r)(\rho x_ix_j\phi )$$
(15)
For the pre-factor of the Hebbian term we take

$$\epsilon (x_i,x_j)=\eta \left(\kappa (2x_i1)(h_i\theta _i)\right)$$
(16)
The constants $`\eta `$ and $`\kappa `$ are positive numbers, the so-called learning rate and margin parameter. Finally, combining the above ingredients and noting that we chose $`\theta _i=0`$, the learning rule reads:

$`\mathrm{\Delta }w_{ij}(t_n)`$ $`=`$ $`\eta [\kappa h_i(t_n)(2x_i(t_n)1)][2x_i(t_n)1]x_j(t_n)`$ (17)
$`+(1r)[\rho x_i(t_n)x_j(t_n)+\phi ]`$
Note that $`x_i(t_n)`$ and $`x_j(t_n)`$ are the activities of neurons in *adjacent* layers, since in our model there are no lateral connections. The constant $`\phi `$ is chosen in such a way that the change in $`_{i,j}w_{ij}`$, where the sum is extended over $`i`$ and $`j`$ in adjacent layers, due to the $`\rho `$-term (not due to the Hebbian term), is independent of $`\rho `$. This can easily be achieved by choosing $`\phi =\rho /P`$, where $`P`$ is the product of the numbers of neurons in two adjacent layers, i.e., $`\phi `$ is equal to either $`\rho /N_\mathrm{I}N_\mathrm{H}`$ or $`\rho /N_\mathrm{H}N_\mathrm{O}`$.

## 3 Numerical Simulations

The network described in the previous section will now first be studied numerically. The numerical details are:

* The initial weights $`w_{ij}(t_0)`$ are randomly chosen with values between $`0.01`$ and $`+0.01`$.
* The punishment rate $`\rho `$ will be kept constant at $`0.02`$. Thus when we vary the $`\eta /\rho `$ ratio, we vary the learning rate $`\eta `$.
* The margin parameter $`\kappa `$, appearing in (16), will be kept at the value $`1`$.
* Whenever the number of neurons in the input, hidden or output layer is fixed, we choose $`N_\mathrm{I}=8`$, $`N_\mathrm{H}=512`$ and $`N_\mathrm{O}=8`$.
* All data are averaged over $`512`$ samples.

The network is confronted with $`p`$ different input patterns $`𝝃_\mathrm{I}^\mu `$, $`(\mu =1,\mathrm{},p)`$, to which correspond equally many target patterns $`𝝃_\mathrm{T}^\mu `$. Learning proceeds as follows. The input layer remains clamped to the first input pattern until the time step at which the target pattern has been found. As soon as this input-output relation $`\mu =1`$ has been learned, we switch to input pattern $`\mu =2`$. After the corresponding target pattern has been found we continue, up to the $`p`$-th input-output relation. Then, we have completed what we will refer to as one ‘learning cycle’.

After this cycle we start the process again, up to the point where the network can recall all $`p`$ input-target relations at once. When that is the case, learning stops. We count the number of learning steps needed to learn all input–output relations.

Chialvo and Bak consider the case of one active neuron in each layer. In section 3.1 we present a numerical experiment with a neural network for which the activities are larger than one, i.e., $`N_\mathrm{I}^{(\mathrm{a})}>1`$, $`N_\mathrm{H}^{(\mathrm{a})}>1`$ and $`N_\mathrm{O}^{(\mathrm{a})}>1`$. In particular we study the total number of learning steps as a function of the ratio $`\eta /\rho `$. In section 3.2 we vary the number of neurons in input and output layer and keep the hidden layer fixed, and vice versa. Finally, in section 4, we interpret our results.

### 3.1 Effect of the Hebbian term

In ‘Learning from mistakes’ Chialvo and Bak  studied the case of *one* active neuron in the input, the hidden and the output layers: $`N_\mathrm{I}^{(\mathrm{a})}=N_\mathrm{H}^{(\mathrm{a})}=N_\mathrm{O}^{(\mathrm{a})}=1`$. We here will study what happens when $`N_\mathrm{I}^{(\mathrm{a})}`$, $`N_\mathrm{H}^{(\mathrm{a})}`$ and $`N_\mathrm{O}^{(\mathrm{a})}`$ all are larger than one.

In our first numerical experiment we take a network with $`p=8`$ input–target relations for which, in each input or target pattern $`\mu `$, the number of active neurons is $`2`$, i.e., $`N_\mathrm{I}^{(\mathrm{a})}=N_\mathrm{H}^{(\mathrm{a})}=N_\mathrm{O}^{(\mathrm{a})}=2`$. In figure 2 the number of learning steps is plotted against the ratio $`\eta /\rho `$ of the two proportionality coefficients related to the Hebbian and the deinforcement term respectively.

From figure 2 we see that, when there is no Hebbian component in the learning rule ($`\eta =0`$), the net needs $`2500`$ learning steps to learn an input-output task. When we add a slight Hebbian component ($`\eta `$ small) the number of learning steps increases, and, hence, the ability of the net to learn diminishes. However, when the Hebbian component becomes more and more important, the number of learning steps starts to decrease dramatically: for $`\eta /\rho `$ between $`0.25`$ and $`0.5`$ the number of learning steps is approximately $`250`$. The Hebbian component, which has the tendency to engrave previously learned patterns, seems to help to *not* forget the old patterns. If $`\eta /\rho `$ exceeds the value $`0.5`$, learning fails. Apparently, the ‘progressive’ $`\rho `$ term, the power of which is to help the network to quickly adapt itself to new input-output relations, cannot conquer the ‘conservative’ power of the $`\eta `$-dependent (i.e., the Hebbian) term. We will come back to these points and consider the effects of the $`\eta `$ and $`\rho `$ terms in some detail in section 4.

### 3.2 Size dependences

In this section we consider the network of figure 1 for varying numbers $`N_\mathrm{I}=N_\mathrm{O}`$ and $`N_\mathrm{H}`$.

#### 3.2.1 Effect of varying the sizes of the input and output layers

In this subsection we test the performance of the network for various sizes of the input, output and hidden layers. In particular, we chose to study the size-dependence for three different values of the learning parameter: $`\eta /\rho =0`$, $`\eta /\rho =0.10`$ and $`\eta /\rho =0.45`$, values which we selected on the basis of the results of the previous subsection.

First, we take a network with the fixed number of $`512`$ neurons in the hidden layer, and only one active neuron per layer. The input and output layers will consist of increasing, equal numbers of neurons, starting with $`N_\mathrm{I}=N_\mathrm{O}=4`$. Moreover, we choose the number of input–output relations $`p`$ to be learned equal to the number of neurons in the input and output layers. The input and output layers will be enlarged in steps of $`4`$ neurons, up to $`N_\mathrm{I}=N_\mathrm{O}=28`$ neurons.

In Figure 3 we give the number of learning steps per pattern for the above mentioned three values of $`\eta /\rho `$. The positive effect of the addition of a Hebbian term to the learning rule becomes more and more convincing with increasing number of input–output relations to be learned.

#### 3.2.2 Effect of varying the size of the hidden layer

Next we consider a network with $`8`$ input neurons, $`8`$ output neurons and $`8`$ subsequent patterns. The number of active neurons is $`2`$ for all input and target patterns. We vary the number of neurons in the hidden layer.

In Figure 4 we have plotted the number of learning steps as a function of the number of neurons in the hidden layer for three values of the quotient $`\eta /\rho `$. Note that, in agreement with figure 2 the number of learning steps is the largest for $`\eta /\rho =0.1`$ (the symbols $`\mathrm{}`$ in figure 4). A suitably chosen value for the coefficient $`\eta `$ of the Hebbian term makes it possible for the network to perform satisfactory with very small number of neurons in the hidden layer (the symbols $``$ in the figure).

## 4 Explanation of the effect of the Hebbian term

The different behavior for different values of $`\eta /\rho `$ is mainly due to its consequences for the effect we call *path interference*, after Wakeling , who studied the critical behavior of the Chialvo & Bak minibrain.

As an example, let us consider the case in which only one neuron is active in each layer. In this case, the ‘path of activity’ from the active input neuron to the corresponding output neuron runs along the outgoing connections with the largest weights. During the learning process, it is possible that an established path (connecting, e.g., the active neuron of pattern $`𝝃_\mathrm{I}^1`$ with the active neuron, in the output layer, of $`𝝃_\mathrm{O}^1`$) is ‘wiped out’ by an attempt to learn one of the other input–output relations. This is likely to happen if the same neuron in the hidden layer becomes active, and, consequently, the connection to the output neuron corresponding to the previously learned pattern is punished by an amount $`\rho `$. This phenomenon of path interference will occur once in a while, irrespective of the values of the parameters $`\eta `$ and $`\rho `$. However, the question whether the memory of the old pattern is wiped out (i.e., whether the connection to the output neuron under consideration is no longer the largest), does depend on the parameters $`\rho `$ and $`\eta `$. To find out how, we should consider the change of this connection compared to the change of the other connections from the same hidden neuron to the other output neurons. For the total relative change, two different learning steps should be taken into account. Firstly, the one at $`t_p`$, at which the right output was found, and the deinforcement term did not contribute, and, secondly, the learning step at $`t_q`$, at which path interference occurred, and the deinforcement term did contribute.

Let $`w_{\mathrm{win}}`$ be the largest outgoing weight from the active neuron in the hidden layer to the output layer, and let $`w_{\mathrm{com}}`$ be a weight value which is representative for one of the other, competing weights connecting the same neuron in the hidden layer to a different neuron in the output layer.

The membrane potentials of each neurons $`i`$ of the output layer are given, according to (4), by $`h_i=w_{ij}`$, where $`j`$ is a neuron of the hidden layer. From (17), with $`x_j=1`$, we find in case of success ($`r=1`$) for the changes of the connections to the winning ($`x_i`$) and the competing ($`x_i=0`$) neurons in the output layer:

$`\mathrm{\Delta }w_{\mathrm{win}}(t_p)`$ $`=`$ $`\eta [\kappa w_{\mathrm{win}}(t_p)]`$ (18)
$`\mathrm{\Delta }w_{\mathrm{com}}(t_p)`$ $`=`$ $`\eta (\kappa +w_{\mathrm{com}}(t_p))`$ (19)
respectively. Similarly, in case of failure ($`r=0`$) these changes are

$`\mathrm{\Delta }w_{\mathrm{win}}(t_q)`$ $`=`$ $`\eta (\kappa w_{\mathrm{win}}(t_q))\rho +\phi `$ (20)
$`\mathrm{\Delta }w_{\mathrm{com}}(t_q)`$ $`=`$ $`\eta (\kappa +w_{\mathrm{com}}(t_q))+\phi `$ (21)
respectively. Only if the increase of $`w_{\mathrm{com}}`$ is larger than the increase of $`w_{\mathrm{win}}`$, the memory path can be wiped out, since then $`w_{\mathrm{com}}`$ may become the largest weight, i.e., if

$$\mathrm{\Delta }w_{\mathrm{com}}(t_p)+\mathrm{\Delta }w_{\mathrm{com}}(t_q)>\mathrm{\Delta }w_{\mathrm{win}}(t_p)+\mathrm{\Delta }w_{\mathrm{win}}(t_q)$$
(22)
We now substitute (18)–(21) into (22). In the resulting inequality we can ignore the values of $`w_{\mathrm{win}}`$ and $`w_{\mathrm{com}}`$ relative to $`\kappa `$ as long as the number of adaptations of $`w_{\mathrm{win}}`$ and $`w_{\mathrm{com}}`$ is small; recall that $`\kappa =1`$, $`\rho =0.02`$, and the initial values of the weights are in the range $`[0.01,0.01]`$ in our numerical experiments. With these approximations, the inequality reduces to $`\rho >4\eta `$. In the opposite case,

$$\rho <4\eta $$
(23)
$`w_{\mathrm{win}}`$ will remain larger than $`w_{\mathrm{com}}`$ and, consequently, path interference will not wipe out learned input–output relations, which explains the decrease of the number of learning steps for $`\eta >\frac{1}{4}\rho `$. For $`\eta >\frac{1}{2}`$ the network is incapable of learning input–output relations. This can be seen as follows. Each time a winning connection is punished, i.e., the output is wrong, it changes approximately by an amount $`\eta \rho `$, whereas the competing connection changes by an amount of $`\eta `$. Hence, only when $`\eta \rho <\eta `$, or, equivalently, when

$$2\eta <\rho $$
(24)
the winning connection decreases more than the competing connection. In the opposite case, $`2\eta >\rho `$, the winning connection remains larger than its competitor, and, at the next learning step, the output will be wrong again.

Combining the inequalities (23) and (24), we find the central result of this article (7), the parameter region for which a Hebbian term in the learning rule is advantageous. This observation is confirmed by the numerical experiment of section 3.1, so, in particular, figure 2.

Note that the reasoning leading to the main results (23) and (24) was based on an assumption regarding the initial values. In particular, it was assumed that the weights were small compared to $`\kappa `$ (which was put equal to $`1`$). In reference  it was shown that the pre-factor (16) of the Hebbian term tends to zero during the learning process:

$$\kappa (2x_i1)(h_i\theta _i)0$$
(25)
implying, that for a small number of active neurons the absolute values of the weights are of the order $`\kappa `$, as follows with (4) and (10). Hence, the assumption that the weights are small compared to $`\kappa `$ ($`\kappa =1`$) is guaranteed to break down at a certain point in the learning process.

## 5 Summary

We have shown, in a particular model, that a Hebbian component in a reinforcement rule improves the ability to learn input–output relations by a neural net.


text:
# Phenomenological implications of brane world scenarios with low tension11footnote 1Talk given by A. Dobado in the X Mexican School of Particles and Fields, Playa del Carmen, México, 2002

## 1 Are extra dimensions real or are they just science fiction?

Most of the motivations for considering extra dimension have a theoretical origin. In the last thirty years most of the new developments in theoretical high energy physics required the introduction of extra dimensions. Well known examples are modern Kaluza-Klein theories (KK), where the isometries of the extra dimensions appear as gauge symmetries, Supersymmetry (SUSY), which in the superfield formulation can be understood as a symmetry involving extra Grassmann dimensions, Supergravity (SUGRA), Superstrings defined consistently only in 10 dimensions and M-theory, which is supposed to be the eleven dimensional theory underlying the five known superstring theories and eleven dimensional SUGRA. One important exception to this rule are the Grand Unified Theories (GUT) but still the most interesting examples from the phenomenological point of view are supersymmetric since they give rise to gauge coupling unification.

The first attempts to extend general relativity to include electromagnetism date back to Theodor Kaluza and Oscar Klein KK . The discovery of 11D SUGRA produced a revival of the KK ideas in the early 80’s. The first string revolution, triggered by the realization that superstrings can provide a theory of gravitation instead of a theory of strong interactions, the cancelation of the anomalies and the discovery of the heterotic string which opened the door for phenomenology, traslated the interest to $`10`$ dimensions with six dimensional compactified spaces (Calabi-Yau, orbifolds…). The second string revolution of the 90’s introduced new ideas such as non-perturbative strings, dualities, branes and the unification of the five known superstring theories in the conjectured M-theory (these two periods of development of string theory can be reviewed in GSW  and P  respectively).

The main phenomenological problem of the old string theories is that they could not be tested since stringy effects were expected to appear at the Plank scale $`M_P`$ which is of the order of $`10^{19}`$GeV. However the new ideas coming from M-theory have inspired new scenarios that could be testable. These scenarios were developed in principle to address the hierarchy problem by putting it in a different setting. The first one was proposed by Arkani-Hamed, Dimopoulos and Dvali (ADD) ADD . The main idea is that our universe is a 3-brane living in a higher $`D=4+N`$ dimensional space (the bulk space) being the extra dimensions compactified to some small volume (Brane World scenario). In this picture the SM particles are confined to the 3-brane but gravitons can propagate along the whole bulk space. Now the fundamental scale of gravity is not the Planck scale any more but another scale $`M_D`$ which is supposed to be of the order of the electroweak scale in order to solve the hierarchy problem. Then it is possible to find the relation $`M_P^2=V_NM_D^{2+N}`$, where $`V_N`$ is the volume of the compatified extra dimension manifold $`B`$. Thus the hierarchy between the Planck and the electroweak (TeV) scale is generated by the large volume of the extra dimensions. The typical size $`R`$ of the extra dimensions ranges from a fraction of mm for $`N=2`$ to about 10 fm for $`N=6`$ (the case $`N=1`$ is already ruled out by the observations in our solar system). The most interesting property of the ADD scenario (in which we will concentrate in the rest of this work) is that it is compatible with the present experimental data, but at the same time it gives rise to many new phenomena that could be tested in the near future.

There are also scenarios where the scale of the extra dimensions $`R`$ is of the order of 1 TeV-1 A . All or some of the SM particles are allowed to propagate along the bulk. This set up is quite appropriate for model building and to deal with gauge coupling unification, SUSY breaking, neutrino spectrum, fermion masses and many other things even if it does not solve the hierarchy problem. In addition there are also models where the hierarchy is generated by the curvature of the extra dimensions as it is the case of the Randall-Sundrum model where the geometry of the space-time is $`AdS_5`$ and thus it cannot be factorized RS .

## 2 Is our Universe a 3-brane?

The most obvious consequence of the the Brane World scenario is the modification of the Newton’s Law at short distances i.e. distances of the order of $`R`$. Some recent experiments are trying to test the Newton Law’s at the submillimeter scale to explore this possibility (see for instance W ). From the point of view of particle physics the main new effects in the ADD scenario are related with the KK mode expansion of the bulk gravitational field

$$g_{\mu \nu }(x,\stackrel{}{y})=\underset{\stackrel{}{k}}{}g_{\mu \nu }^\stackrel{}{k}(x)e^{i\stackrel{}{k}.\stackrel{}{y}/R}$$
(1)
where a toroidal compactification has been assumed for simplicity, $`\stackrel{}{y}`$ represents the $`N`$ extra dimension coordinates and $`\stackrel{}{k}`$ is a $`N`$ dimensional vector with components $`k^a=0,1,2,\mathrm{}`$. Therefore a bulk graviton can be understood as a KK tower of four dimensional massive gravitons with masses of the order of $`k/R`$ with $`k`$ being any natural number (for the $`N=1`$ case) so that the distance in the mass spectrum between two consecutive KK gravitons is of the order of $`1/R`$. This means in particular that the KK graviton spectrum can be considered as approximately continuous for large extra dimensions. In principle we expect two kinds of effects from the KK graviton tower, namely graviton production and virtual effects on other particle production or observables (see for instance HS  and references therein). The rates for the different processes can be computed by linearizing the bulk gravitational field and by coupling the graviton field to the SM energy momentum tensor $`T_{SM}^{\mu \nu }`$. Then expanding the gravitational field in terms of the KK modes one finds the corresponding Feynman rules. As an example one could consider the process $`e^+e^{}\gamma +G`$ where $`G`$ stands for any KK graviton. To compute the total cross section one must sum (integrate for large extra dimensions) over all the KK gravitons. The result is

$$\sigma _{KK}\frac{1}{M_P^2}(\sqrt{s}R)^N\frac{1}{M_D^2}\left(\frac{\sqrt{s}}{M_D}\right)^N$$
(2)
Note that the total cross section is suppressed by powers of $`M_D`$ which is supposed to be of the order of 1 TeV instead of powers of $`M_P`$ which would make the cross section completely negligible. The signature of these events would be one single photon plus missing energy with a continuous spectrum. Virtual effects can be taken into account by considering the KK tower propagator

$$\underset{\stackrel{}{k}}{}\frac{1}{p^2+k^2/R^2}$$
(3)
This propagator produce divergencies for more than one extra dimension even at the tree level that require regularization by introducing some ultraviolet cutoff. This is one of facts that has given rise to the development of the so called deconstructing ACG  or aliphatic CHPW  idea where the extra dimensions are latticed.

## 3 Brane fluctuations

However there is a more physical way to deal with the above problem. So far we have asssumed that the world brane is completely rigid, i.e. it has infinite tension. But rigid objects do not exist in relativistic theories. When brane fluctuations are taken into account two new effects appear. First of all we have to introduce new fields which represent the position of the brane in the bulk space $`(x^\mu ,y^a=\pi ^a(x))`$. These fields are the Goldstone bosons (GB) corresponding to the spontaneous symmetry breaking (SSB) of the translational invariance produced by the presence of the brane (branons). When the effect of the branons is taken into account it is possible to show that the coupling of the SM particles on the brane to any bulk field is exponentially suppressed by a factor $`\mathrm{exp}(M_{KK}^2M_D^2/(8\pi ^2f^4))`$, where $`M_{KK}`$ is the KK mode mass and $`\tau =f^4`$ is the brane tension GB . This effect solves the above mentioned divergency problem. The other important consequence of this effect is that in the case where $`fM_D`$ i.e., when the brane tension scale is much smaller than the fundamental scale of gravity (in the $`D`$ dimensional space), the KK modes decouple from the SM particles on the brane. The conclusion is that for flexible enough branes the only relevant degrees of freedom at low energies in the ADD scenario are the SM particles and the branons. As other GB, branons are expected to be nearly massless 222Some branon mass $`M`$ could be expected from explicit symmetry breaking effects as it happens with pions, which are the GB related to the SSB of the chiral symmetry of strong interactions and weakly interacting at low energies, and their interactions can be described by an effective lagrangian. For homogeneous extra dimension compact space $`B`$, this effective lagragian has the structure of a non-linear sigma model where the coset space is just (is isomorphic to) the space $`B`$ DoMa . In addition the branon couple by pairs to the energy-momentum tensor of the SM, $`T_{SM}^{\mu \nu }`$ through the tensor $`_\mu \pi ^a_\nu \pi ^a`$ quite in the same way as gravitons do. Thus it is possible to obtain the different Feynman rules and amplitudes for branon production starting from SM particles in terms of the SM parameters, the brane tension scale $`f`$, the number of branons $`n`$ ($`n=N`$ for homogeneous spaces $`B`$) and the possible branon mass $`M`$ (ACDM ). The typical signature for branon production would be missing energy with continuous spectrum, but in this case the missing energy spectrum is continuous due to the fact that branons are produced in pairs. In addition by studying branon production one is probing the brane tension scale $`f`$ instead of the fundamental gravitation scale $`M_D`$ which is probed in graviton production. In the following we will consider the branon phenomenology in the ADD scenario with a flexible brane $`fM_D`$ in order to set some bounds on the brane tension and the branon mass from particle physics, astrophysics and cosmology.

## 4 Bounds from Z and W decays

In the Standard Model, the full decay width of the Z boson has three types of contributions, coming from charged leptons, hadrons and neutrinos respectively. The third one corresponds to the so called invisible width, since neutrinos escape without producing any signal in the detectors. In principle, light branons ($`M<M_Z/2`$) being stable and chargeless particles, could also give rise to invisible decays of Z and therefore to deviations with respect to the SM predictions.

The precision measurements done mostly at LEP-I set stringent limits on the Z invisible decay width, which can be translated into strong bounds on the brane tension $`f`$ and the branon mass. Thus, at the 95$`\%`$ confidence level, the variation of the Z invisible width due to decays into branons $`\mathrm{\Gamma }_Z^b(n,f,M)`$ cannot be larger than 2.0 MeV LEP .

We use the results from ACDM  for the first order contribution from branons to this width, which is given by the decay of Z into two branons and two neutrinos: $`Z\overline{\nu }\nu \pi \pi `$. Imposing the LEP-I limit on this decay width we get exclusion plots in the $`fM`$ plane (see Figure 1, left for models with a single branon). For $`n`$ massless branons, the bound on $`f`$ reads: $`f>\mathrm{\hspace{0.17em}11.9}n^{1/8}\text{GeV}`$.

Concerning the $`W^\pm `$ decay, there are mainly two new channels contributing to this process. First, the decay of $`W^\pm `$ into two branons and two leptons: $`W^{}l^{}\overline{\nu }\pi \pi `$. Such leptons can be an electron and an electron antineutrino for the $`W^{}`$ decay, their antiparticles for $`W^+`$, or the analogue pairs of leptons in the rest of families. Second, we have the $`W^\pm `$ decay into two branons and a quark-antiquark pair: $`Wq\overline{q}\pi \pi `$. In principle there should be also a channel including gluons in the final state, however in the SM, such channel is suppressed by a coefficient $`\alpha _S(M_W)/\pi =0.04`$ and for that reason we will also ignore it in the present case. In any case, the result we will obtain will be an strict lower bound to the true decay width into branons.

The bounds on the contribution from new physics to the $`W^\pm `$ width $`\mathrm{\Gamma }_W^b<240`$ MeV, are weaker than the Z case. This limit translates into the following one for massless branons: $`f>\mathrm{\hspace{0.17em}6.9}n^{1/8}\text{GeV}`$. For a single branon, we have plotted the exclusion region for the variables $`f`$ and $`M`$ in Figure 1 (left) ACDM .

## 5 Bounds from direct searches

The differential cross section of the process: $`e^+e^{}\pi \pi \gamma `$ including the two photon polarizations reads ACDM :

$`{\displaystyle \frac{d\sigma _\gamma }{dxd(\mathrm{cos}\theta )}}`$ $`=`$ $`\alpha _{EM}^2{\displaystyle \frac{2s(s(1x)4M^2)^2n}{61440f^8\pi ^2}}\sqrt{1{\displaystyle \frac{4M^2}{s(1x)}}}`$
$`\left[x(33x+2x^2)x^3\mathrm{sin}^2\theta +{\displaystyle \frac{2(1x)(1+(1x)^2)}{x\mathrm{sin}^2\theta }}\right]`$
Such expression could be used in direct searches of branons in colliders. However, we have calculated the contribution to the total cross section of processes involving a single photon in the final state. The total cross section can be calculated analytically, although we will not show it here since the final expression is quite long. However the result for massless branons (M=0) is much simpler:

$`\sigma _\gamma (M=0)`$ $`=`$ $`{\displaystyle \frac{nAs^3}{f^8}}.`$ (5)
where the constant depends on the detection limits. For the LEP-II experiments we get $`A=4.54\times 10^7`$. This has allowed us to obtain new bounds on $`f`$ assuming no observation at LEP-II: $`f138n^{1/8}\text{GeV}`$.

The results in the massive case for one branon assuming no observation is plotted in Figure 1 (right). Thus the interior area limited by the bound curve is potentially excluded by the LEP-II experiment.

Following similar steps to the single-photon case, we compute the limits from the single-Z channel ACDM . In this case the total cross section is obtained by summing over the three Z polarizations. The result for massless branons is now: $`f52n^{1/8}\text{GeV}`$, and the exclusion plot can be found in Figure 1 (right). The excluded area limited by the bound curve is smaller than in the single photon analysis. This is an expected result since the Z coupling to the electron field is smaller than that of the photon, and the Z mass restricts the avalaible phase space. In fact, this restricts the search only to branons with masses: $`M(\sqrt{s}M_Z)/257`$ GeV. In the single-photon channel however, the kinematical range is larger $`M\sqrt{s}/2103`$ GeV. In future colliders, working at very high centre-of-mass energy, the Z mass could be neglected $`\sqrt{s}M_Z`$, and both processes will allow us to get similar bounds.

## 6 Prospects for future linear colliders

Several proposals for the construction of $`e^+e^{}`$ linear colliders in the TeV range are currently under study. The TESLA (TeV Energy Superconducting Linear Accelerator), the NLC (Next Linear Collider) and the JLC (Japanese Linear Collider) are examples of the first generation of these colliders, whereas the CLIC (Compact Linear Collider) would correspond to the second generation. In this section we discuss the sensitivity of these colliders to a hypothetical branon signal.

The physics programme of the new linear collider projects includes the measurement of electroweak parameters with improved precision, such as the invisible Z width or the $`W^\pm `$ width. However, and since the deviations due to the presence of branons increase dramatically with energy, the largest sensitivity to a branon signal is expected in direct searches like single-$`\gamma `$ and single-Z.

For this purpose the LEP-II study of the previous section is extended to higher center-of-mass energies. It is assumed that, at the time of construction of these accelerators, the theoretical and systematic uncertainties on Standard Model processes will be controlled at the level of the femtobarn. Under this assumption, we have estimated the sensitivity limit by scaling the LEP-II estimate by the expected gain in statistics due to the luminosity improvement.

The critical parameter in the analysis is the center-of-mass energy, $`E_{CM}`$. In the single photon channel, this leads to limits on $`f`$ for any branon mass $`M<E_{CM}/2`$. The results of a full study in the $`(f,M)`$ plane and in different experimental contexts are presented in Figure 2. For the single-Z, the bounds are less restrictive, and the study is only applicable to branon masses below $`(E_{CM}M_Z)/2`$, due to kinematic constraints (see also Figure 2).

## 7 Bounds from astrophysics and cosmology

Since branon effects grow strongly with energy, it is not surprising that the bounds from direct searches that we have just obtained with ($`\sqrt{s}200`$ GeV) are more constraining than the indirect ones, in which the energy scale is set by the Z mass ($`M_Z90`$ GeV). In addition the bounds presented here improve the astrophysical ones for massless branons, coming from energy loss in supernova 1987A, from which $`f>10`$ GeV Kugo . For massive branons, the astrophysical limits are much less restrictive (see CDM ).

Concerning cosmology, it is interesting to notice that the allowed range of parameters suggests that branons could have weak couplings and large masses. Since they are stable particles, this makes them natural dark matter candidates. In fact, an explicit calculation shows that their relic abundance can be cosmologically relevant and could account for the fraction of one third of the total energy density of the universe in the form of dark matter presently favoured by observations.

On the other hand, the present limits on dark matter abundance can be used also to restrict the parameters space of the theory. The results are complementary to those of colliders. In the latter case, as we have seen, the unexcluded region corresponds to small cross sections (large $`f`$) and/or large branon masses, whereas having cosmological branon abundances compatible with observations requires either large cross sections (small $`f`$) or small branon masses. These results will be presented elsewhere CDM .

## 8 Summary and conclusions

Brane world scenarios are inspired on modern string (M) theory and offer new insights on many fundamental problems in particle physics. If the fundamental scale of gravity $`M_D`$ is of the order of 1 TeV, gravitons could be produced in the next generation of colliders such as the LHC. Detailed computations of the corresponding cross sections require to take into account the recoil effects of the brane. In the flexible brane limit ($`fM_D`$) the only relevant modes are the SM particles and the branons. Branon production rates and their cosmological and astrophysical implications can be determined in a model independent way in terms of the brane tension and the branon mass.

Acknowledgements: This work has been partially supported by the DGICYT (Spain) under the project numbers PB98-0782, AEN99-0305, FPA 2000-0956 and BFM2000-1326. A.D. acknowledges the hospitality of the LBNL Theory Group, where this work was partially done, economical support from the Universidad Complutense del Amo Program and very especially the organizers of the X Mexican School of Particles and Fields for their kind invitation and congratulates Augusto García and Arnulfo Zepeda for their 60th birthay.


text:
# The origin of Scorpius X-1

## 1 Introduction

Scorpius X-1 was the first discovered and it is the brightest persistent extra-solar celestial X-ray source (Giacconi et al., 1962). It is the prototype low-mass X-ray binary composed of a compact object and a donor star that have masses of 1.4 M and 0.42 M, respectively (Steeghs & Casares, 2002). The orbital period is of 18.9 hr and the compact object is accreting matter from the Roche lobe-filling companion star. Although the mass of the compact object is consistent with that of a neutron star, no thermonuclear explosion (X-ray burst of type I) of the infalling matter on the putative hard surface of a compact object has been reported so far.

The kinematical properties of Sco X-1 can help us in understanding its origin. From the radial velocity, proper motion and distance of the system, the space velocity can be derived. Thereafter, using a mass model for the Galaxy, the galactic path can be tracked. This approach was applied for first time to derive the Galactocentric orbits of the halo black hole binary XTE J1118+480 (Mirabel et al., 2001) and the Galactic disk runaway black hole system GRO J1655-40 (Mirabel et al., 2002).

Scorpius X-1 is a microquasar (for a review see Mirabel & Rodríguez, 1999), namely, a source of relativistic jets with a compact radio core associated to the accreting compact object (Fomalont et al., 2001). The most accurate proper motion of X-ray binaries can be obtained using Very Long Baseline Interferometry (VLBI) techniques at radio wavelengths. In the case of Scorpius X-1, a precise distance of 2.8$`\pm `$0.3 kpc from the trigonometric parallax, and a proper motion of the radio core have been determined using VLBI (Bradshaw et al., 1999). Furthermore, an accurate value for the systemic radial velocity was recently reported by Steeghs & Casares (2002) from spectroscopic observations at optical wavelengths. Because of the unprecedented precision in the distance, proper motion, and systemic radial velocity, Scorpius X-1 is the X-ray binary for which we have obtained the most accurate determination of the galactic path, among the seven X-ray binary systems studied so far (Mirabel & Rodrigues, 2002).

## 2 Space velocity and galactocentric orbit

The position, distance, heliocentric radial velocity, and proper motion are listed in Table 1. These data were used to compute the galactic orbit of Scorpius X-1 using the galactic gravitational potential model by Johnston et al. (1995). The velocity components U, V, and W, directed to the Galactic centre, rotation direction, and North Galactic pole were derived using Johnson & Soderblom (1987)’s equations of transformation, and assuming the sun moves with components (U,V,W) = (9, 12, 7) km s-1 relative to the local standard of rest (LSR) (Mihalas & Binney, 1981). The values of U, V, and W in Table 1 are rather different from the mean values that characterize the kinematics of stars that belong to the thin and thick disk of the Galaxy (Chiba & Beers, 2000).

The parameters and graphic representations of the galactocentric orbit are given at the end of Table 1 and in Figure 1, respectively. Because the distance, proper motion and systemic radial velocity were determined with high precision, the errors in the orbital parameters of the galactocentric orbit of Scorpius X-1 are dominated by the uncertainties in the Galactic potential model. The errors for the orbital parameters quoted in Table 1 are the deviations derived from comparing the prescriptions by Flynn et al. (1996) and de Oliveira et al. (2002) with the one adopted by us (Johnston et al., 1995). The computed orbits that result from the use of the different current mass models for the Galaxy are very similar.

Scorpius X-1 does not participate in the Galactic rotation and its kinematics is unlike those of stars that belong to the thin and thick disk populations. Presently, Scorpius X-1 is moving towards the Galactic plane with a free fall velocity of 74 km s-1, and it is receding from the Galactic Centre at a spatial speed of 125 km s-1.

## 3 Discussion

We first consider the possibility that Scorpius X-1 may be associated with Gould’s belt, which is a spur of the local spiral arm that has a size of a few hundred parsecs. Bradshaw et al. (1997) pointed out that Scorpius X-1 is beyond the Sco-Cen association. The new refined distance of 2.8 kpc by Bradshaw et al. (1999) locates the X-ray source farther beyond the edge of the belt, which is expanding. Since Scorpius X-1 is approaching along the line of sight with a velocity of 113.8 km s-1 (see Table 1), the possibility that Scorpius X-1 has been formed and ejected from Gould’s belt is ruled out.

Two possibilities remain: either Scorpius X-1 was kicked out into such an orbit by the explosion of the massive progenitor of the compact object, or it was formed in the halo itself. 30 Myr ago the binary was at a perigalactic distance of $``$500 pc with a velocity relative to the Galactic Centre of $``$480 km s-1, implying -after subtraction of Galactic rotation- a minimum kick velocity of $``$230 km s-1. This corresponds to a minimum linear momentum of $``$420 M km s-1. Although large, this minimum linear momentum does not exceeds the maximum linear momenta found among runaway solitary neutron stars and millisecond pulsars (Toscano et al., 1999). Besides, Brandt & Podsiadlowski (1995) calculate that after a supernova explosion the resulting low mass X-ray binary can remain bound provided it does not attain a velocity larger than 180$`\pm `$80 km s-1. The minimum kick velocity for Scorpius X-1 is at the limit of the maximum velocity for which the binary can remain bound. Unfortunately, from the theoretical point of view, the timescales after the explosion to achieve persistent Roche lobe overflow by the secondary and circularization of the binary orbit, are both rather uncertain. Therefore, although unlikely, one cannot rule out the possibility that Scorpius X-1 was launched from the Galactic bulge by the explosion of the progenitor of the compact object.

The galactocentric orbit of Scorpius X-1 is consistent with an origin in a globular cluster. The orbit is highly eccentric (e = 0.87), with a perigalactic distance of 0.5 kpc, and a maximum height above the plane of 4.2 kpc. It has comparable parameters to the orbits of globular clusters of the inner halo, such as NGC 6121 and NGC 6712 (Dauphole et al., 1996). However, we find that since the last perigalactic passage 30 Myrs ago, Scorpius X-1 did not intersect any known orbit of a globular cluster.

We point out that the majority of low-mass X-ray binaries of short period are found in globular clusters. Because of the high stellar densities and low velocity dispersions the cores of globular clusters are the best sites in the Galaxy for the formation by close encounters of this class of X-ray binaries. In this context, Scorpius X-1 most probably was formed in the core of a globular cluster in a close encounter of the compact object with a single star or with a binary (Verbunt, 2002), being catapulted out of the core more than 30 Myrs ago.

## 4 Conclusion

Thanks to the unprecedented precision in the distance and proper motion from VLBI measurements by Bradshaw et al. (1999), and the accurate systemic radial velocity measured with optical spectroscopy by Steeghs & Casares (2002), a precise orbital motion of Scorpius X-1 around the galactic centre was computed. The source reached a perigalactic distance of 500 pc at a velocity of 480 km s-1. Presently, it is at a distance from the Sun of 2.8 kpc, falling to the Galactic plane with a vertical speed of 74 km s-1. A possible origin associated to the Gould’s belt can now be ruled out. Scorpius X-1 has an eccentric orbit (e = 0.87) around the Galactic Centre, and the kinematics is similar to that of the most ancient stars and globular clusters of the inner halo. Most probably Scorpius X-1 has been formed in the core of a globular cluster, as the majority of the low mass X-ray binaries of short orbital period. However, we cannot rule out that a natal supernova kick launched Scorpius X-1 into an orbit like this from a birth place in the galactic bulge. In the context of these two possibilities Scorpius X-1 must have been formed more than 30 Myrs ago.

###### Acknowledgements.

We thank Jacques Paul and an anonymous referee for very useful comments on the original manuscript. I.F.M. is a member of the Consejo Nacional de Investigaciones Científicas y Técnicas of Argentina, and I.R. is a fellow of the Conselho Nacional de Desenvolvimento Cienífico e Tecnológico of Brazil. I.F.M. acknowledges support from PIP 0049/98 and Fundación Antorchas.


text:
# Comment on “Quasisaddles as relevant points of the potential energy surface in the dynamics of supercooled liquids” [J. Chem. Phys. 116, 10297 (2002); cond-mat/0203301]

Following theoretical work that emphasized the importance of saddle points on a potential energy surface (PES) in the dynamics of supercooled liquids,Cavagna01  two papers appeared proposing a mapping of instantaneous configurations onto the underlying saddle points of the PES.Angelani00 ; Broderix00  In this mapping the instantaneous configurations are starting points for minimizations on the function $`W=|V|^2`$, where $`V`$ is the potential energy. The motivation for this choice was that the stationary points of the PES, of whatever index, correspond to absolute minima of $`W`$ with $`W=0`$. However, we subsequently pointed out that the vast majority of the minima of $`W`$ do not correspond to stationary points of the PES, but have $`W>0`$.Doye02a  This result was recently confirmed by the authors of Ref. Angelani00,  after reanalysing their data.Angelani02

The minima of $`W`$ must satisfy $`W=2𝐇𝐠=\mathrm{𝟎}`$ where $`𝐇`$ is the Hessian (the second derivative matrix of the potential energy) and $`𝐠`$ is the gradient, $`V`$. There are two types of solution to this equation. First, those where $`𝐠=0`$, which correspond to true stationary points of the PES, and second, those where $`𝐠`$ is an eigenvector of the Hessian with zero eigenvalue. These minima have $`W>0`$ and have a point of inflection in the direction of the gradient on the PES, and we will continue to refer to them as non-stationary points (NSP’s) in order to differentiate them from the first set of minima.

In Ref. Angelani02,  Angelani et al. now argue that it is reasonable to consider the NSP’s as ‘quasi-saddles’: ‘If this number \[the number of inflection directions at the minima of $`W`$\] is small one can think of the minima of W as a “quasi” saddle point (QSP), as they are true saddle points in the subspace orthogonal to the small number of inflection directions. One can conjecture that the properties of a QSP are very similar to that of a true saddle point associate (sic) to an instantaneous configuration’. In this comment we suggest that these properties simply follow from the form of the solutions to $`W=0`$, and do not really indicate any saddle-like character.

We first note that NSP’s are expected to have exactly one zero Hessian eigenvalue (in addition to the three zeros corresponding to translation of the whole system), unless there is some symmetry present in the system giving rise to degenerate eigenvalues, or an unexpected accidental degeneracy, both of which are very unlikely in the present case. This analytic result is confirmed by an analysis of the 8,000 stationary points of $`W`$ of the 256-atom system in our database.Doye02a  Only eight of the 6,920 NSP’s have a non-zero Hessian eigenvalue with magnitude less than $`0.001`$. Using the convergence criteria described in Ref. Doye02a,  the average magnitude of the zero eigenvalue for the NSP’s is $`6\times 10^5`$, and the average magnitude of the smallest non-zero eigenvalue over all 8,000 points is $`0.64`$.units  We therefore suggest that the results in Ref. Angelani02, , where the average number of directions with zero curvature varies from between 1.1 and 4.4 (depending on temperature), are due to incomplete minimization, and that the tolerance used to assign the zero eigenvalues in Ref. Angelani02,  is not tight enough. Instead, the number of extra zero Hessian eigenvalues should be precisely one, independent of temperature.

Secondly, we note that at an NSP the inflection direction is parallel to the gradient. Hence the gradient in the sub-space orthogonal to the inflection direction is zero, and the NSP is a stationary point in this subspace. However, any point in configuration space is a stationary point in the subspace orthogonal to the gradient, and so this result does not imply any additional saddle-like properties.

The authors of Ref. Angelani02,  introduce a further justification for the use of the term ‘quasisaddle’ in a footnote: ‘most of the QSP found are “almost” saddles in the sense that $`W`$ is almost zero numerically’. It is not surprising that the magnitude of the gradient at the NSP’s is significantly smaller than at the instantaneous configuration, since this is what is being minimized. However, analysing our database shows that $`W`$ is not almost zero numerically. The average value of $`W`$ at a NSP varies between 0.28 and 3.25 depending on temperature, compared to a value for the true stationary points of $`2\times 10^8`$, i.e. numerically zero consistent with our convergence criteria.

We therefore disagree with the arguments used in Ref. Angelani02,  to infer saddle-like properties for NSP’s. The most interesting results of Angelani et al. do not depend upon this classification—the observation that the number of negative eigenvalues of the NSP’s extrapolates to zero at the mode-coupling temperature, and is strongly correlated with the magnitude of the diffusion coefficient, is still intriguing. However, the interpretation of these results depends upon a correct understanding of the nature of the NSP’s.


text:
# Strong magnetic field asymptotic behaviour for the fermion-induced effective energy in the presence of a magnetic flux tube

## 1 Introduction

Many authors have dealt with the numerical study of the fermion-induced effective energy, in the presence of magnetic fields of the form of a flux tube, see Refs. . An analytical study of the same problem in two dimensional Euclidean space can be found in Ref. .

A similar topic, is the study of the gluon-induced effective energy in the presence of a colored magnetic flux tube . However, in this case there is an essential difference: the effective energy has nonzero imaginary part and so the system is unstable.

It is worth mentioning the recently developed functional method in Refs. , which is applicable even for magnetic fields which do not exhibit a special kind of symmetry. For example, this method can be applied to a system of two separated magnetic flux tubes.

In Ref.  we presented an asymptotic formula (Eq. (49)) for the effective energy when the characteristic magnetic field strength $`B_m=\stackrel{}{B}(\stackrel{}{x})𝑑\stackrel{}{S}/\pi d^2`$ (or the magnetic flux $`\mathrm{\Phi }`$ of the tube) tends to infinity and the spatial size $`d`$ of the magnetic flux tube is kept fixed. However, there are some points which were not clearly explained. In fact, the arguments, which led us to the asymptotic formula, are based on a numerical study of the integral of Eq. (10), as we will see in the main part of this paper. The aim of this work is to present this study in detail. The final asymptotic formula is given by Eq. (23).

In addition, in Sec. 4 we compared the asymptotic formula of Eq. (23) with the strong field limit of the derivative expansion and find agreement.

## 2 Effective energy in the presence of a magnetic flux tube

The 3+1 dimensional fermion-induced effective energy in the presence of a magnetic field is given by the equation:

$`E_{eff(3+1)}={\displaystyle \frac{1}{2iT}}Tr\mathrm{ln}(\overline{)}D^2+m_f^2)`$ (1)
where $`\overline{)}D=\gamma ^\mu D_\mu `$ $`(\mu =0,1,2,3)`$ and $`D_\mu =_\mu ieA_\mu `$. The gamma matrices satisfy the relationship $`\{\gamma ^\mu ,\gamma ^\nu \}=2g^{\mu \nu }`$ and $`T`$ is the total length of time.

We have shown in Ref.  that a magnetic field with strength independent of $`z`$ coordinate, directed towards $`z`$-axis, has renormalized effective energy per unit length $`L_z`$222$`L_z`$ is the length of the space box towards $`z`$ direction. equal to 333For the case of massless fermions ($`m_f=0`$) the renormalized effective energy is given by appendix A in Ref. .

$`E_{eff(3+1)}^{(ren)}={\displaystyle \frac{1}{4\pi }}{\displaystyle \underset{\{n\}}{}}(E_{\{n\}}+m_f^2)\mathrm{ln}({\displaystyle \frac{E_{\{n\}}+m_f^2}{m_f^2}}){\displaystyle \frac{1}{4\pi }}{\displaystyle \underset{\{n\}}{}}E_{\{n\}}`$ (2)
where $`E_{\{n\}}`$ are the eigenvalues of the planar operator $`(\gamma ^mD_m)^2`$ $`(m=1,2)`$, and $`\{n\}`$ is a set of quantum numbers. Note that $`E_{\{n\}}0`$ since the operator $`(\gamma ^mD_m)^2`$ is positive definite as the square of a hermitian operator. For the sake of simplicity we will drop the index $`(ren)`$ in the rest of this paper.

If the magnetic field is cylindrically symmetric with finite magnetic flux, the spectrum of the operator $`(\gamma ^mD_m)^2`$ is continuous 444Besides continuous spectrum, zero modes may exist according to the Aharonov-Casher theorem .. Thus we can set $`E_{\{n\}}=k^2`$.

In order to perform the summation over $`\{n\}`$ that appears in Eq. (2), we need a density of states, given by the equation:

$$\rho (k)=\underset{l,s}{}\left(\rho _{l,s}^{(free)}(k)+\frac{1}{\pi }\frac{d\delta _{l,s}(k)}{dk}\right)=\underset{l,s}{}\rho _{l,s}^{(free)}(k)+\frac{1}{\pi }\frac{d\mathrm{\Delta }(k)}{dk}$$
(3)
where $`_{l,s}\rho _{l,s}^{(free)}(k)`$ is the density of states for free space, and $`\delta _{l,s}(k)`$ is the phase shift which corresponds to $`l^{th}`$ partial wave with momentum $`k`$ and spin $`s`$. The function $`\mathrm{\Delta }(k)`$ is defined by the equation:

$$\mathrm{\Delta }(k)=\underset{L+\mathrm{}}{lim}\underset{s,l=L}{\overset{L}{}}\delta _{l,s}(k)$$
(4)
From Eqs. (2) and (3), if we drop the field independent term $`_{l,s}\rho _{l,s}^{(free)}(k)`$ and integrate by parts, we obtain

$$E_{eff(3+1)}=\frac{1}{2\pi ^2}_0^+\mathrm{}k\mathrm{ln}\left(\frac{k^2+m_f^2}{m_f^2}\right)(\mathrm{\Delta }(k)c)𝑑k$$
(5)
where $`c=lim_{k+\mathrm{}}\mathrm{\Delta }(k)`$. Our numerical study, in Ref. , shows that $`c=\pi \varphi ^2`$ ($`\varphi =e\mathrm{\Phi }/2\pi `$). This means that $`c`$ is independent of the special form of the magnetic field we examine and depends only on the total magnetic flux of the field $`\mathrm{\Phi }`$. Also, in Ref. , we see that the function $`\mathrm{\Delta }(k)`$, when $`k+\mathrm{}`$, tends to $`c`$ fairly rapidly, and thus, the integral over $`k`$ in Eq. (5) is convergent.

Eq. (5) is suitable for the numerical computation of the 3+1 dimensional fermion-induced effective energy in the presence of a magnetic flux tube. The calculation of the phase shifts is performed by solving an ordinary differential equation. For details see Refs. .

## 3 The asymptotic behaviour for the effective energy when $`B_m+\mathrm{}`$ and $`d`$ and $`m_f`$ are kept fixed

The function $`\mathrm{\Delta }(k)`$ depends on three quantities: $`k`$, $`B_m=2\varphi /d^2`$ and $`d`$. From these we construct the following two dimensionless quantities555We have assumed the rescaling $`B_meB_m`$.: $`kB_m^{1/2}`$ and $`B_md^2`$ (=$`2\varphi `$). For dimensional reasons the function $`\mathrm{\Delta }(k)`$ can be put into the form $`\mathrm{\Delta }(k)=G(kB_m^{\frac{1}{2}},B_md^2)`$. Setting $`\mathrm{\Delta }(k)=G(kB_m^{\frac{1}{2}},B_md^2)`$ and making the change of variable $`y=kB_m^{\frac{1}{2}}`$ in Eq. (5) we obtain

$`E_{eff(3+1)}`$ $`=`$ $`{\displaystyle \frac{B_m}{2\pi ^2}}{\displaystyle _0^+\mathrm{}}y\mathrm{ln}\left({\displaystyle \frac{y^2+m_f^2/B_m}{m_f^2/B_m}}\right)\left(G(y,B_md^2)c\right)𝑑y`$ (6)
$`=`$ $`{\displaystyle \frac{B_m}{2\pi ^2}}{\displaystyle _0^+\mathrm{}}y\left(\mathrm{ln}(y^2+m_f^2/B_m)\mathrm{ln}(m_f^2/B_m)\right)\left(G(y,B_md^2)c\right)𝑑y`$ (7)
In the strong magnetic field case ($`B_m/m_f^2>>1`$) we obtain:

$`E_{eff(3+1)}`$ $``$ $`{\displaystyle \frac{B_m}{\pi ^2}}{\displaystyle _0^+\mathrm{}}y\mathrm{ln}y\left(G(y,B_md^2)c\right)𝑑y`$ (8)
$``$ $`{\displaystyle \frac{B_m}{2\pi ^2}}\mathrm{ln}(B_m/m_f^2){\displaystyle _0^+\mathrm{}}y\left(G(y,B_md^2)c\right)𝑑y`$
There are many ways to achieve a large ratio $`B_m/m_f^2`$ by changing the three independent variables $`B_m`$, $`m_f`$ and $`d`$. However, the most interesting case is when $`B_m+\mathrm{}`$ and $`m_f`$ and $`d`$ are kept fixed (see also Refs. ). This means that we keep the spatial size of the magnetic field configuration fixed and we increase the characteristic magnetic field strength $`B_m`$ or the magnetic flux $`\varphi =B_md^2/2`$.

From Eq. (8) we see that the asymptotic behaviour of the effective energy when $`B_m+\mathrm{}`$ ($`m_f`$ and $`d`$ are kept fixed) depends on the asymptotic behaviour of the integrals:

$`I_1(B_md^2)`$ $`=`$ $`{\displaystyle _0^+\mathrm{}}y\left(G(y,B_md^2)c\right)𝑑y`$ (9)
$`I_2(B_md^2)`$ $`=`$ $`{\displaystyle _0^+\mathrm{}}y\mathrm{ln}y\left(G(y,B_md^2)c\right)𝑑y`$ (10)
In the case of the integral of Eq. (9) it was possible to find an analytical expression:

$$_0^+\mathrm{}y(G(y,B_md^2)c)𝑑y=\frac{1}{12B_m}d^2\stackrel{}{x}B^2(\stackrel{}{x})$$
(11)
In order to give a proof for the above equation we will use Eq. (13) in Ref.

$$Tr(\gamma ^mD_m)^2=\underset{\{n\}}{}E_{\{n\}}=\frac{1}{6\pi }d^2\stackrel{}{x}B^2(\stackrel{}{x})$$
(12)
By using Eq. (3) for the density of states we obtain

$`{\displaystyle \underset{\{n\}}{}}E_{\{n\}}={\displaystyle _0^+\mathrm{}}k^2\rho (k)𝑑k={\displaystyle _0^+\mathrm{}}k^2\left({\displaystyle \underset{l,s}{}}\rho _{l,s}^{(free)}(k)+{\displaystyle \frac{1}{\pi }}{\displaystyle \frac{d\mathrm{\Delta }(k)}{dk}}\right)𝑑k`$ (13)
If we drop the field independent term $`_{l,s}\rho _{l,s}^{(free)}(k)`$ and set again $`\mathrm{\Delta }(k)=G(kB_m^{\frac{1}{2}},B_md^2)`$ and $`y=kB_m^{1/2}`$ we find

$`{\displaystyle \underset{\{n\}}{}}E_{\{n\}}`$ $`=`$ $`{\displaystyle _0^+\mathrm{}}k^2{\displaystyle \frac{1}{\pi }}{\displaystyle \frac{d\mathrm{\Delta }(k)}{dk}}𝑑k`$ (14)
$`=`$ $`{\displaystyle \frac{2}{\pi }}{\displaystyle _0^+\mathrm{}}k(\mathrm{\Delta }(k)c)𝑑k`$
$`=`$ $`{\displaystyle \frac{2}{\pi }}B_m{\displaystyle _0^+\mathrm{}}y(G(kB_m^{1/2},B_md^2)c)𝑑y`$
Eq. (11) is a straightforward result of Eqs. (14) and (12).

For convenience we write

$$B(\stackrel{}{x})=B_m\stackrel{~}{B}(\stackrel{}{x}/d)$$
(15)
where the magnetic field $`\stackrel{~}{B}(\stackrel{}{x})`$ has a characteristic magnetic field strength $`\stackrel{~}{B}_m=1`$ and range $`d=1`$, Eq. (11) can then be written in the form

$$_0^+\mathrm{}y(G(y,B_md^2)c)𝑑y=\frac{B_md^2}{12}d^2\stackrel{}{x}\stackrel{~}{B}^2(\stackrel{}{x})$$
(16)
From the above relationship it is obvious that the integral $`_0^+\mathrm{}y(G(y,B_md^2)c)𝑑y`$ is proportional to $`B_md^2`$.

I could not find an analytical expression, like that of Eq. (16), for the integral $`I_2(B_md^2)`$ of Eq. (10). However, a numerical study, which is presented in Fig. 1, suggests a linear 666Our numerical study can not exclude the possibility of a weak logarithmic growing of $`I_2(B_md^2)/B_md^2`$, as we will see at the end of this section asymptotic behaviour for this integral 777This implies that the asymptotic behaviour of the integral of Eq. (10) is the same with that of integral $`_0^+\mathrm{}y(G(y,B_md^2)c)𝑑y`$., when $`B_m+\mathrm{}`$, independently of the particular shape of the magnetic field. Thus we can write

$$_0^+\mathrm{}y\mathrm{ln}y\left(G(y,B_md^2)c\right)𝑑yc1B_md^2$$
(17)
where the constant $`c1`$ is independent of $`B_md^2`$, and its value depends only on the special form of the magnetic field configuration.

As we see in Fig. 1 we have performed numerical calculations for four magnetic field configurations with quite different kinds of inhomogeneity.

$`B_1(r)`$ $`=`$ $`B_m{\displaystyle \frac{1}{(r^2/d^2+1)^2}}`$ (18)
$`B_2(r)`$ $`=`$ $`B_m\mathrm{exp}({\displaystyle \frac{r^2}{d^2}})`$ (19)
$`B_3(r)`$ $`=`$ $`B_m\theta (dr)`$ (20)
$`B_4(r)`$ $`=`$ $`B_m\mathrm{\hspace{0.25em}3}(1{\displaystyle \frac{r}{d}})\theta (dr)`$ (21)
We expect that our numerical work, presented in Fig. 1, convinces the reader that the asymptotic behaviour of $`I_2(B_md^2)`$, as given by Eq. (17), does not depend on the kind of inhomogeneity of the magnetic flux tube.

Thus, when $`B_m+\mathrm{}`$ and $`m_f`$ and $`d`$ are kept fixed, from Eqs. (8),(16) and (17) we obtain

$`E_{eff(3+1)}{\displaystyle \frac{c1B_m^2d^2}{\pi ^2}}{\displaystyle \frac{B_m^2d^2}{24\pi ^2}}\mathrm{ln}(B_m/m_f^2){\displaystyle d^2\stackrel{}{x}\stackrel{~}{B}^2(\stackrel{}{x})}`$ (22)
We see that the logarithmic term dominates, so we can write

$`E_{eff(3+1)}^{asympt}={\displaystyle \frac{B_m^2d^2}{24\pi ^2}}\mathrm{ln}(B_m/m_f^2){\displaystyle d^2\stackrel{}{x}\stackrel{~}{B}^2(\stackrel{}{x})}`$ (23)
If, instead of the characteristic magnetic field strength $`B_m=\mathrm{\Phi }/\pi d^2`$ we use the parameter $`\varphi =e\mathrm{\Phi }/2\pi `$ we see that the fermion-induced effective energy, for $`\varphi +\mathrm{}`$, is proportional to $`\varphi ^2\mathrm{ln}\varphi `$.

$`E_{eff(3+1)}^{asympt}\varphi ^2\mathrm{ln}\varphi `$ (24)
This asymptotic behaviour is correct, independently of the particular shape of the magnetic flux tube, as it was shown in the above numerical study (See Fig. 1).

Note that the above asymptotic formula, of Eq. (23), is valid and for other cases of large ratio $`B_m/m_f^2`$. For example we note two: a) for $`m_f0`$ when $`B_m`$ and $`d`$ are kept fixed and b) for $`d0`$ when $`m_f`$ and $`\varphi =B_md^2/2`$ are kept fixed. For these two cases, the asymptotic formula, of Eq. (23), is obtained straightforwardly from Eqs. (8) and (16), and thus knowledge of the asymptotic behaviour of the integral $`I_2(B_md^2)`$ is not necessary.

Finally, we note that our numerical study in Fig. 1 is not sufficient to exclude a logarithmic growing of $`I_2(B_md^2)/B_md^2`$. However this logarithmic dependence, even if it really exists, should be weak compared to the logarithmic dependence which appears in the asymptotic formula of Eq. (23), otherwise it would be clearly shown in Fig. 1 and it is not.

In order to observe a possible weak logarithmic growing we ought to perform computations in an exponentially large range of values of $`B_m`$. We could not execute computations for so wide a range, as the differential equation for the computation of the phase shifts (Eq. (69) in Ref. ) is stiff for large $`B_m`$. However, we tried an estimate of this logarithmic dependence by fitting a curve of the form $`c_1+c_2\mathrm{ln}(B_md^2)`$ to our data, as they are presented in Fig. 1. We find that $`c_2=0.0047\pm 0.0005`$ for $`B_1`$, $`c_2=0.0100\pm 0.0005`$ for $`B_2`$, $`c_2=0.0062\pm 0.0005`$ for $`B_3`$ and $`c_2=0.0120\pm 0.0006`$ for $`B_4`$. These results imply that, even if this logarithmic term $`c_2\mathrm{ln}(B_md^2)`$ really exists, it would not contribute to the effective energy (as it is given by the asymptotic formula of Eq. (23) above) more than 5 per cent.

## 4 Derivative expansion for $`B_m+\mathrm{}`$ when $`d`$ and $`m_f`$ are kept fixed

An approximate way for the computation of the fermion-induced effective energy is the derivative expansion . This method gives accurate results for smooth magnetic field configurations. Note that the derivative expansion fails for the magnetic fields of Eqs. (20) and (21), as the magnetic field of Eq. (20) is discontinuous and the magnetic field of Eq. (21) has discontinuous first derivative.

Our aim is to compare the strong field limit of the derivative expansion with the asymptotic formula of Eq. (23).

If we keep the first two terms of the derivative expansion, the effective energy, in the case of 3+1 dimensions, is given by the equation

$$E_{eff(3+1)}^{(der)}=E_{(3+1)}^{(0)}[B]+E_{(3+1)}^{(1)}[B,(B)^2]$$
(25)
where

$`E_{(3+1)}^{(0)}={\displaystyle d^2\stackrel{}{x}\frac{B^2(\stackrel{}{x})}{8\pi ^2}_0^+\mathrm{}\frac{1}{s^2}(\mathrm{coth}(s)\frac{1}{s}\frac{s}{3})e^{sm_f^2/B(\stackrel{}{x})}𝑑s}`$ (26)
$`E_{(3+1)}^{(1)}=({\displaystyle \frac{1}{8\pi }})^2{\displaystyle d^2\stackrel{}{x}\frac{(B)^2}{B(\stackrel{}{x})}_0^+\mathrm{}\frac{1}{s}(\frac{d}{ds})^3[s\mathrm{coth}(s)]e^{sm_f^2/B(\stackrel{}{x})}𝑑s}`$ (27)
For $`B_m+\mathrm{}`$ from Eqs. (26) and (27) we obtain

$`E_{(3+1)}^{(0)}{\displaystyle \frac{B_m^2d^2}{24\pi ^2}}\mathrm{ln}(B_m/m_f^2){\displaystyle d^2\stackrel{}{x}\stackrel{~}{B}^2(\stackrel{}{x})}`$ (28)
$`E_{(3+1)}^{(1)}B_m({\displaystyle \frac{1}{8\pi }})^224\zeta ^{}(2){\displaystyle d^2\stackrel{}{x}\frac{(\stackrel{~}{B})^2}{\stackrel{~}{B}(\stackrel{}{x})}}`$ (29)
Eq. (29) can be derived straight from Eq. (27) if we assume that $`B_m>>m_f^2`$ (for details see Ref. ). In order to derive Eq. (28) from Eq. (26) we have used a similar way with that presented in Ref. .

Thus Eq. (28), which gives the strong field limit of the derivative expansion, is in agreement with the asymptotic formula of Eq. (23).

## 5 Conclusions

We emphasize that the arguments, which led us to the asymptotic formula of Eq. (23), are based on a numerical study. In particular, we computed numerically $`I_2(B_md^2)`$ for four magnetic field configurations with quite different kinds of inhomogeneity (see Eqs. (18), (19), (20) and (21)). Our results in Fig. 1 suggest that $`I_2(B_md^2)`$ is proportional to $`B_md^2`$ (for $`B_md^2>>1`$) independently of the specific form of the magnetic field configurations. The consequence of this result is that also the final formula, which states that $`E_{eff(3+1)}^{asympt}\varphi ^2\mathrm{ln}\varphi `$ for $`\varphi +\mathrm{}`$, is valid independently of the specific form of the magnetic field.

Finally we showed, comparing with the asymptotic formula of Eq. (23), that the leading term of the derivative expansion gives the same asymptotic behaviour for the fermion-induced effective energy in the strong magnetic field regime even in the case (like that of the magnetic fields of Eqs. (20) and $`(21)`$) when the derivative expansion is invalid, because the next to the leading term diverges.

## 6 Acknowledgements

I am grateful to Professor G. Tiktopoulos for numerous valuable discussions, as well as Professors K. Farakos and G. Koutsoumbas for useful comments and suggestions.


text:
# Mott transition and transport crossovers in the organic compound 𝜅-(BEDT-TTF)2Cu[N(CN)2]Cl

## Abstract

We have performed in-plane transport measurements on the two-dimensional organic salt $`\kappa `$-(BEDT-TTF)2Cu\[N(CN)2\]Cl. A variable (gas) pressure technique allows for a detailed study of the changes in conductivity through the insulator-to-metal transition. We identify four different transport regimes as a function of pressure and temperature (corresponding to insulating, semi-conducting, ”bad metal”, and strongly correlated Fermi liquid behaviours). Marked hysteresis is found in the transition region, which displays complex physics that we attribute to strong spatial inhomogeneities. Away from the critical region, good agreement is found with a dynamical mean-field calculation of transport properties using the numerical renormalization group technique.

The Mott metal-insulator transition (MIT) is a key phenomenon in the physics of strongly correlated electron materials. It has been the subject of extensive experimental studies in transition metal oxides, such as (V1-xCrx)2O3 or chalcogenides such as NiS2-xSex (for a review, see Ref. Imada et al. (1998)). In contrast to chemical composition, hydrostatic pressure allows in principle to sweep continuously through the transition. For these materials however, the appropriate range of pressure is several kilobars. For this reason, many fundamental aspects of the MIT are yet to be studied in detail. This issue is particularly important in view of recent theoretical predictions for e.g. spectroscopy and transport close to the transition, which should be put to experimental test Georges et al. (1996). Layered charge-transfer salts of the $`\kappa `$-(BEDT-TTF)$`{}_{2}{}^{}\mathrm{X}`$ family (where X is a monoanion) offer a remarkable opportunity for such a study. Indeed, these compounds are known to display a great sensitivity to hydrostatic pressure Ito et al. (1996). The $`\kappa `$-(BEDT-TTF)2Cu\[N(CN)2\]Cl compound in particular (abbreviated $`\kappa `$-Cl below) displays a very rich phase diagram with paramagnetic insulating, antiferromagnetic insulating, superconducting and metallic phases when pressure is varied over a range of a few hundred bars McKenzie (1997); Lefebvre et al. (2000); Ito et al. (1996).

In this work, we report on an extensive experimental study of the in-plane resistivity of the $`\kappa `$-Cl compound for a range of pressure spanning both the insulating and the metallic phases. In contrast to previous studies Ito et al. (1996), pressure is varied continuously using a Helium gas cell, at constant temperature. By analyzing the pressure and temperature dependence of the measured resistivity, we have identified important crossover lines, which are summarized on the phase diagram of Fig.1. These crossovers separate four different regimes of transport (to be described below) within the paramagnetic phase, corresponding to an ”insulator”, a ”semiconductor”, a ”bad metal” and a Fermi-liquid metallic regime.

Our measurements were performed both with increasing and decreasing pressure sweeps, yielding a determination of the two spinodal lines $`P_1^c(T)`$ and $`P_2^c(T)`$ shown in Fig. 1 (see also Fig. 3). A critical comparison is made to dynamical mean-field theory (DMFT Georges et al. (1996)), which describes successfully the observed transport crossovers, while the critical region itself reveals more complex behaviour. Numerical Renormalization Group (NRG) calculations of the resistivity within DMFT are reported here for the first time, and compare favorably to the experimental data.

This paper is based mainly on the data set displayed on Fig. 2. The in-plane conductivity was measured using a standard four- terminal method as a function of pressure between 1 bar and 1 kilobar. The temperature range 13K-52K was investigated, with pressure sweeps performed every Kelvin (for clarity, only selected temperatures are displayed on Fig. 2). The data clearly demonstrate how pressure drives the system from a low-conductivity insulating regime at low pressure (below 200 bar) to a high-conductivity metallic regime at high pressure (above 350 bar), with a more complex transition region in the 200-350 bar range.

The pressure interval in which the normalized difference between these two measurements exceeds our experimental precision (inset of Fig. 3) corresponds to the region of coexistence between the insulating and metallic phases. This interval is found to be rather large ($`100`$ bar) at $`35`$ K, and still measurable (although the signal is close to our error bars) at $`39`$ K. Note that a much narrower interval would be found from the width at half-maximum. We emphasize that varying pressure rather than temperature is a much more accurate experimental technique for the investigation of the coexistence region, because of the shape of the spinodal lines in the $`(P,T)`$ plane. This also ensures a better accuracy on pressure (of order 1 bar) as compared to a cooling experiment. Our measurements yield a determination of the critical point in the $`T_c39`$ K to $`T_c42`$ K range ($`P_c280`$ bar), somewhat higher than previously reported from NMR Lefebvre et al. (2000) or sound velocity Fournier et al. (2002) experiments. The data in Fig. 3 are far from the ”ideal” situation of a sharp discontinuity in $`\sigma (p)`$ for $`T<T_c`$, and present rather a wide region of phase coexistence, in which the measured $`\sigma (P)`$ displays sizeable noise. This suggests a complex physics in this regime, with the likely appearance of insulating and metallic domains in the material, and pinning by disorder.

We now focus on the low-pressure, insulating regime. In the inset of Fig. 4, the conductivity data in this pressure range have been replotted as $`\mathrm{log}\rho `$ vs $`1/T`$. At lower temperature (below 50 K), an activation law $`\rho \mathrm{exp}(\mathrm{\Delta }/2T)`$ describes the data quite well. The temperature range over which this activation law applies is limited both from below and from above. Indeed, at the lowest temperatures a transition into the antiferromagnetic insulating phase is found (around 25 K at ambient pressure, see Fig. 1), while a crossover to a different insulating regime is observed at higher temperatures, as discussed below. Despite this limited range, a reasonable determination of the insulating gap $`\mathrm{\Delta }`$ in the ”paramagnetic insulator” region of Fig. 1 can be achieved, since $`\mathrm{\Delta }T`$. An approximately linear pressure dependence (see Fig. 4) is found: $`\mathrm{\Delta }7402P_{\mathrm{bar}}`$ (in Kelvin). Hence, the activation gap is still a large scale close to the coexistence region ($`\mathrm{\Delta }400\text{K}`$ at $`150`$ bar), of the order of ten times the critical temperature associated with the transition. $`\mathrm{\Delta }`$ cannot be reliably determined above $`150`$ bar, but a rough extrapolation would lead to a gap closure around $`370`$ bar, on the high-pressure side of the coexistence region of Fig. 1, suggesting that the finite-temperature Mott transition in this material is not driven by the closure of the Mott gap. At higher temperature (above $`T50`$ K at ambiant pressure), the $`T`$-dependence of the resistivity deviates from the above activation law (inset of Fig. 4). This crossover into a ”semi-conducting” regime, which extends over much of the high-T part of the phase diagram, is indicated on Fig. 1. Interestingly, the crossover scale is an order of magnitude smaller than the insulating gap.

The high pressure regime above $`300`$ bar is now considered. As demonstrated in the right inset of Fig. 5, the resistivity in this regime has a quadratic, Fermi liquid dependence upon temperature at low temperatures: $`\rho =\rho _0(P)+A(P)T^2`$. The quality of the fit, obtained by replotting the data set of Fig.2 as a function of $`T^2`$, illustrates the high precision of the variable pressure technique.

The coherence temperature $`T^{}`$ above which this law is no longer valid (of the order of $`35`$ K at $`500`$ bar) defines a crossover into a ”bad metal” regime, as indicated on Fig. 1. The prefactor $`A(P)`$ of the $`T^2`$ dependence is found to depend strongly on pressure, as displayed on Fig. 5, and the product $`A(P)\times (T^{})^2`$ remains approximately pressure-independent. These findings correspond to a strongly correlated Fermi liquid regime at low temperature. We note that the residual resistivity has a weaker pressure-dependence than A, but does increase close to the coexistence region. While $`A(P)`$ cannot be determined precisely below $`280`$ bar, the data suggest a divergency of $`A`$ at $`P220`$ bar, significantly smaller than the pressure at which the extrapolated insulating gap would vanish (of order $`370`$ bar, see above). This suggests that the closure of the Mott-Hubbard gap and the loss of Fermi liquid coherence are two distinct phenomena, associated with very different energy scales, as is also clear from the fact that the coherence scale $`T^{}`$ (a few tens of Kelvin) is much smaller than the insulating gap $`\mathrm{\Delta }`$ (several hundreds Kelvin).

In order to better characterize the crossover into the ”bad metal” as temperature is increased above $`T^{}`$, we have performed measurements in a wider temperature range, up to $`300`$ K, as displayed on Fig. 6.

We confirm the general behaviour reported by other authors Ito et al. (1996). At moderate pressures (a few hundred bars), the resistivity changes from a $`T^2`$ behaviour below $`T^{}`$, to a regime characterized by very large values of the resistivity (exceeding the Ioffe-Regel-Mott limit by more than an order of magnitude) but still metallic-like ($`\mathrm{d}\rho /\mathrm{d}T>0`$). For pressures in the transition region, this increase persists until a maximum is reached, beyond which a ”semi-conducting” regime is entered ($`\mathrm{d}\rho /\mathrm{d}T<0`$). This regime is continuously connected to that found on the insulating side, at temperatures above $`50`$ K. Fig. 1 summarizes the four regimes of transport which we have described so far.

We now compare these experimental findings to the DMFT description of the Mott transition. Similarities between DMFT and some physical properties of BEDT organics have been emphasized previously McKenzie (1998); Merino and McKenzie (2000). One of the key qualitative outcomes of our measurements is the separation of energy scales found in the transition region. This is also a distinctive feature of DMFT, in which Hubbard bands are well separated from the quasi-particle peak in the correlated metal, and a coexistence region of insulating and metallic solutions delimited by two spinodal lines is found. The three transport regimes found on the metallic side of the transition are well accounted for within DMFT Rozenberg et al. (1994); Majumdar and Krishnamurthy (1995); Georges et al. (1996); Merino and McKenzie (2000). In this theory, a Fermi- liquid regime with long lived quasiparticles applies below the coherence scale $`T^{}`$, with $`\rho =AT^2`$ and $`A1/(T^{})^2`$. As $`T`$ approaches $`T^{}`$, the lifetime decreases sharply while the low-energy local density of states (d.o.s) is strongly depleted. This leads to a large resistivity, increasing rapidly with $`T`$, characteristic of a bad metal. For $`TT^{}`$, the local d.o.s displays a pseudogap and a ”semi-conducting” transport regime is found, with $`d\rho /dT<0`$. We have performed DMFT calculations of the resistivity for a simple Hubbard model at half-filling. The NRG technique was used Bulla et al. (2001), which is crucial to capture correctly the enhancement of $`A`$, and to determine accurately the resistivity maximum. The coupling was fixed at $`U=4000`$ K, consistent with the measured critical temperature (within DMFT, $`T_cU/100`$). A semi-circular band was used, with a bandwidth $`W`$ adjusted for each pressure, as indicated on Fig. 6. Comparison to the measured $`\rho (T)\rho _0`$ involves also a global scale factor (the same for all curves). We find a good agreement between our calculations and the data, up to $`T150`$ K, for pressures ranging from the transition (300 bar) to $`10`$ kbar. The fitted bandwidth increases by about $`40\%`$ in that range. This value, together with the magnitudes of $`U`$ and $`W`$ is consistent with theoretical Rahal et al. (1997) and experimental Klehe et al. (2000) estimates for these materials. We attribute the discrepancy observed above $`150`$ K to the strong thermal dilatation of the materials, leading to a smaller effective bandwidth, an effect which could be taken into account by allowing for a T-dependence of $`W`$ in the calculation. The same effect might also contribute to the observed crossover on the insulating side at $`T50`$ K, even though DMFT also leads to a similar prediction of a purely electronic crossover at $`T\mathrm{\Delta }`$ Rozenberg et al. (1994); Georges et al. (1996).

Finally, we observe that our data in the critical regime around $`(T_c,P_c)`$ do not appear to obey the critical behaviour established within DMFT, namely that of a liquid-gas transition in the Ising universality class Kotliar et al. (2000). The observed transition is much too smooth, and the coexistence region too wide, to be described in that manner. Besides the low-dimensionality of these compounds, this could be attributed to strong spatial inhomogeneities or disorder in the transition region.

To conclude, we have performed detailed transport measurements on the $`\kappa `$-Cl compound using a variable pressure technique. Four transport regimes have been identified as a function of temperature and pressure, separated by crossover lines. Hysteresis experiments have confirmed the first-order nature of the transition and allowed for the first detailed determination of the coexistence region. DMFT provides a good qualitative description of these crossovers, and new DMFT-NRG calculations of transport are in good agreement with the data. Further theoretical work is needed to understand the critical regime, however. We note that recent ultrasound velocity measurements Fournier et al. (2002) have revealed anomalies along crossover lines very similar to the two high-temperature lines reported here from transport measurements. This calls for further theoretical work aiming at connecting these two effects.

###### Acknowledgements.

We are grateful to R. Chitra, V. Dobrosavlevic, G. Kotliar, R. McKenzie, M. Poirier, M. Rozenberg and A.M. Tremblay for useful discussions. T.A.C., S.F. and A.G. are grateful to KITP-UCSB for hospitality during the final stages of this work (under NSF grant PHY99-07949). T.A.C. acknowledges support from the SFB 195 of the DFG.