text:
# Aging and effective delays in colloidal glasses under shear

## I Introduction

It has been known for long that many physical systems such as molecular liquids, metallic liquids, polymer melts when cooled at low temperature exhibit a dynamical transition without structural changes. Their typical relaxation time increases by more than 10 decades and can no longer be measured when the temperature is low enough. This phenomenon is known as the glass transition. When the typical relaxation times becomes superior to any reasonable experimental time scale, the systems displays a non-equilibrium behavior. This out-of-equilibrium behavior is revealed by a continuous slowing down of the microscopic dynamics of the system. This phenomena is referred as aging.

Aging processes have been extensively studied in the past decades. The basic experiments to test aging consist in quenching the system from the high temperature phase into the glassy phase and then measuring the time evolution of its physical properties such as response or correlation functions. The non equilibrium behaviour appears in the dependance of theses functions with the time $`t_w`$ that the system had spend in the glassy phase. For these simple thermal histories $`t_w`$ is called the age of the system. The slowing down of the dynamics following such a rapid quench is coined as simple aging (see eg:Struik (1978); Vincent et al. (1996)). It often happens in a self similar manner such that the responses or correlation functions can be rescaled on a master curve when plotted against the reduced variable $`\frac{t}{t_w^\mu }`$. $`\mu `$ is then called the aging exponentBouchaud et al. (1998). This past decades a great number of experiments showed the universality of this kind of behavior for an increasing variety of systemsKityk et al. (2002); Angel et al. (2000); Kircher and Bohmer (2002).
The glassy state can also be reached by increasing the density. This can be realized on colloidal glasses van Megen and Underwood (1994) that also exhibit a glass transition upon an increase of volume fraction Knaebel et al. (2000); Ramos and Cipelletti (2001). When rapidly driven into their glassy region, this systems also display classical aging with an aging exponent usually inferior to 1. However temperature and density are not practical control parameters for colloidal glasses. It has nonetheless been observed Derec et al. (2000); Cloitre et al. (2000) that theses type of systems can be fluidized under mechanical sollicitation (oscillatory shear strain or constant shear rate). Moreover, it has been recently shown that the slowing down of the dynamics following the stop of the mechanical drive is that of a simple aging experiment. Hence the cessation of high shear is often referred as a rheological quench.

Simple aging experiments are the simplest protocoles that can be implemented on glassy systems. However, more elaborate temperature histories gave deeper insight of the organisation of the microscopic dynamics of theses amorphous materials. A common temperature history that has been tested on this systems is a transient temperature shift. After a quench into the glassy phase, the system is transiently driven from a non equilibrium state into another non equilibrium state. This protocole implemented on spin Dupuis et al. (2001), molecular Kityk et al. (2002), polymer Bellon et al. (2002), ferroelectric Doussineau and Levelut (2002) glasses cast into light new phenomena such as rejuvenation, memory effects or Kovacs effects Kovacs and al (1979). The two first effects revealed a very peculiar organisation of the phase space of theses systems (hierarchical organization Lederman et al. (1991) or chaotic behavior Fisher and Huse (1986); Bray and Moore (1987)). The last effect has been recently interpreted as a non monotonic change of the dynamics at all length scales Berthier and Bouchaud (2002). In this paper we transpose this more elaborate temperature protocole to a colloidal glass using transient mechanical shearing periods instead of transient temperature shifts. The observed phenomenology following this two types of perturbations can thus be compared.

We use MultiSpeckle Diffusing Wave Spectroscopy (MSDWS) to probe the transient change of the microscopic dynamics upon shear application. In the first part of this paper we will recall the technique. We will secondly confirm that colloidal glasses can be driven into a steady reproducible state by high oscillatory shear. We also found that the systems has a full aging behavior after the shear cessation. We then show that the system can keep aging under low enough shear. We thus determine the shear amplitude that drive the system into a steady regime. Thirdly, we show the change of the dynamics after a transient mechanical perturbation can be described by three regimes. High amplitudes perturbations lead to a rejuvenation of the system whereas low amplitude ones help the system to age (overaging). Intermediate amplitudes lead to a dynamical state where both signs of overaging and rejuvenation are simultaneously present depending on the considered time scales. We discuss a description of the two first regimes in terms of effective delay. We then explain why effective delay is an insufficient notion to describe our result and propose a model independent interpretation of our data in terms of changes of shape of the relaxation time distribution for the system. We show that our results are consistent with the existence of a weak long term memory and that they do not give any hint of a hierarchical organisation of the relaxation times. We eventually propose a spatial scenario to explain the observed phenomenology.

## II Experimental methods

In this section we describe the sample we used as well as the MSDWS technique.

### II.1 Sample preparation

The sample is a commercial suspension of polystyrene spherical beads of diameter 162 nm copolymerized with acrylic acid (1%) that creates a charged corona stabilizing the microspheres. The corona prevents from aggregation by steric and electrostatic repulsions. A concentrated suspension of theses colloidal particules is prepared by dialysis from a dilute suspension. Great care is taken that no air bubble is introduced in the dense suspension during sample preparation. The volume fraction of the sample is fixed at $`\varphi =50\%\pm 0.5\%`$ and was determined by drying. For that volume fraction the sample is pasty and displays a weak yield stress. No sign of crystallization could be noticed over several month. The volume fraction was chosen to fulfill two criteria:

* The sample has to be concentrated enough so that the structural relaxation can be observed ($`\alpha `$ modes). The $`\alpha `$ modes must also display a measurable aging behavior over all the experiment time scale.
* The volume fraction has to be low enough so that the strain is homogeneous, at least on macroscopic length scale all through the sample. We used a high resolution camera and a x20 zoom to check the shear strain was homogeneous at a macroscopic level over the whole sample. Effect like shear banding or wall slippage were not detected.

The sample is then carefully loaded through a syringe into a custom made shear cell. The cell is sealed with silicon oil which viscosity is chosen so that it forms a stable meniscus around the cell for several months.

### II.2 Diffusing wave spectroscopy set-up

Since we extensively described the Multispeckle Diffusing Wave Spectroscopy (MSDWS) in a previous article Viasnoff et al. (2002), we will simply recall the principle of this technique as well as our experimental set-up.
MSDWS is a dynamical light scattering technique that enables the monitoring of particles motions in a turbid non ergodic media. We use an Argon-ion laser (Spectra Physics BeamLok 2Watts). The beam is expanded to a diameter of 1 cm and is incident on the plate-plate shear cell with a 3mm gap. Multiply scattered light is collected by a 50-mm Nikon camera lens. It is focused onto an iris diaphragm. The lens is set up so that a one-to-one image of the scattered light from the output plane of the sample cell is formed on the diaphragm. The scattered light emerging from the diaphragm is detected by a Dalsa CCD camera, model CAD1-128A, which is placed approximately 15 cm behind the diaphragm. It is an 8-bit camera which we run at 500 frames per second. The images are transferred to a computer by a National Instruments data acquisition board, model PCI1422. Figure 1 shows the experimental set-up.

The interference state of the light collected onto the CCD chip is a random pattern called a speckle pattern. We analyse the temporal fluctuations of this pattern via the two times intensity autocorrelation :

$$g_2(t_w,t)=\frac{I_i(t_w)I_i(t+t_w)_i}{I_i(t_w)_iI_i(t+t_w)_i}$$

where $`I_i(t_w)`$ is the intensity recorded on pixel $`i`$ at time $`t_w`$, and $`\mathrm{}_i`$ is the average over all the pixels of the considered frames. The intensity autocorrelation function can be related to the normalized electric field autocorrelation function $`g_1(t_w,t)`$ via the Siegert relation:

$$|g_1(t_w,t)|^2=\frac{1}{\beta }(1g_2(t_w,t))$$

where $`\beta `$ is a normalizing factor. This technique does not require any time average. It thus allows a precise study of non ergodic systems in transient regimes (see ref Viasnoff et al. (2002)). Although the direct extraction of the scatters motion from $`|g_1(t_w,t)|^2`$is far from being established in the case of dense sample, we will show that interesting informations can be extracted from the changes of correlation function shape.

## III Experimental results

In this section, we present the different dynamical behavior displayed by our colloidal glass when submitted to an oscillatory shear. We first show that a complete rejuvenation of the system can be achieved under a high shear strain. After the shear cessation the colloidal glass displays a generic aging behavior. We then demonstrate that the aging process survives to a low amplitude shear. Under such conditions, not only the system keeps on aging but the aging rate is made faster by the shear application. The dynamics thus appears older. We call this phenomenon: ”overaging”. For intermediate shear amplitude, the system presents both behavior depending on the different probed time scales.

### III.1 Complete rejuvenation and simple aging

We first establish that the microscopic dynamics of the system can be driven to a steady state by the application of a high shear. Whatever the past history of the sample, its microscopic dynamic can resume an equilibrium behaviour in a finite amount of time if the sample undergo a mechanical shearing of high amplitude. Indeed, as presented on figure 2 (a), we submitted our sample to various histories. They consist in a series of both different periods of shearing at various rates and amplitudes and of periods of rest for various duration. The sample is then submitted to an oscillatory shear of frequency f=1Hz and amplitude A=30% during D=100s. After the shear cessation the sample is left at rest and its microscopic dynamics is probed by measuring the two times intensity autocorrelation functions $`|g_1(t_w,t)|^2`$. The results are presented on figure 2 (b).

Whatever the history before the oscillations may be, the set of curves obtained after the shear cessation is always the same. Further more these curves remain unchanged, within our experimental precision, when the duration D varies from 40 sec to 10000 sec, and/or when the amplitude $`A`$ ranges within 20% and 35% (that corresponds to the upper limit of our shear cell). Hence applying an oscillatory strain for 100s at 30% erases all the memory of the past history of the sample. It is then completely rejuvenated. As demonstrated on figure 2(b), the system starts aging when the shear is stopped. The auto correlation functions depend both on t, the time elapsed since the beginning of the measurement, and $`t_w`$, the time elapsed since the shear cessation. All theses curves can be rescaled and their characteristic decorrelation times $`\tau _\alpha `$ (define as the times needed for the functions to decorrelate by half their initial values) vary like $`t_w^{1.06\pm 0.08}`$ as demonstrated in Viasnoff and Lequeux (2002). The cessation of the high shear can thus be view as a rheological quench.

The system is now submitted to a shear strain of amplitude 30%, 1Hz for 100s. The amplitude is then lowered to 3%, the frequency is changed to 0.1Hz and the correlation functions are recorded while the shear is On. This history is sketched on figure 3(a). Whereas previously, the decorrelation was due only to thermally driven motions of the scatterers, we now monitor together thermal and forced motions of the latex spheres due to the external shear. The characteristic shape of $`|g_1(t_w,t)|^2`$ is given by fig 3(b). The curves displays a pseudo-periodic structure. For each period of shear the latex particles go back close to their reference position they occupied at time $`t_w`$. This explains why the correlation fonctions rises back for each period. This shape is usually called an echoHebraud et al. (1997). However, the thermal and shear induced activations lead to a certain degree of irreversibility in the scatterers motion Hebraud et al. (1997); R.Höhler (1997); Petekidis et al. (2003). This is reflected by a change in the height of the echo. A complete description of the technique can be found in Hebraud et al. (1997).

Figure 3(b) shows that despite the application of the shear the systems keeps on aging over a long period of time. Indeed the height of the echoes keeps evolving even after 1000s. This shows that a 5.5% shear is not able to drive the system out of its glassy phase, since no equilibrium could be found over the time range 1s-1000s. We do not claim, that there is no equilibrium state that could be reach after a longer time spent under shear; however, for the time range we studied our system presents a typical aging behavior. We would now like to evaluate the amplitude above which the system can reach a steady state. A straightforward way to do so would be to increase the amplitude of the shear and look when the dynamics under shear reaches a steady state. This amplitude could then be considered as the limit of the glassy region. However the echo technique fails before such an amplitude can be found. Indeed for too high shear amplitude the echo technic suffers from sampling problems: the camera acquisition rate being too small, the height of the echoes is then measured via one single data point whose value do not necessarily corresponds to the echo’s real height. Another scheme has to be found.

We proceed as described on figure 4(a). The sample is fully rejuvenated by a high shear. On technical grounds the shear is stopped for 10 s after which the oscillation is switched back on at a frequency of 1Hz for different amplitudes $`a`$ for various durations $`d`$. The oscillation is then stopped and an intensity correlation function is measured 0.1s after the shear cessation. Its half decorrelation time $`\tau _\alpha (0.1s)`$ is then evaluated as a fonction of both $`a`$ and $`d`$. The results are displayed on figure 4 (b).

It clearly shows that for amplitudes inferior to 7.9%, $`\tau _\alpha (0.1s)`$ does not reach any steady state for the duration we studied. However for amplitudes superior to 8.9% a steady state is reached and $`\tau _\alpha (0.1s)`$ does not depend any more on $`d`$. In addition, $`\tau _\alpha (0.1s)`$ does dependent anymore on $`a`$ for amplitudes superior to 20% as explained in the first paragraph. The limit of the glassy phase lies thus between 7.9% and 8.9%. Actually, for amplitudes lying between these two values the determination of $`\tau _\alpha (0.1s)`$ is made difficult by large fluctuations of its values from experiment to experiment.

### III.2 high shear behavior

We evaluated the boundary of the glassy phase by considering only one point ($`\tau _\alpha (0.1s)`$) of the first correlation function following the shear cessation. In this section, we now focus on the full set of curves that derive from this shear history. The duration of the second burst is now fixed to 100 seconds. We concentrate on the dynamics following the shear cessation for different amplitudes. As in the previous section the system is fully rejuvenated by a high shear, then it is left at rest for 10s, before it undergoes another burst of oscillation at 1Hz during 100s for an amplitude of 14.9%. Hence the second burst drives the system out of its glassy region. The set of intensity correlation functions is recorded after the low intensity burst cessation. We call $`t_w`$ the time elapsed since the cessation of the high shear (this notation is changed compared to previously published papers). We compare the shape of the curves after the shear $`|g_1(t_w,t,\gamma =14.9\%)|^2`$ to the curves obtained with no second burst: $`|g_1(t_w,t,\gamma =0\%)|^2`$. As observed on figure 5(a), it appears that $`|g_1(t_w,t,\gamma =14.9\%)|^2`$ has the same general shape as $`|g_1(t_w,t,\gamma =0\%)|^2`$. However, for the same $`t_w`$, $`|g_1(t_w,t,\gamma =14.9\%)|^2`$ decorrelates more quickly than $`|g_1(t_w,t,\gamma =0\%)|^2`$. Hence the set of sheared curves seems younger than the unperturbed one. This result is not surprising and is a direct consequence of the fluidisation of the system under shear. The analysis of the curves can be pushed further.
Since the global shape of the two sets of curves seems similar it is interesting to see if one can be deduced from the other by a simple shift in time. Thus we computed the value of the effective delay defined as:

$$|g_1(t_w+t_{eff},t,\gamma =0\%)|^2|g_1(t_w,t,\gamma =14.9\%)|^2$$

for different $`t_w`$. For all the $`t_w`$ that we studied a negative value of $`t_{eff}`$ could be found that led to a good overlay of the perturbed and unperturbed curves (see figure 5(b)). We found a unique effective delay $`t_{eff}=70s\pm 5s`$ for all $`t_w`$.

Notice also that the perturbed and reference curves at the same age seems to merge on a log scale for large enough $`t_w`$. This results from the fact the perturbed and unperturbed set of curves can be deduced one from the other by a simple shift in time. At long enough $`t_w`$, the differences between $`t_w`$ and $`t_w+t_{eff}`$ can no longer be seen.

### III.3 Low shear behavior

We now focus on the behavior following a small amplitude shear (2.9%). From figure 4 it appears that the systems keeps on aging under the mechanical oscillation. One could think that the spontaneous aging is nonetheless hindered by the shear. Actually, this is the opposite that happens: the sollicitation helps the system to age as demonstrated on figure 6(a). We compare the shape of the curves after the shear $`|g_1(t_w,t,\gamma =2.9\%)|^2`$ to the curves obtained with no second burst: $`|g_1(t_w,t,\gamma =0\%)|^2`$. It appears that $`|g_1(t_w,t,\gamma =2.9\%)|^2`$ decorrelates more slowly than $`|g_1(t_w,t,\gamma =0\%)|^2`$. Hence the set of sheared curves seems older than the unperturbed one.

It clearly shows that the effect of low shear is to age the system faster. We called this effect overaging. Since the global shape of the two sets of curves seems similar it is interesting to see if one can be deduced from the other by a simple shift in time. Thus we computed the value of the effective delay. A positive value of $`t_{eff}`$ could be found that led to a good overlay of the perturbed and unperturbed curves (see figure 6(b)) for all studied $`t_w`$. However the value of $`t_{eff}`$ now varies with $`t_w`$ as it is shown on figure 7. It thus clearly demonstrates that even though each perturbed curve can be deduced from an unperturbed curve by a simple time translation, the full set of sheared curves cannot be inferred that way from the unperturbed one.

It thus differs from the high shear case. However, $`t_{eff}`$ seems to plateau for large waiting times. For long waiting time the ratio $`\frac{t_{eff}}{\tau _\alpha }`$ tends to zero. Indeed $`\tau _\alpha `$ increases as $`t_w`$ whereas $`t_{eff}`$ tends to plateau. Hence $`t_{eff}`$ becomes increasingly difficult to determine as $`t_w`$ is getting larger. It also implies that on a log scale the perturbed correlation functions and the unperturbed ones tend to merge. This behavior is thus the same as the behavior for high shear.

### III.4 Intermediate shear behavior

If now the amplitude of the second burst is set at a higher value, the system behavior differs from the one we have just described. We set the amplitude $`a`$ of the secondary oscillation to the value of $`5.9\%`$ without changing the other parameters. As previously explained the set of perturbed curved recorded after the shear cessation is compared to that of the unperturbed case. The results are displayed on figure 8(a). The main difference from the preceding case lies in the fact that the shape of $`|g_1(t_w,t,\gamma =5.9\%)|^2`$ is not the same as the one of $`|g_1(t_w,t,\gamma =0\%)|^2`$.

The change in shape can be described as follows:
\- For short $`t_w`$ (1s), the function $`|g_1(t_w,t,\gamma =5.9\%)|^2`$, decreases faster than $`|g_1(t_w,t,\gamma =5.9\%)|^2`$ at short time t but slower at long time t. The dynamics thus seems partially rejuvenated at short t and overaged at long t. It follows that the two curves for a fixed $`t_w`$ eventually cross one another.
\- For longer $`t_w`$ (1212s), the rapid decrease at short time disapears and the function $`|g_1(t_w,t,\gamma =5.9\%)|^2`$ decorrelates more slowly than the unperturbed one. One could notice that the shape $`|g_1(t_w,t,\gamma =5.9\%)|^2`$ in the limit of large $`t_w`$ is close to that of the unperturbed case. Indeed in that limit one can again find a positive effective delay $`t_{eff}`$ that allows the overlay of an unperturbed curve and a sheared one (fig 8(b)). This delay also depends on $`t_w`$. As in the former case the ratio $`\frac{t_{eff}}{t_w}`$ seems to tend towards 0 and the perturbed and reference curves merge on a log scale in the limit of large $`t_w`$.

## IV Discussion

In some previous articlesViasnoff and Lequeux (2002); Viasnoff et al. (2003), we explained how the observed phenomenology can be described by the Soft Glassy Rheology model Sollich et al. (1997); Sollich (1998). In this section we try to extract some physical ingredients that characterize the change of the microscopic dynamics of our system under shear. We analyze our data with no reference to any specific theoretical model in order to make our conclusions as general as possible. We first show that our results are fully compatible with a weak long term memory (WLTM)Bouchaud et al. (1998), and that they can be consistently understood considering the full distribution of relaxation times in the system. We then comment upon the ability to rejuvenate totally the sample with a shear of finite amplitude applied for a finite duration. We also emphasize the differences between our observed rejuvenation and the phenomenon called rejuvenation for spin glasses. Eventually we propose a spatial scenario that could explain the phenomenology.

### IV.1 Weak long term memory

In the long waiting time limit the three regimes display a common behavior: the perturbed and unperturbed curves seems to merge, at least on a log scale. In this section we will show that it can be interpreted in terms of Weak Long term Memory Bouchaud et al. (1998).
In a simple aging experiment, the correlation fonctions $`|g_1(t_w,t,\gamma =0\%)|^2`$ rescale onto a master curve. It thus shows that their evolution during such a type of experiment is autosimilar in time: $`|g_1(t_w,t,\gamma =0\%)|^2`$ can be written as a simple function $`f(\frac{t}{t_w})`$ of the variable $`\frac{t}{t_w}`$. However, it is not straightforward that any shear history would lead to correlation functions that can be written that way. Nonetheless, the experimental section shows that in the long time limit the perturbed correlation functions $`|g_1(t_w,t,\gamma )|^2`$ approximately write as $`|g_1(t_w+t_{eff},t,\gamma =0\%)|^2`$. In the case of a large amplitude, one can define a constant negative effective delay, in the case of a small sollicitation $`t_{eff}`$ varies with $`t_w`$, however $`t_{eff}`$ seems to plateau at large $`t_w`$. For intermediate $`\gamma `$ the situation is slightly different; no $`t_{eff}`$ can be found for most values of $`t_w`$. The shapes of the correlation fonctions clearly differ from the shape of the simple aging one. However, in the limit of large $`t_w`$ the shape of the perturbed curves tend towards the reference one. Hence for any strain history and in limit of large $`t_w`$, the perturbed correlation functions can be written as:

$$|g_1(t_w,t,\gamma )|^2|g_1(t_w[1+\frac{t_{eff}}{t_w}],t,\gamma =0\%)|^2$$

which also reads, using the master function $`f`$:

$$|g_1(t_w,t,\gamma )|^2f\left(\frac{t}{t_w+t_{eff}}\right)f\left(\frac{t}{t_w+o(t_w)}\right)$$

It physically means that, whatever its initial dynamical state, the system will always spontaneously evolves towards a unique asymptotic dynamical state. The effect of any finite perturbation applied for a finite amount of time will thus eventually fade away.
This observation is reminiscent of the concept of weak long term memory (WLTM)Bouchaud et al. (1998). WLTM means that the system is able to explore its whole phase space whatever its initial state and even though some regions may take a very long time to be reached. As a consequence, the disturbance of the dynamics by a finite perturbation will disappear in the limit of large waiting times. Even though this principle has been rigourously observed in many theoretical models, it is hardly discussed in experiments. Our results are good hints that such a principle applies to soft glassy materials, at least from the point of view of their microscopic dynamics.

### IV.2 Modifications of the relaxation time distribution

In the preceding sections, we discussed the shape of the correlation fonctions. An equivalent way of describing them is to consider their distributions $`P(\tau _{dec})`$ of decorrelation times. Theses distributions are linked to the distributions of relaxation times $`P(\tau _{rel})`$ in the sample. This point of view is able to give us a coherent image of the observed phenomenology.

The aging process viewed from the point of view of $`P(\tau _{rel})`$corresponds to a slow drift of the distribution towards slower and slower relaxation times. Furthermore, the rescaling of all the correlation functions implies that this slow drift happens in a self similar manner. Hence, the distribution must rescale with the waiting time. It thus must write:

$$P(\tau _{rel},t_w)=P_m(\frac{\tau _{rel}}{t_w})$$

where $`P_m`$ is the master function for P. Similarly, the weak long term memory behavior can be understood as the convergence of any initial distribution towards a unique asymptotique master function $`P_m`$ in the limit of large waiting times. At that point, the observed phenomenology can be understood as follows:

* When the perturbation drives the system out of its glassy phase, the distribution is shifted towards faster relation times. It is, within the experimental precisions, similar to a distribution through which the system passes during its reference evolution. In fact a simple time translation allows a good overlay of the whole set of correlation curves.
* Oppositely, a small amplitude perturbation leads to an overaging of the system. It means that the microscopic dynamics is slowed down after the perturbation. As a result the distribution is shifted towards longer relaxation time. This result is not intuitive, and the nature of the shift is not trivial. One could think that it corresponds to a speeding up of the autosimilar drift of the distribution. If it were so, the full set of the perturbed correlation functions would corresponds to a time translation of the reference one via a unique effective delay $`t_{eff}`$. However, we demonstrated, in section III.3, that $`t_{eff}`$ actually depends on $`t_w`$. It implies that the perturbation enriched the distribution P with slow relaxation time in a way that differs from that of a spontaneous evolution. This point will be further discussed in section IV.3.
* For intermediate amplitudes, the shape of the distribution P is greatly modified compared to the reference case. This modification entails the drastic change in shape of the correlation function. The distribution is enriched simultaneously with short and long relaxation times. This dynamical state is reflected by the faster decrease of $`|g_1(t_w,t,\gamma )|^2`$ for small $`t`$ and its slower decrease for large $`t`$. This holds for time $`t_w`$ shortly following the shear cessation. However, for larger $`t_w`$, we have already mentioned that the system shows a weak long term memory. It means that the dynamical state that we have just described eventually converges towards the distribution $`P_m`$. Our results further show that this convergence occurs in the following fashion: the short time part of P reaches quickly retrieved its unperturbed shape. It is reflected by the rapid disappearance of the fast decrease of the perturbed correlation function for small $`t`$. The large relaxation times, however, remain overpopulated compared to the reference history. The correlation function thus show a slowest decay at large $`t`$. The rejuvenated part of the distribution has disappeared, and the overaged part is left alone. A positive effective delay can thus be found. For theses values of $`t_w`$, the intermediate and small amplitude cases are equivalent.

Theses conclusions and interpretations are illustrated on figure 9(left) and do not depend on any underlying model. However, in a previous paper we showed that such a phenomenological scenario can be achieved within the frame of the Soft Glassy Rheology ModelViasnoff et al. (2003).

### IV.3 relevance of the effective delay description

The literature of glasses often describes the effect of a small positive or negative shift in temperature in term of effective age or effective delay. The two usual ways to evaluate the effective delay for these systems consist in either measuring one time quantities such as the susceptibility or the dielectric constant at a fixed frequency $`(\chi \mathrm{"}(\omega ),\epsilon \mathrm{"}(\omega ))`$ or two times quantities such as relaxation functions $`(M(t_w,t),G(t_w,t))`$. The effective delay is determined by shifting the curves after the perturbation onto a reference curve. Apart from some deviations shortly after the perturbation cessation, a good overlay can be found. In this section we emphasize the fact that some misinterpretations can arise if the full evolution of a two times quantity is not considered to evaluate the effective delay.

If the description of the change in the dynamics in terms of effective delay would be an exact result, it would then mean that the perturbation drives the system into a dynamical state through which it would have passed during a simple aging history. Hence an effective age would be fully sufficient to characterize this state. In this event, the relaxation time distribution after the perturbation at time $`t_w`$ would be truly identical to the reference distribution at time $`t_w+t_{eff}`$. This is schematically presented on fig 9 a). As a consequence the evolution of the perturbed distribution would be the exact evolution of the reference one but shifted by the same delay $`t_{eff}`$.

In this event the determination of an effective delay by the overlay of one single correlation curve with a reference one, or the overlay of the evolution of a single typical point of each correlation function with the evolution of that same point for a reference history, would lead to the same result. This is what happens in the high shear regime (see fig 5) where $`t_{eff}`$ is found constant with $`t_w`$. However, in the low shear regime, $`t_{eff}`$ depends on $`t_w`$. This implies that looking at the evolution of a single point value is thus no longer equivalent to looking at the system at a fixed $`t_w`$. The shape of the distribution is changed. However, if the change is slight for relaxation times corresponding to the time window over which the measurement is performed, an apparent $`t_{eff}`$ can be found. However, if the evolution of the distribution depends on the entire spectrum of relaxation times, the evolution of the probed part of the distribution may differ from the reference one. This situation is illustrated on figure 9 b). Although the distribution seems to overlay over a certain time range a t a fixed time, their evolution differs, since the evolution takes into account the whole distribution. Hence the value $`t_{eff}`$ has to be varied with $`t_w`$ in order to keep a good overlay of the reference and perturbed distributions.
We remark that this variation in time of $`t_{eff}`$ cannot be detected by the usual methods described in the first lines of this section. This effect can be measured here since our experiment gives easy access to the evolution of a two times quantity:$`|g_1(t_w,t)|^2`$. This is the reason why we think this notion of effective delay should be considered carefully.

### IV.4 High shear rejuvenation

Macroscopic observations clearly shows that the rheological properties of such systems can reach a steady state under high shearDerec et al. (2000); Cloitre et al. (2000); Bonn et al. (2002); van Megen and Underwood (1994). The shear cessation acts then as a thermal quench from a high temperature state. However, the type of shear that are usually used are either continuous ones at a fixed shear rate, or oscillatory ones at very high amplitude (typically 1000%). We now show that, in the case of a sample with a volume fraction close to that of its dynamical transition, the microscopic dynamics can be entirely renewed applying an oscillatory shear of amplitude around 20%. This shear sollicitation is able accelerate the rate of the structural relaxations down to the shortest time that our set up can measure $`0.05s`$. A careful look at the correlation functions displayed on figure 2 shows a fast decay of all the correlation functions in their very first points. This decay corresponds to the tail of the fast $`\beta `$ decay existing in all glasses. They correspond to individual fast relaxation modes that are reminiscent of the relaxation modes in the liquid phase. The long time decay that we studied here is often referred as the $`\alpha `$ modes that arise from cooperative structural rearrangements. Theses modes are the signature of the vicinity of the glass transition. Notice that on the correlation function taken 0.1s after the shear cessation, the $`\alpha `$ and $`\beta `$ modes are approximately separated by only 1 decade. It thus means that under high shear theses modes are hardly separated in time. The relaxation processes are thus very close to that of a fluid. This is the logical reason why higher or longer shear periods do not modify the dynamics following the shear cessation.

We also want to point out a real difference between the rejuvenation phenomenon we observe here and the rejuvenation described for spin glasses or polymers. As demonstrated by figure 4, rejuvenation occurs if the shear amplitude $`a`$ is superior to 8.9%. The condition for our sample to be rejuvenated is that it is rheologically quenched from a steady state into a glassy state. In spin glasses the observed rejuvenation (restarting of the dynamics after a negative temperature step) occurs when the system is thermally quenched from a glassy state into another glassy state. In those systems, the dynamics at lower temperature does no depend on the time spent at a higher temperature. These effects can be described among other descriptions by a hierarchical organisation of the phase space Hammann et al. (1992). In our case, no such descriptions are needed to explain our results since rejuvenation occurs only if the colloidal glass is driven out of its glassy region.
In order to prove that there is no rejuvenation effects in the sens of spin glasses, we performed the following experiment. We compare the correlation function of a system under shear to the correlation fonction of the same system left at rest. The shear amplitude is fixed to 5.9% such that it lies in the region where the system keeps on aging despite the shear application see fig 4. The shear history is displayed on figure 10 (a).

Figure 10 (b) shows that the height of the echoes decays faster than the reference curve. Hence, the forced dynamics results in an additional degree of irreversibility in the scatterers motion. It clearly appears that the correlation function decreases faster under shear than when left at rest. Though this experiments is the exact analogous to the negative temperature shift experiment, its results are clearly at odd with the rejuvenation effect for spin glasses Vincent et al. (1996). We thus conclude that the similarity between temperature shift and oscillating shear only lies in the fact that they transiently modify the rearrangement rate.

We finally make the following remark: In any case (rejuvenation, overaging or crossover regime) the microscopic dynamics is accelerated while the system is under shear. This is fully consistent with the shear thinning behavior of such systems. The different behaviors described in the previous sections result from the spontaneous evolution of the structured formed when the shear is stopped. In that prospect, overaging is a result of the shear thinning behavior
On the other hand shear thickening material may present a different phenomenology: the dynamics under shear may appear slower than the spontaneous one. Thus the shear cessation can be followed by an acceleration of the dynamics. Even though such results may look phenomenologically identical to the spin glass rejuvenation effect, they may arise for totally different reasons such as anisotropy and not because of hierarchical organisation of relaxation times.

### IV.5 Real space proposed scenario

In the preceding sections we described our results in terms of distribution of relaxation times. However, we do not know the very nature of the relaxation modes associated with theses times. In this section, we propose a picture in real space of a possible scenario that can occur in such materials. We are aware that this scenario is certainly to simplistic, and can be view as a basis for more elaborated ones.
On macroscopic length scales, the sample is affinely deformed by the macroscopic shear. However, the heterogeneous structure of the rheological properties of the material on smaller length scale can result in stress and/or strain concentration on certain regions. If the stress is locally superior to the local yield stress, a local rearrangement is induced by the macroscopic deformation.
$``$ At small shear amplitude, one could think that these rearrangements events occur in the most fragile regions that are randomly distributed in the sample. They may involve few particles. The external shear thus acts as a ”boost” for the rearrangement rate. and increased probability is provided to the system to find more stable configuration. The effect of mechanical deformation may thus be to add to the thermal noise a random mechanical noise due to the increased number of rearrangements.
$``$ At high amplitudes the shear strain is sufficient to disrupt all structures in the sample. Shear is high enough to overcome all local yield stress. The system is fluidized, an equilibrium dynamics is found that resemble that of a liquid (see fig 4(b)). This corresponds to the usual scenario of shear thinning and total rejuvenation in such systems.
$``$ At intermediate amplitude one can naively imagine that two types of zones coexist: the shear amplitude is such that shear induced rearrangements occur in localized mesoscopic zones of the sample. Though on macroscopic length scale the system seems to be homogeneously shear (no macroscopic shear banding), the strain is localized in fluidized areas that thus have a rapid dynamic under and after shear cessation. Theses rapidly rearranging zones may delimit stiff and thus slow dynamics areas that undergo a small straining and mechanical noise from their boundaries. The small perturbation applied onto theses zones may thus lead to an overaging of their dynamics as in the small shear amplitude scenario. This is would explain the simultaneous appearance of the two trends (rejuvenation and overaging) within the same sample. The volume fraction of the overaged zones would highly depend on the shear amplitude, ranging from 0 for a high amplitude to a macroscopic size for very small amplitudes.
We also want to point out that such a scenario occurs only if no macroscopic shear banding occur. In that case a total rejuvenation of the sample by a stationary shear is not possible. This is why theses effect are only observed in a narrow range of volume fraction as discussed in Viasnoff et al. (2003).

## V conclusion

In this article we discussed in detail the effect of shear upon the microscopic dynamics of a colloidal glass. We showed theses effect can be understood as a change in the relaxation time distribution. Small amplitude shear result in drift of the distribution towards longer relaxation time. This drift happen in a different fashion than the spontaneous aging one. Intermediate amplitude broadens the distribution both in the short and long relaxation time direction. High amplitude shear allow a steady state to be reached that entirely rejuvenates the system. We further show that no matter the initial distribution shape, it always converges towards a unique distribution for long waiting times. We also made a point out of the ambiguity of the determination of an effective age for glassy systems submitted to small perturbations.

The authors thank J.P Bouchaud and L Berthier for enlightening discussions.


format: text

text:
# TUM-HEP-513/03 Can four neutrinos explain global oscillation data including LSND & cosmology?11footnote 1Talk given by T.S. at the NOON 2003 workshop, February 10-14, 2003, Kanazawa, Japan

## 1 Introduction

The neutrino oscillation interpretations of the solarsolar , sno  and KamLANDkamland  neutrino experiments, atmosphericskatm , macro  neutrino data, and the LSND experimentlsnd  require three neutrino mass-squared differences of different orders of magnitudePakvasa:2003zv . Since it is not possible to obtain this within the Standard Model framework of three active neutrinos it has been proposed to introduce a light sterile neutrinosterile  to reconcile all the experimental hints for neutrino oscillations. Here we present an analysis of the global neutrino oscillation data in terms of four-neutrino mass schemes, including data from solar, KamLAND and atmospheric neutrino experiments, the LSND experiment, as well as data from short-baseline (SBL) experimentsKARMEN , CDHS , bugey  and reactor experimentslbl  reporting no evidence for oscillations. We find that for all possible types of four-neutrino schemes different sub-sets of the data are in serious disagreement and hence, four-neutrino oscillations do not provide a satisfactory description of the global oscillation data including LSND. The details of our calculations can be found in Refs. 4nu01 , 4nu02 , solat02 .

## 2 Notations and approximations

Four-neutrino mass schemes are usually divided into the two classes (3+1) and (2+2), as illustrated in Fig. 1. We note that (3+1) mass spectra include the three-active neutrino scenario as limiting case. In this case solar and atmospheric neutrino oscillations are explained by active neutrino oscillations, with mass-squared differences $`\mathrm{\Delta }m^2sol`$ and $`\mathrm{\Delta }m^2atm`$, and the fourth neutrino state gets completely decoupled. We will refer to such limiting scenario as (3+0). In contrast, the (2+2) spectrum is intrinsically different, as there must be a significant contribution of the sterile neutrino either in solar or in atmospheric neutrino oscillations or in both.

Neglecting CP violation, in general neutrino oscillations in four-neutrino schemes are described by 9 parameters: 3 mass-squared differences and 6 mixing angles in the unitary lepton mixing matrix. Here we use a parameterisation introduced in Ref. 4nu01 , which is based on physically relevant quantities: the 6 parameters $`\mathrm{\Delta }m^2sol`$, $`\theta sol`$, $`\mathrm{\Delta }m^2atm`$, $`\theta atm`$, $`\mathrm{\Delta }m^2lsnd`$, $`\theta lsnd`$ are similar to the two-neutrino mass-squared differences and mixing angles and are directly related to the oscillations in solar, atmospheric and the LSND experiments. For the remaining 3 parameters we use $`\eta _s,\eta _e`$ and $`d_\mu `$. Here, $`\eta _s(\eta _e)`$ is the fraction of $`\nu _s(\nu _e)`$ participating in solar oscillations, and ($`1d_\mu `$) is the fraction of $`\nu _\mu `$ participating in oscillations with $`\mathrm{\Delta }m^2atm`$ (for exact definitions see Ref. 4nu01 ). For the analysis we adopt the following approximations: (1) We make use of the hierarchy $`\mathrm{\Delta }m^2sol\mathrm{\Delta }m^2atm\mathrm{\Delta }m^2lsnd`$. This means that for each data set we consider only one mass-squared difference, the other two are set either to zero or to infinity. (2) In the analyses of solar and atmospheric data (but not for SBL data) we set $`\eta _e=1`$, which is justified because of strong constraints from reactor experimentsbugey , lbl .

Due to these approximations the parameter structure of the four-neutrino analysis gets rather simple. The parameter dependence of the four data sets solar, atmospheric, LSND and NEV is illustrated in Fig. 2. The NEV data set contains the experiments KARMENKARMEN , CDHSCDHS , Bugeybugey  and CHOOZ/Palo Verdelbl , reporting no evidence for oscillations. We see that only $`\eta _s`$ links solar and atmospheric data and $`d_\mu `$ links atmospheric and NEV data. LSND and NEV data are coupled by $`\mathrm{\Delta }m^2lsnd`$ and $`\theta lsnd`$.

## 3 (2+2): ruled out by solar and atmospheric data

The strong preference of oscillations into active neutrinos in solar and atmospheric oscillations leads to a direct conflict in (2+2) oscillation schemes. We will now show that thanks to the latest solar neutrino data (especially from SNOsno ) in combination with the KamLAND experimentkamland , and the improved SK statistic on atmospheric neutrinosskatm  the tension in the data has become so strong that (2+2) oscillation schemes are essentially ruled out.222Details of our analyses of the solar, KamLAND and atmospheric neutrino data can be found in Refs. solat02 , KL . For an earlier four-neutrino analysis of solar and atmospheric data see Ref. concha4nu .

The thin lines in the left panel of Fig. 3 show the $`\mathrm{\Delta }\chi ^2`$ of latest solar data as a function of $`\eta _s`$, the parameter describing the fraction of the sterile neutrino participating in solar neutrino oscillations. The 99% CL bounds $`\eta _s0.44`$, when the SSM constraint on the 8B-flux is included, and $`\eta _s0.61`$ for free 8B-flux reflect the strong preference for active neutrino oscillations of solar data. Recently, the outstanding results of the KamLAND reactor experimentkamland  confirmed the LMA solution of the solar neutrino problemKL , sol+KL . Apart from this very important result the sensitivity of KamLAND to an admixture of a sterile neutrino is rather limited (see think lines in Fig. 3). The combined analysis leads to the 99% CL bound $`\eta _s0.5`$ for 8B free, whereas in the SSM constrained case the bound is unchanged.

In contrast, in (2+2) schemes atmospheric data imply $`\eta _s0.65`$ at 99% CL, in clear disagreement with the bound from solar data. In the right panel of Fig. 3 we show the $`\chi ^2`$ for solar data and for atmospheric combined with SBL data as a function of $`\eta _s`$. Furthermore, we show the $`\chi ^2`$ of the global data defined by

$$\overline{\chi }^2(\eta _s)\mathrm{\Delta }\chi ^2sol(\eta _s)+\mathrm{\Delta }\chi _{\mathrm{atm}+\mathrm{sbl}}^2(\eta _s).$$
(1)
From the figure we find that only if we take both data sets at the 99.95% CL a value of $`\eta _s`$ exists, which is contained in the allowed regions of both sets. This follows from the $`\chi ^2`$-value $`\chi _{\mathrm{PC}}^2=12.2`$ shown in the figure. In Refs. 4nu02 , PG  we have proposed a statistical method to evaluate the disagreement of different data sets in global analyses. The parameter goodness of fit (PG) makes use of the $`\overline{\chi }^2`$ defined in Eq. (1). This criterion evaluates the GOF of the combination of data sets, without being diluted by a large number of data points, as it happens for the usual GOF criterion (for details see Ref. PG ). We find $`\chi _{\mathrm{PG}}^2\overline{\chi }_{\mathrm{min}}^2=23.5`$, leading to the marginal PG of $`1.3\times 10^6`$. We conclude that (2+2) oscillation schemes are highly disfavoured by the disagreement between the latest solar and atmospheric neutrino data. This is a very robust result, independent of whether LSND is confirmed or disproved.333Sub-leading effects beyond the approximations adopted here should not affect this result significantly. Allowing for additional parameters to vary might change the ratio of some observablesPas:2002ah , however, we expect that the absolute number of events relevant for the fit will not change substantially.

## 4 (3+1): strongly disfavoured by SBL data

It is known for a long time3+1early , okada  that (3+1) mass schemes are disfavoured by the comparison of SBL disappearance dataCDHS , bugey  with the LSND result. In Ref. cornering  we have calculated an upper bound on the LSND oscillation amplitude $`\mathrm{sin}^22\theta lsnd`$ resulting from NEV and atmospheric neutrino data. From Fig. 4 we see that this bound is incompatible with the signal observed in LSND at the 95% CL. Only marginal overlap regions exist between the bound and global LSND data if both are taken at 99% CL. An analysis in terms of the parameter goodness of fit4nu02  shows that for most values of $`\mathrm{\Delta }m^2lsnd`$ NEV and atmospheric data are compatible with LSND only at more than $`3\sigma `$, with one exception around $`\mathrm{\Delta }m^2lsnd6`$ eV2, where the PG reaches 1%. These results show that (3+1) schemes are strongly disfavoured by SBL disappearance data.

## 5 Comparing (3+1), (2+2) and (3+0) hypotheses

With the methods developed in Ref. 4nu01  we are able to perform a global fit to the oscillation data in the four-neutrino framework. This approach allows to statistically compare the different hypotheses. Let us first evaluate the GOF of (3+1) and (2+2) spectra with the help of the PG method described in Ref. PG . We divide the global oscillation data in the 4 data sets SOL, ATM, LSND and NEV. Following Ref. 4nu02  we consider

$$\begin{array}{ccc}\overline{\chi }^2& =& \mathrm{\Delta }\chi ^2sol(\theta sol,\mathrm{\Delta }m^2sol,\eta _s)+\mathrm{\Delta }\chi ^2atm(\theta atm,\mathrm{\Delta }m^2atm,\eta _s,d_\mu )\hfill \\ & +& \mathrm{\Delta }\chi ^2nev(\theta lsnd,\mathrm{\Delta }m^2lsnd,d_\mu ,\eta _e)+\mathrm{\Delta }\chi ^2lsnd(\theta lsnd,\mathrm{\Delta }m^2lsnd),\hfill \end{array}$$
(2)
where $`\mathrm{\Delta }\chi _X^2=\chi _X^2(\chi _X^2)_{\mathrm{min}}`$ ($`X`$ = SOL, ATM, NEV, LSND). In Tab. 1 we show the contributions of the 4 data sets to $`\chi _{\mathrm{PG}}^2\overline{\chi }_{\mathrm{min}}^2`$ for (3+1) and (2+2) oscillation schemes. As expected we observe that in (3+1) schemes the main contribution comes from SBL data due to the tension between LSND and NEV data in these schemes. For (2+2) oscillation schemes a large part of $`\chi _{\mathrm{PG}}^2`$ comes from solar and atmospheric data, however, also SBL data contributes significantly. This comes mainly from the tension between LSND and KARMENChurch:2002tc , which does not depend on the mass scheme and, hence, also contributes in the case of (2+2). Therefore, the values of $`\chi _{\mathrm{PG}}^2`$ in Tab. 1 for (2+2) schemes are higher than the one given in Sec. 3, where the tension in SBL data is not included.

The parameter goodness of fit is now obtained by evaluating $`\chi _{\mathrm{PG}}^2`$ for 4 DOFPG . This number of degrees of freedom corresponds to the 4 parameters $`\eta _s,d_\mu ,\theta lsnd,\mathrm{\Delta }m^2lsnd`$ describing the coupling of the different data sets (see Eq. (2) and Fig. 2). The best GOF is obtained in the (3+1) case. However, even in this best case the PG is only 0.56%. The PG of $`1.6\times 10^6`$ for (2+2) schemes shows that these mass schemes are essentially ruled out by the disagreement between the individual data sets.

Although we have seen that none of the four-neutrino mass schemes can provide a reasonable good fit to the global oscillation data including LSND, it might be interesting to consider the relative status of the three hypotheses (3+1), (2+2) and the three-active neutrino scenario (3+0). This can be done by comparing the $`\chi ^2`$ values of the best fit point (which is in the (3+1) scheme) with the one corresponding to (2+2) and (3+0). First we observe that (2+2) schemes are strongly disfavoured with respect to (3+1) with a $`\mathrm{\Delta }\chi ^2=17.8`$, implying a CL of 99.87% for 4 DOF. Further we find that (2+2) is only slightly better than (3+0), which is disfavoured with a $`\mathrm{\Delta }\chi ^2=20.0`$ with respect to (3+1).

## 6 Adding information from cosmology

Recently a precision measurements of the cosmic microwave background (CMB) radiation has been published by the WMAP collaborationwmap . Under certain assumptions about the cosmological model and in combination with other cosmological data a stringent upper bound on the sum of the neutrino masses $`m_\nu \stackrel{<}{}\mathrm{\hspace{0.25em}0.7}`$ eV is obtainedwmap , Elgaroy:2003yh . Since this bound is of the same order as $`\sqrt{\mathrm{\Delta }m^2lsnd}`$ it will lead to further restrictions in four-neutrino mass schemeswmap3+1 . To combine this information with data from oscillation experiments we use the results of an analysis performed by S. Hannestadhannestad . There, WMAPwmap  results have been combined with other CMB dataCMB , data from the large scale structure survey 2dFGRSLSS , the HST Hubble key projectHST , and Type Ia supernovae observationsSN . It has been found that a non-zero $`m_\nu `$ can be compensated partially by an increase in the effective number of neutrino species. Therefore, the original bound obtained in Ref. wmap , which was calculated assuming three neutrinos, cannot be applied in a straightforward way. In the case of four neutrinos a somewhat weaker bound of $`m_\nu <1.38`$ eV (95% CL) is obtainedhannestad . To combine this result with data from oscillation experiments we use the $`\chi ^2`$ from cosmology for four neutrinos, shown in Fig. 6 of Ref. hannestad .

In the following analysis we consider (3+1) schemes of the hierarchical type (first two spectra in Fig. 1) and neglect the masses of the three lower mass states, such that $`m_\nu \sqrt{\mathrm{\Delta }m^2lsnd}`$.444Note that all other four neutrino schemes will be much more disfavoured, since more mass states will contribute at the LSND mass scale. The $`\mathrm{\Delta }\chi _{\mathrm{cosmology}}^2`$ from Ref. hannestad  under this assumption is shown in the right panel of Fig. 5 as a function of $`\mathrm{\Delta }m^2lsnd`$. From the left panel of this figure one finds that cosmology excludes the region $`\mathrm{\Delta }m^2lsnd\stackrel{>}{}\mathrm{\hspace{0.25em}3}`$ eV2. However, the islands at $`\mathrm{\Delta }m^2lsnd0.91`$ eV2 and 1.7 eV2 are practically unaffected by cosmological data. Performing a GOF analysis shows that adding cosmological data gives $`\chi _{\mathrm{PG}}^2=17.1`$ and 19.6 for the two islands. This should be evaluated for 5 DOFPG , leading to a PG of 0.43% and 0.15%, respectively, which is comparable to the value given for (3+1) in Tab. 1. We note, however, that a modest improvement of the bound on $`m_\nu `$ by future cosmological data will drastically worsen the fit of the last remaining islands in the four-neutrino parameter space.

A further cosmological constraint on four-neutrino schemes comes from Big Bang nucleo synthesis (BBN)okada , BGGS . In the absence of any non-standard mechanism, like a large lepton number asymmetry, (3+1) as well as (2+2) schemes will lead to an effective number of neutrino species $`N_\nu =4`$BBNeq . This is in conflict with the upper bound on $`N_\nu `$ derived in recent analyseshannestad , Crotty:2003th , Barger:2003zg  including CMB and BBN data. Taking into account the uncertainty related to the primordial abundances of He and D, upper bounds in the range $`N_\nu <3.13.3`$ ($`2\sigma `$)Barger:2003zg  are obtained.

## 7 Conclusions

Performing a global analysis of current neutrino oscillation data we find that four-neutrino schemes do not provide a satisfactory fit to the data:

$``$ The strong rejection of non-active oscillation in the solar+KamLAND and atmospheric neutrino data rules out (2+2) schemes, independent of whether LSND is confirmed or not. Using an improved goodness of fit method especially sensitive to the combination of data sets we obtain a GOF of only $`1.6\times 10^6`$ for (2+2) schemes.

$``$ (3+1) spectra are disfavoured by the disagreement of LSND with short-baseline disappearance data, leading to a marginal GOF of $`5.6\times 10^3`$. If LSND should be confirmed we need more data on $`\nu _e`$ and/or $`\nu _\mu `$ SBL disappearance to decide about the status of (3+1).

$``$ Recent cosmological data lead to further restrictions of the allowed parameter space. Only two islands aground $`\mathrm{\Delta }m^2lsnd0.91`$ eV2 and 1.7 eV2 in hierarchical (3+1) schemes survive.

Our analysis brings the LSND hint to a more puzzling status, and the situation will become even more puzzling if LSND should be confirmed by the up-coming MiniBooNE experimentminiboone . Recently an explanation in terms of five neutrinos has been proposedSorel:2003hf . In such a (3+2) scheme solar and atmospheric data are explained dominantly by active oscillations, similar to the (3+1) case. The fit of SBL data is improved with respect to (3+1) by involving two large mass splittings $`\mathrm{\Delta }m_{41}^20.9`$ eV2 and $`\mathrm{\Delta }m_{51}^220`$ eV2. Alternative explanations of LSND by anomalous muon decayBabu:2002ic  seem to be in disagreement with the negative result of KARMEN, or require a very radical modification of standard physics, like CPT violationBarenboim:2002ah .

Acknowledgements: T.S. would like to thank the organisers for the very interesting workshop and for financial support. Furthermore, we thank S. Hannestad for providing us with the $`\chi ^2`$ of cosmological data, and S. Pastor for discussions. This work was supported by Spanish grant BFM2002-00345, by the European Commission RTN grant HPRN-CT-2000-00148 and by the ESF Neutrino Astrophysics Network. M.M. was supported by Marie Curie contract HPMF-CT-2000-01008, T.S. was supported by the “Sonderforschungsbereich 375-95 für Astro-Teilchenphysik” der Deutschen Forschungsgemeinschaft and M.A.T. was supported by the M.E.C.D. fellowship AP2000-1953.


text:
# Thermal Stabilization of the HCP Phase in Titanium

## Abstract

We have used a tight-binding model that is fit to first-principles electronic-structure calculations for titanium to calculate quasi-harmonic phonons and the Gibbs free energy of the hexagonal close-packed (hcp) and omega ($`\omega `$) crystal structures. We show that the true zero-temperature ground-state is the $`\omega `$ structure, although this has never been observed experimentally at normal pressure, and that it is the entropy from the thermal population of phonon states which stabilizes the hcp structure at room temperature. We present the first completely theoretical prediction of the temperature- and pressure-dependence of the hcp-$`\omega `$ phase transformation and show that it is in good agreement with experiment. The quasi-harmonic approximation fails to adequately treat the bcc phase because the zero-temperature phonons of this structure are not all stable.

preprint: LA-UR-03-0239

The experimentally observed ground state of titanium at room temperature is the hexagonal close-packed (hcp) structureSikka et al. (1982), even though zero-temperature first-principles electronic-structure calculationsGyanchandani et al. (1990); Ahuja et al. (1993); Jomard et al. (1998); Greeff et al. (2001); Mehl and Papaconstantopoulos (2002) predict the omega ($`\omega `$) structure to be the energetically favored ground state. There are two possible causes for this discrepancy: (1) this is another failure of the local-density approximation (LDA) for calculating electronic-structure of solids, or (2) thermal effects stabilize the hcp phase. In this paper we show that the second argument is correct: the thermal occupation of phonon states favors the hcp structure above a pressure-dependent transition temperature. Only by including the vibrational entropy in the first-principles calculations can theory be brought into agreement with experiment.

Titanium is a good case for testing our ability to predict solid-solid phase transformations with first-principles methods. It is an element that displays a rich phase diagram with several recently discovered high-pressure phasesVohra and Spencer (2001); Akahama et al. (2001). When alloyed with other elements, it leads to materials of major technical importance for the airline, space, and other industries.

It is possible to use first-principles methods to calculate the zero-temperature internal energy, $`\mathrm{\Phi }_0(V)`$, as well as the free energy contributions from the ions, $`F_I(T,V)`$, and the electrons, $`F_E(T,V)`$, to give the complete equation of state

$$F(T,V)=\mathrm{\Phi }_0(V)+F_I(T,V)+F_E(T,V).$$
(1)
Because of the large computational effort required to calculate the free energy of Eq. 1 over a significant range of volumes and temperatures, we have chosen to use a tight-binding model instead of conventional first-principles electronic-structure methods. We will demonstrate the high accuracy of our model and believe that it should be essentially identical with the best available first-principles methods.

The model’s accurate emulation of first-principles calculations (in this case the Wien97 codeBlaha et al. (1999)) requires optimization of a single tight-binding parameterization such that the model closely reproduces the relative total energies of selected crystal structures. The structures are chosen to be both simple and representative: simple enough to make the first-principles calculations feasible, and representative of all microscopic environments the model is aimed at describing. A comparison of total energies as a function of volume with structures fulfilling these requirements are shown in Fig. 1. Not shown but included in the fit as independent structures are fcc crystals distorted to correspond to the longitudinal and transverse phonons at the reciprocal-space high-symmetry point X. These additional structures fine-tune the model’s ability to correctly describe lattice vibrationsRudin et al. (2002). Crucial for convergence of the fitting procedure is to fit the cubic-structure energy bands at high-symmetry points in reciprocal spaceJones and Albers (2002); Papaconstantopoulos (1986).

The electronic contribution to the free energy depends on the electronic density of states (DOS) as a function of volume, $`n(E,V)`$. The occupation of these states, given by the Fermi distribution $`f(E,T)=[e^{(EE_f)/(k_BT)}+1]^1`$, determines their entropyPathria (1972),

$$S_{el}(T,V)=k_B[f\mathrm{ln}f+(1f)\mathrm{ln}(1f)]n(E,V)𝑑E,$$
(2)
and hence the electronic contribution to the free energy $`F_E(T,V)=TS_{el}(T,V)`$. This contribution is small compared to that of the ions but is included for completeness.

The ionic contribution to the free energy depends on the phonon DOS, $`g(\omega ,V)`$, through the zero-point energy,

$$U_{zero}(V)=\frac{1}{2}_\mathrm{\Omega }\mathrm{}\omega g(\omega ,V)𝑑\omega ,$$
(3)
as well as the temperature-dependent, quasi-harmonic free energy,Pathria (1972)

$$F_H(T,V)=k_BT_\mathrm{\Omega }\mathrm{ln}[1e^{\mathrm{}\omega /k_BT}]g(\omega ,V)𝑑\omega .$$
(4)
This differs from the volume-independent harmonic approximation in that the quasi-harmonic approximation uses phonons calculated for small displacements (i.e., harmonic) at each volume. The resulting phonon frequencies are therefore volume dependent but temperature independent (they are calculated at zero temperature). The approximation requires that the phonons for any given crystal structure and volume are stable at zero temperature, which may not be true if the crystal can lower its total energy by spontaneously and continuously distorting into another crystal structure. At elevated temperatures (near melting) this approximation also neglects some important anharmonic contributions to the free energy that can be responsible for stabilizing some high-temperature phases.

We have used the tight-binding model to calculate the electronic and the phonon DOS by evaluating the relevant energy eigenvalues and dynamical matrices on a fine mesh of wave vectors in the first Brillouin zone. In both cases the spectrum is smeared with a Gaussian, chosen to be as small as possible while keeping the DOS smooth and continuous. The mesh is refined until the free energy converges.

The dynamical matrix at a given wave vector q is the Fourier transform of the force constants, which we calculate from the tight-binding model by the direct-force methodKunc and Martin (1982); Wei and Chou (1992); Frank et al. (1995); Parlinski et al. (1997). This method requires large simulation cells consisting of repeated unit cells transposed by vectors $`\mathrm{}`$, the forces on all atoms are calculated in response to the displacement of the atom(s) in one unit cell. Our simulation cells contain 54 (hcp), 81 ($`\omega `$), and 128 atoms (bcc); for more details we refer the reader to Ref.Rudin et al. (2002).

Figure 2 compares our calculated phonon dispersion for the hcp structure to experimental dataStassis et al. (1979). The overall agreement is good, in particular near the zone center (around $`\mathrm{\Gamma }`$) and along the vertical direction (parallel to the c-axis). Frequencies of modes in horizontal directions (perpendicular to the c-axis) away from the zone center are somewhat high compared to experiment. The main characteristics, which the tight-binding model captures, give a reliable DOS, the relevant entity used in Eqs. 3 and 4 to calculate the phonon contribution to the free energy.

The phonon density of states is calculated for the hcp111The $`c`$ over $`a`$ ratio of the hcp structure is varied to give the lowest zero-temperature energy, i.e., in particular at high pressures the calculation is done for a hexagonal crystal with a $`c/a`$ larger than the ideal value., $`\omega `$, and bcc structures at a dozen volumes corresponding to pressures in the range from zero to 170 GPa and is then used to evaluate the Gibbs free energyPathria (1972), $`G(T,P)=F(T,V)+PV`$.

Comparison of the Gibbs free energies as a function of pressure and temperature for the hcp and $`\omega `$ crystal structures allows us to identify the phase transformation line. This approach is extremely demanding because tiny errors in the relative free energies can dramatically alter the phase boundary. Figure 3 shows good agreement between our theoretically calculated $`\alpha \omega `$ phase transformation line with the experimental phase diagramYoung (1991). The calculated transition temperature at zero pressure is approximately 280 K.

At normal pressure we know of no experimental data that shows that the $`\omega `$ phase is the correct zero-temperature ground state. However, it has recently been shown that energy barriers exist for this phase transformationTrinkle et al.  that may make it extremely difficult to verify at low temperatures. Thus the kinetics may be just too sluggish for the system to bring itself into true thermodynamic equilibrium.

We were not able to similarly calculate and include in Fig. 3 the transition into the bcc structure. The open bcc structure appears to be dynamically stabilized by entropy and anharmonic effects. At zero temperature some phonon modes are unstable due to a spontaneous instability of the bcc crystal structure into the $`\omega `$ phase. There is no energy barrier to this continuous transformation. In particular, this causes the longitudinal phonon mode with wave vector $`𝐪=\frac{2}{3}(1,1,1)`$, which distorts the bcc crystal into the $`\omega `$ phase, to be unstable below pressures of roughly 40 GPa. The mechanical instability of this phonon mode at low temperatures is also seen experimentally; the martensitic phase transition from the bcc phase to the hcp phase prevents the high-temperature bcc phase from being quenched to room temperaturePetry et al. (1991). One branch of transverse phonon modes with wave vectors $`(\xi ,\xi ,0)`$, which are related to the bcc to hcp transition, also appears unstable in our treatment. Both of these instabilities are consistent with strictly first-principles calculations of zero-temperature phononsPersson et al. (2000).

Figure 4a shows the relevant segment of the $`(\xi ,\xi ,\xi )`$ branch. The energy of the corresponding pathway from the bcc phase to the $`\omega `$ phase, shown in Fig. 4b, develops a local minima for the bcc phase at pressures of roughly 40 GPa, but at lower pressures the crystal can deform from bcc to $`\omega `$ without encountering an energy barrier. This pathway can be viewed as planes moving in accordance with the wave vector $`𝐪=\frac{2}{3}(1,1,1)`$ or as chains of atoms along the body diagonal moving with respect to each other. Calculations of the energy needed to move such a single chain of atoms indicate that for small displacements the movement is energetically favored, whereas the movement of a single atom is not.

If only the stable phonons are included in the evaluation of Eq. 4, the bcc structure is calculated to be lowest in Gibbs free energy at pressures that are an order of magnitude higher than experimentally observed. A correct calculation requires a more sophisticated treatment that includes phonon-phonon interactions, i.e., anharmonic contributions to the free energy.

In conclusion, we have for the first time used a completely theoretical approach, without any experimental input, to calculate the Gibbs free energies of competing crystal phases in titanium. Comparison of the hcp and $`\omega `$ Gibbs free energies results in a pressure- and temperature-dependence of the structural phase transition in good agreement with experiment and serves to explain how finite-temperature entropy effects stabilize the experimentally found room-temperature hcp crystal structure over the theoretically found zero-temperature $`\omega `$-phase ground state.

We thank Matthias Graf, Carl Greef, and Dallas Trinkle for helpful and encouraging discussions. This research is supported by the Department of Energy under contract W-7405-ENG-36. All FLAPW calculations were performed using the Wien97 packageBlaha et al. (1999).


text:
# Synchronization of delayed systems in the presence of delay time modulation

## Abstract

We investigate synchronization in the presence of delay time modulation for application to communication. We have observed that the robust synchronization is established by a common delay signal and its threshold is presented using Lyapunov exponents analysis. The influence of the delay time modulation in chaotic oscillators is also discussed.

, , , , , and

Synchronization phenomenon was discovered by Huygens in two coupled pendulum clocks. The first description of synchronization in chaotic systems by Fujisaka and Yamada  gave new rise to scientific attention to the phenomenon and, over the past decade, synchronization in chaotic oscillators has been rigorously investigated not only for understanding of its role in nonlinear systems but for its potential applications in various fields . Synchronization in chaotic oscillators is characterized by the fact that transverse variable $`𝐗=𝐱_1𝐱_2`$ converges to zero in the limit of $`t\mathrm{}`$, while longitudinal variable $`𝐱=𝐱_1+𝐱_2`$ is still in a chaotic state. Up to now the kinds of synchronization  that have been observed include complete synchronization in identical systems and phase synchronization (PS), lag synchronization, and frequency synchronization in slightly detuned systems, while generalized synchronization and generalized phase synchronization have been described in systems with different dynamics.

One of the most important fields where chaos has been practically applied is in communication. The major concern in this field was that an encoded message is vulnerable to extraction by nonlinear dynamic forecasting , when it is masked by the signal from the low-dimensional chaotic system. Then it becomes essential to develop high-dimensional chaotic systems having a multiple number of positive Lyapunov exponents to implement a communication system. In this regard, one time-delayed system received a lot of attention: that is,

$$\dot{x}=f(x(t),x(t\tau _0)),$$
(1)
where $`\tau _0`$ is a delay time . The time-delayed system exhibits intriguing characteristics with increasing $`\tau _0`$. Despite a small number of physical variables, the embedding dimension and number of positive Lyapunov exponents increases as the delay time increases, and the system eventually transits to hyperchaos. For this reason, a time-delayed system was considered to be the most suitable candidate for implementing communication systems.

However, it was newly discussed the communication based on a time-delayed system is still vulnerable because the delay time $`\tau _0`$ can be exposed by several measures, e.g., autocorrelation , filling factor, and one step prediction error, etc . If the delay time is known, the time-delayed system collapses to the simple manifold in the $`(\dot{x},x(t),x(t\tau _0))`$ space . So the time-delayed system becomes quite a simple one to someone who knows the delay time $`\tau _0`$, and that message encoded by the chaotic signal can be extracted by the common attack methods .

Accordingly, the most crucial questions to address regarding to the application of time-delayed system to communications come to be ”Is there any scheme which hides the delay time from an eavesdropper”, and ”Is it possible to find out a robust synchronization regime under the scheme”. In this paper, we consider a generalization of time-delayed systems in which the delay time is not constant but modulated in time, which was introduced to study the general property of the delay differential system . Applying the generalization to two coupled chaotic systems, we focus our discussion on answering the two questions above. Thus, we would show that there exists a robust synchronization regime even if the delay time is modulated and that the imprints of delay time in time series are completely wiped out by the modulation of delay time.

First, to explain the delay time modulation (DTM), we consider a single logistic map with time-delayed feedback, and observe what an effect the modulated delay time $`\tau `$ has on the characteristics of the time series: $`x(t+1)=\lambda \stackrel{~}{x}(t)(1\stackrel{~}{x}(t))`$, where $`\stackrel{~}{x}(t)=(1\beta )x(t)+\beta x(t\tau )`$ with $`\tau (t)=g(t,\tau _0)`$. The modulation function of the delay time is given by:

$$\tau =g(t,\tau _0).$$
(2)
If we take $`g(t,\tau _0)`$ as a constant, the system returns to the usual time-delayed system. As we shall see, the nontrivial choice of the modulation function not only has a crucial effect on the characteristics of synchronization, but also remarkably transforms the characteristics of the time series. Figure 1 shows autocorrelations and filling factors  for two different modulations. In the case of the constant delay time (Fig. 1 (a) and (c)), one can easily notice that the peaks indicated by arrows exactly correspond to the multiples of the delay time $`\tau _0`$. On the contrary, when the modulation is sinusoidal (Fig. 1(b) and (d)) the imprints of delay time are completely wiped out, and the detection of delay time seems to be impossible. We emphasize that the disappearance of the imprints of delay time is caused by the fact that the delay time is changing in time. Accordingly, Figure 1 clearly shows that the results do not depend on the method of the detection.

To analyze synchronous behaviors, we consider a system of two coupled Lorenz oscillators with common time-delayed feedback:

$$\begin{array}{ccc}\dot{x}_1\hfill & =\hfill & \sigma (y_1\stackrel{~}{x}_1),\hfill \\ \dot{y}_1\hfill & =\hfill & \stackrel{~}{x}_1z_1+r\stackrel{~}{x}_1y_1,\hfill \\ \dot{z}_1\hfill & =\hfill & \stackrel{~}{x}_1y_1bz_1,\hfill \\ & & \\ \tau \hfill & =\hfill & g(t,\tau _0),\hfill \end{array}\}\text{Transmitter}$$
(3)
$$\begin{array}{ccc}\dot{x}_2\hfill & =\hfill & \sigma (y_2\stackrel{~}{x}_2),\hfill \\ \dot{y}_2\hfill & =\hfill & \stackrel{~}{x}_2z_2+r\stackrel{~}{x}_2y_2,\hfill \\ \dot{z}_2\hfill & =\hfill & \stackrel{~}{x}_2y_2bz_2,\hfill \end{array}\}\text{Receiver }$$
(4)
where $`\sigma =10.0`$, $`r=28.0`$, and $`b=8/3`$ and $`\stackrel{~}{x}_1=(1\beta )x_1(t)+\beta y_1(t\tau )`$ and $`\stackrel{~}{x}_2=(1\beta )x_2(t)+\beta y_1(t\tau )`$. Here we take the modulation function as $`g(t,\tau _0)=\frac{0.95\tau _0}{2}\mathrm{sin}(\omega t)+\frac{\tau _0}{2}`$ and $`\tau _0=60`$. We note that the systems are coupled by common signal $`\beta y_1(t\tau )`$. Figure 2 shows the trajectories in the phase space with time-delayed feedback. One can see that the shape of the attractor is so seriously deformed that it is impossible to identify the original shape of the Lorenz attractor. This deformation of the attractor is the typical characteristic of the time-delayed system which stems from the transition to hyperchaos by time-delayed feedbacks .

The synchronization transition can be easily detected by observing the relative motion of the two oscillators. The equation of motion describing the relative motion is obtained from Eqs. (3)-(4):

$`\dot{X}`$ $`=`$ $`\sigma (Y(1\beta )X),`$
$`\dot{Y}`$ $`=`$ $`(1\beta )(x_1Z+Xz_1XZ)`$ (5)
$`\beta y_1(t\tau )Zr(1\beta )XY,`$
$`\dot{Z}`$ $`=`$ $`(1\beta )(x_1Y+Xy_1XY)`$
$`+\beta y_1(t\tau )ZbZ,`$
where $`X=x_1x_2`$, $`Y=y_1y_2`$, and $`Z=z_1z_2`$. The above equations are nonautonomous in themselves and, so, by iterating the above equation with Eqs. (3)-(4), we obtain the transverse Lyapunov exponents (TLEs)  which enable us to determine the synchronization threshold.

Figure 3 shows two largest TLEs as a function of $`\beta `$ and $`\omega `$. One can clearly see the threshold to CS on $`(\beta ,\omega )`$ plane of Fig. 3 (a) where the largest TLE becomes negative. Also the threshold to PS, where the vanishing TLE becomes negative, is plotted on $`(\beta ,\omega )`$ plane of Fig. 3 (b). The synchronization regimes are presented in Fig. 4 as a contour plot. Contrary to a common intuition, the CS regime becomes wider as the modulation frequency $`\omega `$ increases, even though the regime is smeared with an amplitude death (AD) . The amplitude death phenomenon due to DTM in our study is similar to the transition to periodic orbit in Ref.  in the sence that in both cases the systems lose their exponential instability by the interaction. However different apect is that our system converges to the fixed point without oscillation instead of the periodic orbit. These effects show that DTM can play a role of stabilizing scheme such as chaos control. The appearance of the wider synchronization regime in Fig. 4 implies that the scheme will be effective to implement real communication systems, since the higher modulation frequency widens the synchronization regime and it would enable us to construct more of communication channels.

Figure 5 shows the temporal behaviors of the difference between $`x_1`$ and $`x_2`$ and the modulated delay time $`\tau `$. The temporal behaviors below (point A in Fig. 4) and above (point B in Fig. 4) the threshold are presented in Fig. 5 (a) and Fig. 5 (b), respectively. They clearly show the threshold behaviors of CS at two points. The delay time $`\tau `$ is modulated as shown in Fig. 5 (c). It is interesting to see the autocorrelations and filling factor for different modulations. For the constant delay time, the delay time is exactly indicated by the position of the dips pointed by arrows in Fig. 6 (a) and (c). In the case of DTM (Fig. 6 (b) and (d)), the dips disappear.

In general introducing a $`d`$-dimensional chaotic oscillator, we can consider a chaotic modulation of delay time: $`\tau =g(𝐱,\tau _0)`$ where $`\dot{𝐱}=\frac{1}{ϵ}𝐅(𝐱)`$ and $`ϵ`$ is a time scaling parameter to control the rate of modulation like $`\omega `$ in sinusoidal modulation. We can expect that the chaotic modulation will not work a fundamental change in contrast with the sinusoidal modulation, because the delayed feedback $`x(t\tau )`$ hardly correlated with $`x(t)`$ already even in the sinusoidal modulation as shown in Fig. 6 (b). However, in the application to communication, the chaotic modulation would enhance the security of the system because it may increase the complexity of the attractor.

It is worth discussing the physical implementation of the systems described by Eq. (3)-(4). We think that such a situation can be observed in optical systems. For example, one can consider two laser systems coupled by reflected laser beam from a vibrating mirror. In that case, Lorenz equation describes a laser system and the vibrating mirror plays a role of modulation function.

In conclusion, we have studied the effects of DTM in chaotic oscillators. We have observed that the imprints due to the constant time delay are completely wiped out by the modulation in the chaotic map and flow. We have also demonstrated that the robust synchronization can be established in the presence of DTM. We expect this effects to be useful in constructing real communication systems.

We thank M. S. Kurdoglyan for helpful discussions. This work is supported by Creative Research Initiatives of the Korean Ministry of Science and Technology.


text:
# Temperature dependence of magnetism near defects in SrB6.

## Abstract

The T-dependence of magnetic moments in SrB6 is studied through spin-polarized band calculations for a supercell of Sr27B156 containing a B6 vacancy. The magnetic moment decays rather quickly with T despite the fact that only electronic Fermi-Dirac effects are included. This result and the T-dependence of moments near a La impurity can hardly explain the reports of a very high Curie temperature in hexaborides. but suggest that the magnetism is caused by some other type of impurity.

Early reports of weak ferromagnetism (FM) in lightly doped CaB6 or SrB6, associated with an unusually high Curie temperature (TC) , have been followed by several theoretical and experimental studies. For instance, it has been speculated that the FM could be of excitonic origin , be a spontaneous polarization of the dilute electron gas , or be caused of special conditions around impurities and thereby lead to band magnetism . The band structure of the pure hexaboride , which is metallic when calculated within the density functional (DF) approach, has been questioned, since one GW-calculation found a gap at the Fermi energy, $`E_F`$ . These results are in contrast to a recent independent GW-calculation by Kino et al, where the ground state is metallic . Also, there are conflicting results concerning the interpretation of the experimental results on the weak magnetism. Quantum oscillations indicate that the FM state is not intrinsic to the bulk, but that impurities might be important . This conclusion was supported by Mori and Otani  who found that the FM signal could be removed by a treatment of the crystals, and suggested that that the FM is due to iron impurities. Other defects such as vacancies have been suggested as well . Possible vacancies of whole B6-clusters were anticipated in calculations by Monnier and Delley  who showed that a sizable FM moment surrounded such a vacancy in a supercell of CaB6. They suggest that the observation of weak ferromagnetism in dilute doped hexaborides is caused either by a surface or interface effect, since the B6 vacancy is modeling the natural cleavage planes in polycrystals, or by some amount of B6 vacancies in pure crystals.

In the present letter, by use of the linear Muffin-Tin orbital (LMTO) band method in the local spin-density approximation, we study the variation of the magnetic moment around a B6 vacancy in Sr6 as function of temperature. This is done by self-consistent spin-polarized calculations, in which T-dependent electronic excitations are modeled by the Fermi-Dirac distribution, but no thermal disorder of the lattice is taken into account. As in ref.  we consider a 3$`x`$3$`x`$3 supercell with one B6 vacancy in SrB6 so that the total cell is Sr27B156. The missing B6 atoms form empty MT-spheres, and the calculations include additional empty spheres in the open part structure. The cell contains 19 nonequivalent MT-spheres, 4 Sr, 8 B and 7 empty ones. Other details of the calculation are as in ref. . The density-of-states (DOS) functions near $`E_F`$ for the non-polarized case, are shown in Fig. 1 for two sets of k-points (4 and 10) in the irreducible part of the Brillouin zone. It is seen that the essential features of the DOS are developed already when using 4 k-points, with only small differences in the fine details. In particular, the DOS peak at $`E_F`$ is similar in the two sets of k-points, so that the criterion for Stoner magnetism should come out well already when using only 4 k-points. The Stoner factor is calculated to be 1.05, just above the limit for ferromagnetism, 1.0, indicated by the horizontal line in fig. 1. The broken line shows the partial DOS from the six equivalent B-atoms closest to the B6-vacancy. It is seen that roughly 1/3 of the total DOS at $`E_F`$ comes from these sites although they represent only about 3 percent of the total number of atoms in the cell. In the spin-polarized calculations (at low T), when using 4 k-points, one finds a magnetic moment, $`m=`$2.1 $`\mu _B`$ per cell, in agreement with ref. . This moment decreases to about 1.8 $`\mu _B`$ when the number of k-points is increased to 10 for the same temperature. This type of magnetism is a standard Stoner magnetism , where the large DOS on the 6 B-atoms nearest to the vacancy contribute most to the Stoner instability, cf. Fig 1. Similarly, about 40 percent of the total moment in the spin-polarized calculations comes from these six B-sites closest to the vacancy.

Spin-polarized calculations at different temperatures in the Fermi-Dirac distribution, show that the moment goes to zero when T is increased to about 700 K. The two sets of k-points have been used at the lowest and highest T, and both sets give a vanishing moment at T=700K. Calculations for intermediate T, using 4 k-points, show first a gradual decrease of $`m`$ at low T, while the curve drops sharply as T is approaching 700 K, as shown in Fig 2. This T-dependence is understood from the fact that $`E_F`$ falls on a narrow peak in the DOS, so that thermal smearing, due to electronic excitations, can reduce the effective DOS and the Stoner factor. Since TC in the doped hexaborides is of the same order, 700 K, it is tempting to assign the observed TC to the mechanism of Stoner magnetism due to the B6 vacancy. However, the typical behavior for ferromagnets is that the temperature at which S goes below the critical value for ferromagnetism, is generally much larger than the Curie temperature (TC) for real materials . By taking into account additional smearing effects on the DOS coming from thermally induced disorder of the structure one obtains more realistic values of TC in mean-field calculations . The atomic displacements from thermal vibrations at 700 K can be of the order 0.25 a.u. (corresponding to 7-8 percent of the B-B distance) for a material with a Debye temperature of 350 K, while large smearing effects on DOS peaks appear already at a lower degree of disorder . The present calculations include only the T-dependence via the Fermi-Dirac distribution and the effects from thermal disorder are expected to reduce the moment further. Therefore, it is difficult to understand the high TC in the real materials from the T-dependence of the moments within the Stoner mechanism. One possibility would be that the local structure around a B6 vacancy somehow can escape effects of thermal disorder, in which case it is justified to consider only the effect Fermi-Dirac smearing. The position of $`E_F`$ relative to the peak appears to be at the optimal place for the largest possible Stoner factor for this supercell, and an increased height of the peak at larger T would be opposite to the normal T-dependence of the DOS.

Another mechanism leading to magnetism is due to a gain of Coulomb energy in addition to the exchange energy, by means of charge transfers. This mechanism can be activated if there is a large derivative $`N^{}`$ of the DOS of an impurity band, leading to very small moments near the impurity . This difference with the normal Stoner mechanism, suggests that this weak magnetism could resist better to temperature and broadening of the DOS by disorder, since $`N^{}`$ can be large within a wider energy interval compared to that of a peak in the DOS. However, the results from calculations for a La impurity in a 3$`x`$3$`x`$3 cell give only a partial support to this expectation . The DOS of this cell is such that $`N^{}`$ is maximal a bit below $`E_F`$, and even if there is a small moment at 600 K, it tends to zero at low T. It is difficult to have a sizable moment both at low and high T. By optimizing the position of $`N^{}`$ relative to $`E_F`$ (and the amplitude of $`N`$) by a combination of the virtual crystal approximation and a non-conventional choice of linearization energies in the band calculation, it is possible to have a weak moment within some range of T. Such band results are not ab-initio, but they are useful to show a correlation of moment and charge transfers. They also show that magnetism within this mechanism depend on fine details of the electronic structure near the impurity, which in turn suggest that such magnetism should be a result of circumstances. Furthermore, it is expected that thermal disorder should be as important for this type of magnetism as it is for the standard Stoner magnetism.

In conclusion, the band results confirm that weak magnetism is possible, either as standard Stoner magnetism around a B6 vacancy or assisted by Coulomb energies near a La-impurity, but it is more difficult to understand the high $`T_C`$ values from both of these mechanisms. These results and recent experimental reports  motivate similar electronic structure studies of the properties near iron impurities in these materials. However, the exact state of iron contamination, impurity sites or clustering on surfaces, is yet not clear.