text:
# Spallation Neutron Production by 0.8, 1.2 and 1.6 GeV Protons on various Targets.

## I Introduction

Large numbers of neutrons can be produced through spallation reactions induced by an intermediate energy (around 1 GeV) proton accelerator on a heavy element target. With the progress in high intensity accelerators it is now possible to conceive spallation sources that could compete with high flux reactors. Several spallation sources for solid state and material physics are under construction or study in the USA (SNS SNS ), in Europe (SINQ SINQ , ESS ESS ) and in Japan (NSP NSP ). Spallation neutrons can also be used in Accelerator Driven Systems (ADS) to drive sub-critical reactors, in which long-lived nuclear waste could be burnt Bow ; tak  or energy produced Rub . All these systems have in common a spallation target made of a heavy, either solid (W, Pb) or liquid (Hg, Pb, Pb-Bi eutectic) metal in a container (generally steel) which is separated from the vacuum of the accelerator by a thin window.

A detailed engineering design of a spallation target needs a precise optimisation of its performances in terms of useful neutron production and a proper assessment of specific problems likely to occur in such systems, like induced radioactivity, radiation damage in target, window or structure materials, additional required shielding due to the presence of high energy neutrons, etc… This can be done by using Monte-Carlo transport codes describing the interaction and transport of all the particles created in nuclear reactions occurring inside the system to be designed. Generally, a high energy transport code (often based on the HETC code from HETC ), in which elementary interactions are generated by nuclear physics models, is coupled below 20 MeV to a neutron transport code like MCNP MCNP  that uses evaluated data files. Although the spallation mechanism has been known for many years, the models used in such codes, Intra-Nuclear Cascade followed by evaporation-fission, have never been really validated on experimental data and large discrepancies remain both between experimental data and model predictions and between different models. This was particurlarly obvious from the OECD/NEA intercomparisons ICC1 ; ICC2 ; ICC3  of these codes, regarding neutron and residue production. This led to the conclusion that many improvements of the models are still needed but also that there was a lack of experimental data to make a good validation, especially above $`800`$ MeV. Among the needed data, the energy and angular distributions of spallation produced neutrons are essential for model probing: the high energy part of the spectrum allows the testing of the intra-nuclear cascade while the low energy part of the spectrum is sensitive to the excitation energy at the end of the intra-nuclear cascade stage and the evaporation model. They are also important to optimize the target geometry since secondary particles contribute to the propagation of the inter-nuclear cascade in a thick target and high energy neutrons are responsible for radiation damage in target and structural materials.

During the last years, a wide effort has been made in several laboratories to measure spallation data regarding neutron multiplicity distributions Pie ; Led ; Lot ; Hil , light charged particles NESSI  and heavy residues Mic ; Wla ; TE ; GSI-Au  in order to establish a base for the test and validation of the spallation physics models. We report in this paper on neutron production double-differential cross-sections, measured at the SATURNE synchrotron, induced by $`0.8`$, $`1.2`$ and $`1.6`$ GeV protons impinging on different targets. The experimental setup, although already discussed in detail in Bor ; Mar , is presented in section 1 while the results are displayed in section 2. The aim was to measure at one energy (1.2 GeV) neutron spectra on nuclei representative of different parts of the periodic table of elements and at the same time corresponding to materials used in targets or structures of accelerator driven systems: Al, Fe, Zr, W, Pb and Th. At the two other energies only Fe and Pb targets were studied. Angular distributions covering $`0^{}`$ to $`160^{}`$ were obtained and allowed the determination of average neutron multiplicities per reaction above 2 MeV. The Pb data have already been published in a letter XL  but are again reported here for the sake of completeness. Section 3 is devoted to the comparison of the data to different intra-nuclear cascade models.

## II Experimental apparatus

The slow extraction of the beam delivered by the Saturne synchrotron did not allow conventional time-of-flight measurement using the HF signal of the accelerator. Therefore, the time-of-flight had to be measured between the incident proton passing through a thin scintillator placed in the beam and the detected neutron. For the highest neutron velocities this method becomes highly imprecise due to the limited available flight path, and a measurement of the (n,p) charge exchange by a spectrometer technique was adopted. These two methods had previously been tested in a first experiment in which neutrons were measured at 0 and have been respectively described in detail in two previous papers Bor ; Mar . Then, a new experimental area was designed to allow complete angular distribution measurements. The scheme of this experimental area is shown in Fig. 1. The beam comes from the left, hits the studied target and is deflected to a beam stopper (composed of lead and tungsten blocks) by a dipole magnet. To ensure the detection of the neutrons emitted by the target and not those coming from the beam stop, a large shielding of heavy concrete was built around the target. Because the most energetic neutrons are emitted in the forward direction this shielding is thicker between 0 and 90 (3.5 meters) than at backward angles (2.5 meters). The shielding is pierced by 12 circular holes at 0 and every 15 from 10 to 160. Each channel is composed of two consecutive cylinders of respectively $`110mm`$ and $`136mm`$ in diameter with a length of respectively $`1600mm`$ and $`1800mm`$ at forward angles and $`1100mm`$ and $`1400mm`$ at backward angles. The solid angle is determined by the size of the neutron detectors, not the collimator. The two detection methods were used in dedicated runs due to the fact that they required different beam intensities: between $`10^9`$ and $`10^{11}`$ particles per second for the spectrometer (because of the low neutron-to-proton charge exchange efficiency) while the time-of-flight was limited to less than $`10^6`$ particles per second by the use of the in-beam scintillator. The spectrometer was not used at angles larger than 85 since very high energy neutrons are expected only at forward angles.

### II.1 The time-of-flight method.

The Saturne synchrotron was delivering continuous beam during spills of $`500ms`$ with a repetition frequency depending on the beam energy (typically 1.5 seconds for 800 MeV protons in this experiment). A classic time-of-flight measurement between the neutron detector and a HF signal of the accelerator was thus not possible. We just recall here the principle of the method and the modifications performed since the first experiment at $`0^{}`$ Bor .

The time-of-flight is measured between the incident proton, tagged by a thin plastic scintillator (SC), and a neutron sensitive NE213 liquid scintillator (see Fig 2). The beam intensity is fixed at a maximum of 106 particles per second so that individual incident protons could be counted by SC. The target-detector distance depends on the angle but is about 8.5 meters. Up to ten angles can be explored simultaneously using neutron detectors composed of a cylindrical cell of NE213 liquid scintillator coupled to a photomultiplier. Six of them are cells of the multi-detector DEMON Til  and the other four (called DENSE) are smaller detectors, optimized for low energy measurements. The latter detectors allow energy measurements with a reasonable precision from 2 to 14 MeV, while the DEMON cells are used between 4 and 400 MeV. The characteristics of the DEMON and DENSE detectors are given in table 1.

The DEMON detectors are placed in a shielding of paraffin loaded with Borax and Lithium to reduce the background. The DENSE detectors, smaller and thus less sensitive, do not need such a protection. The energy threshold of the detectors is adjusted using the Compton edge of the gamma rays delivered by a 22Na and a 137Cs radioactive source. The detection thresholds of the DENSE and DEMON are tuned to 1.0 and 1.9 MeV respectively. This allows measurements with a sufficiently well defined efficiency above 2 and 4 MeV respectively.

These detectors are sensitive to neutrons, $`\gamma `$-rays and charged particles. A plastic scintillator NE102 placed in front of each counter (see Fig 2) tags events induced by a charged particle. The neutron-gamma discrimination is performed by a pulse shape analysis done as follows Til : the charge delivered by the photomultiplier is measured by a QDC 1612F during two different time intervals, a prompt one (125 ns long) and a delayed one (185 ns long delayed by 65 ns) giving two charge values (Qf) and (Qs) respectively. A bidimensional spectra Qf vs Qs allows the separation of neutrons and gammas (see Fig. 3).

The neutron detector starts a gate of 500 ns, longer than the time-of-flight of 2 MeV neutrons (440 ns to cover the 8.5 meters from the target to the detector). The signal from the scintillator (SC) in the beam is delayed by 500 ns and is the stop of the time measurement. Due to the beam intensity, although limited to less than $`10^6p/s`$, several protons can be detected during the gate, but only one has induced the detected neutron. For a common start (the neutron) up to ten stop signals are converted and stored with a multistop module (LECROY TDC 3377). The sum of these time spectra contains the real and accidental events. The accidental contribution is determined by the time measurement of unmatched start and stop signals. The background effect is taken into account by a measurement with an empty frame at the location of the target.

The knowledge of the neutron detector efficiency being crucial for this experiment, measurements and calculations have been performed to determine it over the whole energy range. From 2 to 17 MeV, measurements were made at the Bruyères-le-Châtel Van de Graff accelerators as described in Bor . The quasi mono-energetic neutrons are produced by 7Li(p,n)7Be, 3H(p,n)3He, 3H(d,n)4He and 2H(d,n)3He reactions and the efficiency is determined by comparison with a standard detector (full triangles in fig.4 and fig. 5).

At higher energies (30 $``$ E $``$ 100 MeV), experiments have been performed at the TSL Uppsala facilities in Sweden Thun . A neutron beam is produced by 7Li(p,n)7Be reaction in the \[100-180 MeV\] range. The neutron detector efficiency is then measured using n-p elastic scattering and the simultaneous detection of the correlated n and p; the so called associated particle method. For energies from 150 to 800 MeV, the $`d+Be`$ break-up reaction is used at Saturne to produce quasi monoenergetic neutrons. The deuteron beam intensity is measured by activation of a carbon foil Que  and the neutron flux is deduced from $`d+Be`$ break-up cross-sections Lec .

The results are displayed in Fig. 4 and Fig. 5. For the Saturne measurement, the neutron flux attenuation along the 7 to 8.5 meters of flight path in air and in the veto scintillator is the same during the efficiency measurements and the real experiment: therefore it has not to be corrected for. On the contrary, experimental points measured at Bruyères and Uppsala have to be corrected to take into account the difference in neutron flux attenuation due to the difference in distance and the absence of the veto detector. The corrected efficiencies used in the data analysis are represented by the solid line in Fig. 4 and Fig. 5 for DEMON and DENSE detectors respectively. The efficiency calculations performed with O5S O5S  (diamonds) and a modified version of the KSU (triangles) Monte-Carlo codes, agree fully with experimental results. In the original KSU code KSU  some reactions are not taken into account and at 90 MeV the sum of the cross sections over all the inelastic processes gives 165 $`mb`$, that is 90 $`mb`$ less than the global cross section measured by Kellog Kellog . In our modified KSU code, the total inelastic cross-section has been normalised to the Kellog measurements by an appropriate increase of the $`(n,\alpha )`$ light response. This is justified by the fact that the missing reactions in the KSU model produce essentially recoil nuclei like deuterons, tritons, alpha, lithium or boron whose light production is close to that of alphas. The line presented in Fig. 4 is the efficiency used for all DEMON detectors. The DENSE detectors being used only between 2 and 14 MeV, efficiency has been determined by measurements only in the 2-17 MeV range (see Fig. 5). The efficiency was assumed to be the same for all the detectors of the same type.

The beam is monitored by the start scintillator located in front of the target. Uncertainties on the cross-section determination are due to statistical and systematic errors. Systematic errors come mainly from the subtraction of accidental coincidences, gamma rejection and efficiency determination. The error on the latter mainly depends on the absolute calibration procedure used at the different accelerators and are estimated to be 10% for the Bruyeres measurement (knowledge of the standard detector), 4% at Uppsala ($`np`$ cross-section) and 10% at Saturne ($`d+Be`$ cross-section). The values are summarized in table 2.

The neutron energy resolution depends on a time and a geometrical component and is given by :

$$\frac{\sigma _E}{E}=\gamma \left(\gamma +1\right)\left[\left(\frac{\sigma _l}{l}\right)^2+\left(\frac{\sigma _t}{t}\right)^2\right]^{1/2}$$
(1)
with $`\gamma `$ the Lorentz factor
$`\frac{\sigma _t}{t}`$ the time resolution
$`\frac{\sigma _l}{l}`$ the geometrical component
The time resolution ($`\sigma _t=1.5ns`$ and 2.0 ns for DEMON and DENSE detectors respectively) is estimated by the measurement of the FWHM of the prompt gamma peak on the time-of-flight spectra. The geometrical component comes from the target thickness (1 to 3 cm) and from the size of the sensitive area of the detector. The interaction probability being constant as a function of the depth the standard uncertainty is $`\sigma _l=\frac{L}{2\sqrt{3}}`$ ISO . Thus $`\sigma _l=6cm`$ and 1.5 cm for DEMON and DENSE detectors respectively.

The energy resolution is plotted as a function of the energy in figure 6. It appears clearly that above 400 MeV, this time-of-flight method doesn’t allow neutron energy measurement with a resolution better than $`12\%`$. Therefore another complementary detection system has been developed and is described in the next section.

### II.2 Proton recoil spectrometer.

High energy (i.e. above 200 MeV) neutrons are detected using $`(n,p)`$ scattering on a liquid hydrogen converter and detecting the emitted proton in a magnetic spectrometer (Fig 7). We present here only a comprehensive description of the measurement with a special emphasis on the modifications compared to the the first experiment at $`0^{}`$ which is detailed in ref Mar .

The spectrometer is composed of the dipole magnet VENUS, which generates a 0.4 T field, and of 3 multiwire proportional chambers, C1, C2 and C3, of respective active area 20x20cm2, 80x40cm2 and 100x80cm2. Each chamber is composed of 2 sets of wires equipped with PCOS2 electronics allowing the localization in the $`XY`$ plane. The wires are spaced by 1.27mm for C1 and 2mm for C2 and C3.

The acquisition is triggered by the coincidence between the plastic scintillator S1 and the large wall of NE102 plastic scintillators behind VENUS. This wall is made of 20 horizontal slats with a photomultiplier on each side.

A second dipole magnet, CHALUT, deviates in the horizontal plane the charged particles created in the target, in the concrete or in air. The field integral of this magnet was $`0.4Tm`$ during our experiment. The thin plastic scintillator SA tags the possible charged particles remaining in front of the liquid hydrogen target. The spectrometer (hydrogen target and detectors) and CHALUT are placed on platforms which could rotate from 0 to $`85^{}`$ only since, for larger angles, very few neutrons with energies higher than 400 MeV are expected. A single setting of the magnetic field in VENUS ($`0.4T`$) was chosen during the standard measurements. The center of the chamber C3 was shifted by $`40cm`$ to the left with respect to the beam axis in order to optimize the detection of deflected protons.

The trajectory reconstruction of the charged particles emitted by the hydrogen target is deduced from the impact coordinates in C1, C2 and C3. The well known magnetic field inside VENUS gives their momentum. The geometric calibration of the multiwire chambers was done using a low intensity $`800MeV`$ proton beam without the hydrogen target and successive magnetic fields of $`0.2`$, $`0.`$, $`0.2`$ and $`0.4`$ Tesla.

The liquid hydrogen target is a cylinder of 12.8 cm in diameter. It has a useful thickness of $`0.94g/cm^2`$ and is located 8.45 meters from the production target, covering an angular aperture of $`0.43^{}`$. The entrance and outgoing windows are titanium foils, 100 $`\mu m`$ thick .

Several types of charged particles are created in the H2 target through the following processes : np$``$np, np$``$np$`\pi ^o`$, np$``$pp$`\pi ^{}`$, np$``$d$`\pi ^o`$, np$``$2n$`\pi ^+`$. The charged particles with different masses are identified using the biparametric representation of the S1-Wall time-of-flight versus the momentum measured with the wire chambers and VENUS (Fig. 8).

The incident proton beam is monitored by two telescopes viewing a 50 $`\mu m`$ mylar foil placed upstream in the beam ($`20m`$ from the target). The absolute calibration of these telescopes is obtained by a comparison with the activation of a carbon sample Que  measured in a dedicated run. A calibration with activation of Al foils was also done and gives a very comparable result.

The response function of the spectrometer which takes into account the contribution of elastic and inelastic processes arising in the hydrogen target has been measured with quasi-monoenergetic neutrons produced by the break-up of deuterons or 3He beams on a Be target. From various beam energies, neutrons between 0.2 and 1.6 GeV were produced for this calibration. The neutron flux is obtained from the known (n,p) elastic scattering cross-sections CEX  and the normalized response functions are then used to unfold the measured proton spectra. This procedure gives the normalized neutron energy distribution. It is described in detail in ref Mar . The maximum energy available at Saturne (1.6 $`GeV/A`$ for $`{}_{}{}^{3}He`$) and the growing importance of the inelastic processes set a limit to this unfolding procedure.

Systematic errors in this method arise mainly from the beam calibration, the spectrometer response function and the unfolding procedure. The estimations of these 3 errors are given in table 3. They are less than $`10\%`$ at 800 and 1200 MeV but reach $`15.5\%`$ at 1600 MeV due to the increasing inelastic contribution. Error bars on the results presented in this paper take into account only statistical uncertainties except at the very low energies where the increase of uncertainty associated with the proximity of the detection threshold is added.

### II.3 Targets

The same targets were used for both the time-of-flight and spectrometer methods. Because of the low beam intensity imposed by the detection of the incident proton with the SC scintillator, we had to use rather thick targets in order to keep a significant counting rate. They were 3cm diameter cylinders made of natural material with thickness in the centimeter range, shown in Table 4, depending on the elements.

## III Double-differential cross-sections

### III.1 Experimental results

In order to show the consistency of the three different sets of detectors, figs. 10 and 10 display details of the double-differential cross-sections obtained for a Fe target at 1600 MeV. In fig. 10 data obtained at 10 and 25 with the DEMON detectors (filled circles) and the spectrometer (squares) are shown. It can be seen that in the overlap regions, i.e.between 200 and 400 MeV, the data are compatible within the error bars. Actually, 1600 MeV corresponds to the worst case since, as mentioned above, the spectrometer unfolding procedure is approaching its limits of reliability. Other examples of the good agreement between both methods regarding Pb at 800 and 1200 MeV were shown in XL . In fact, for all the measurements, the data from the DEMON detector and spectrometer always agree within less than 15$`\%`$ at 1600 MeV and 10$`\%`$ at lower energies. Concerning the comparison of spectra obtained between 4 and 12 MeV with the DENSE and DEMON detectors, measurements at the same angles were done only for Pb at 1200 MeV and also shown in XL . However, since below 15 MeV neutrons mostly come from an evaporation process which is practically isotropic in the laboratory system, results obtained at near angles can be compared. This has been done in fig. 10 where 25 and 145 DEMON spectra are plotted together with respectively 40 and 160 DENSE ones. Both sets of data are consistent and agree within less than 10%. This appears to be verified whatever the target and the energy. Therefore, a single set of data merging the different measurements by taking their mean values in the overlap regions has been processed and will be shown in the following.

Fig. 11 shows comparisons of our data with previously obtained ones by Amian et al. Ami  and Nakamoto et al. Naka  using time-of-flight. In each case we compare the data at the closest possible angles, adding when appropriate results from DEMON and DENSE detectors at two different angles. Actually, our 10, 25-40, 55, 85-100, 115-130 and 145-160 data are displayed together with the previous ones at respectively 15, 30, 60, 90, 120 and 150. Our Pb measurements at 800 MeV (left) fully agree with Amian ones, as already noticed in XL , while Nakamoto cross-sections are systematically lower at low neutron energy and higher at energies between 10 and 100 MeV. For Fe (center), only data from Ami  are available and we observe a slightly less good agreement between the two works. However, it should be stressed that contrary to Ami  our (and ref Naka ) targets are not really thin and secondary reactions increase the number of low energy neutrons (as discussed below). This effect is visible only below 4 MeV and appears to be larger for a 3cm thick Fe target than for a 2cm thick Pb one, as shown by the simulations in fig. 12. This could explain why we measure more neutrons than Amian et al. at low energies in the case of iron. At 1600 MeV for Pb (fig. 11 right), we have compared our data to ref Naka  data obtained at 1500 MeV. This is possible since, from our 1200 and 1600 MeV measurements, we could infer that cross-sections should differ by less than 5% between 1500 and 1600 MeV, apart from the high energy part of the spectra at very forward angles. As observed at 800 MeV, we get higher cross-sections below 7 MeV and lower ones at intermediate energies. Since the thickness of the targets are the same in both cases, this can be understood only by differences in the neutron detector efficiency determination. As mentioned earlier, at intermediate energies, our experimentally determined efficiency Thun  is higher than the one calculated using the standard KSU code and thus we used a modified version. As far as we know, it was the standard version of KSU that was used to determine the detector efficiency in Naka . This could explain the discrepancy. At high energies, with the spectrometer, we obtain a much better energy resolution than in ref. Naka  (because of their limited flight path) that allows us to distinguish structures due to direct reactions at forward angles.

As already mentioned, the thickness of our targets induces some distortion in the double-differential spectra because of the slowing down of the incident proton and the probability that some of the energetic emitted particles may undergo secondary reactions. The first point was discussed in XL  and results in a slight shift and a broadening of the quasi-elastic and inelastic peaks at very forward angles. The second effect is expected to lead to a depopulation of the high and intermediate energy parts of the spectra and an increase of the number of low energy emitted neutrons. Calculations using the LAHET high-energy transport code system LAHET  (using Bertini as intra-nuclear cascade model and pre-equilibrium) were performed for both a target with the actual geometry and an infinitely thin one in order to assess the order of magnitude of the effect. In fig. 12, results are shown for Pb and Fe at 800 MeV. It can be seen that the difference is very small for the 2cm lead target, for which it is significant only between 2 and 3 MeV, and a little larger for the 3cm iron one. In the latter case, the disappearance of intermediate energy neutrons is also perceptible. Similar results are found with the other targets and at other energies, the effect being maximum for Fe and W.

All the measured angular distributions are presented together with the model calculations in the next section. Data have been taken at 0.8, 1.2 and 1.6 GeV on Pb and Fe targets and at 1.2 GeV on Th, Pb, W, Zr, Fe and Al targets at 0, 10, 25, 55, 85, 130 and 160.

### III.2 Comparison with models

Spallation reactions are generally described by a two step mechanism: a first stage in which successive hard collisions between the incident particle and the individual nucleons of the target nucleus lead to the emission of a few fast nucleons, then, the decay of the excited remnant nucleus by emission of low energy particles or, sometimes for heavy nuclei, by fission. The first step is generally described by Intra-Nuclear Cascade models while evaporation-fission models are used for the second one. Some authors introduce a pre-equilibrium stage between intra-nuclear cascade and de-excitation. In high energy transport codes, the most widely used intra-nuclear cascade model is the old Bertini Ber  one dating from 1963. However, several other models are available, such as the Isabel Isa  and the Cugnon Cugn ; Cugn1  INCL models, which have brought some improvements in the physics. The most widely used evaporation model in the domain of spallation reactions is the Dresner model Dres , usually associated with the Atchison Atch  fission model.

The high energy part of the neutron spectra enables one to directly probe the intra-nuclear cascade models. Low energy neutrons, which are the majority of the neutrons produced in spallation reactions, are emitted during the evaporation process. However, their number mainly depends upon the intra-nuclear cascade stage since the cascade determines the initial excitation energy of the decaying hot residue and, therefore, the number of evaporated particles. Actually, evaporation neutron spectra for a given excitation energy are not expected to depend very much on the evaporation-fission model contrary to light charged particle spectra or to residual nuclei production for which emission barriers and competition between the different decay modes are not so well established. This is why we have made comparisons with different intra-nuclear cascade models in order to test their validity and understand their differences, using always the same evaporation-fission model.

All the calculations discussed in the following have been done with high energy transport codes in which the actual thickness and diameter have been taken into account. In order to have sufficient statistics, calculations were done for angular bins of 5 degrees (except at 0 degrees where it is only 2.5 degrees). The INCL model does not predict a correct total reaction cross-section mainly because the diffuseness of the nuclear surface is not taken into account. Therefore, the INCL calculations were renormalized to the total reaction cross-sections given by the Bertini model which appear to be in very good agreement with experimental values from Bar .

In XL , for lead, we presented calculations performed with the TIERCE TIERCE  high-energy transport code system developed at Bruyères-le-Châtel (which is very similar to LAHET) using either the Bertini or the Cugnon INCL model with the same evaporation-fission model (based on the Dresner-Atchison model). It was shown that, at the three measured energies, the Bertini model was largely overpredicting the experimental data while INCL was giving a rather good agreement. This was ascribed to the higher excitation energy, $`E^{}`$, obtained at the end of the cascade stage with the Bertini calculation than with INCL. This assumption can be verified in fig. 14 where the $`E^{}`$ distribution obtained with Bertini (solid line) is shown to extend to much higher values than INCL (dashed-dotted line) and gives also a higher average value (265 versus 213 MeV). These calculations were performed for thin targets at 1 GeV. The same observations were also made in NESSI  and DF  where a similar plot was shown for p+Au reactions and INCL was found to give the best agreement with the excitation energy deduced from the neutron multiplicity distributions. Several reasons can explain the difference in $`E^{}`$ between the two models : first, INCL leads to the emission of more pions than Bertini. However, the difference in average $`E^{}`$ due to the energy carried away by the pions is only 30 MeV. Second, as mentioned in XL , the Pauli blocking is treated in a different way. In Bertini, only collisions of nucleons with momentum larger than the Fermi momentum are allowed while, in INCL, the actual phase space occupation rate is taken into account. This leads to a less stringent condition, therefore more cascade particles can escape and make the energy remaining in the nucleus lower. This is illustrated in fig. 14 where the $`E^{}`$ distributions with a strict Pauli blocking (as in Bertini) and the standard one are shown. The decrease is obvious.

In the LAHET code system LAHET  it is possible to add after the intra-nuclear cascade stage a pre-equilibrium Preq  stage which is expected to reduce the excitation energy of the nucleus by emission of intermediate energy particles prior to the evaporation. Besides, this is also the recommended option by the LAHET authors Prael . Also available is the Isabel model which can be used only up to 1 GeV. As can be seen in fig. 14, both models lead to excitation energy distributions close to the one found with INCL. Isabel is used with the partial Pauli blocking (recommended) option, which, as in INCL, is supposed to take into account the depletion of the phase space due to the emission of cascade particles. Here, we show calculations performed with both models and the same Dresner-Atchison evaporation-fission at 800 MeV for the Pb and Fe target. In the following, Bertini plus pre-equilibrium will be referred to as BPQ. Fig. 15 presents the calculated neutron spectra compared to the experimental data. It can be observed that, for Pb, the BPQ calculation reproduces very well the data, except at very forward angles and high neutron energies where the peak corresponding to the excitation of the $`\mathrm{\Delta }`$ resonance appears much too high. This is a deficiency of the Bertini intra-nuclear cascade model, already pointed out in Delta  as due to a bad parameterisation of the $`NNN\mathrm{\Delta }`$ reaction angular distribution. The problem does not exist with Isabel. Both models correctly predict the low energy part of the spectra. This can be understood by the respective excitation energy distributions being similar in their extension to the one found with INCL (see fig. 14). The high energy neutrons above 85 are also well reproduced by both calculations but Isabel underestimates cross-sections at backward angles in the intermediate energy region. For iron, Isabel presents the same features as for lead while BPQ now overpredicts low and intermediate energy neutron production at forward angles, indicating that the angular distribution of pre-equilibrium neutrons is probably too much forward-peaked.

At 1200 MeV, the use of Isabel in LAHET being limited to 1 GeV, we compare the data with only BPQ and INCL calculations in figs. 16, 17 and 18, for all the targets. We also performed calculations, which are not shown here, using Bertini without pre-quilibrium. Whatever the target, this model yields too many low energy neutrons, emphasizing that it leads to too high excitation energies. For Th, Pb and W, both BPQ and INCL models give a reasonable agreement with the data, although BPQ tends to slightly overestimate the production of intermediate energies neutrons. As the target becomes lighter, this trend is amplified and BPQ begins to also overpredict low energy cross-sections. This is an indication that the addition of a pre-equilibrium stage after intra-nuclear cascade to decrease the too large excitation energy found in Bertini may not be the proper solution: in fact, it seems difficult to obtain the correct evaporation neutron production without overestimating intermediate energy cross-sections (which are enhanced by pre-equilibrium). On the contrary, INCL reproduces quite well the results for all the targets, proving that the model has a correct mass dependence. Only for the light targets at very backward angles, the high energy neutron production is underpredicted.

At 1600 MeV, fig. 19 displays the results for the Pb and Fe targets. For BPQ the tendencies noticed at 1200 MeV are growing worse: even for Pb, the agreement is not very good between 10 and 40 MeV. Since the high energy part of the spectra is always rather well reproduced (except at 0), this seems to point out a wrong dependence of the pre-equilibrium emission also with incident energy. Here again, INCL gives a satisfactory agreement with the data for both targets.

In summary, we can conclude that INCL is able to globally reproduce the bulk of our data, with some slight discrepancies in the angular distributions. The Bertini model followed by pre-equilibrium, although it is found to be an improvement compared to Bertini alone, works well for Pb at 800 MeV but fails as the energy is increased and the target gets lighter.

## IV Average multiplicities per reaction neutron

Since our double-differential cross-sections nearly cover the full angular range with sufficiently close measurements, it has been possible to infer from the data average neutron multiplicities per reaction above the energy threshold of our detectors. This has been done by interpolating between the measured angles, integrating over 4$`\pi `$ and then dividing by the reaction cross-section taken from Bar ) and the result is shown in tables 5, 6 and 7 for 2 different energy bins, corresponding roughly to evaporation and cascade neutrons respectively. Since we divide by the reaction cross-section, the multiplicities obtained are numbers of neutrons per primary reaction and therefore contain the effect of secondary reactions. The interpolation between angles was done assuming that the angular dependence of the cross-sections was the same as the one calculated by the TIERCE code. This was necessary in particular between 2 and 4 MeV where we have only a few points from the DENSE detectors. The uncertainty on these interpolations was assessed by using different intra-nuclear cascade models in TIERCE and different interpolation procedures. The values given in the tables take into account this uncertainty plus the systematic errors discussed in section 2. The experimental values are compared with the neutron multiplicities given by the two codes, TIERCE-INCL and LAHET-BPQ. For the calculations, 0-2 MeV and total multiplicities are also given. Results are also shown in fig. 20 as a function of the mass of the targets studied at 1200 MeV and in fig. 21 versus incident energy for the iron and lead targets. The comparison of the experimental to calculated average neutron multiplicities confirms in a rather concise way what has been observed in the preceeding section. In all cases, INCL agrees with the data within the error bar while BPQ tends to overpredict 2-20 MeV neutron multiplicities, i.e. evaporation neutron production, especially at 1200 and 1600 MeV. For high energy neutrons (above 20 MeV) the sensitivity to the models is less important. This is most likely because of compensating effects, BPQ predicting more intermediate energy neutrons because of pre-equilibrium while INCL spectra often extend to higher energies. However, a significant deviation from the experiment is found at the two highest energies for iron with BPQ.

Neutron multiplicity distributions have been measured by the NESSI collaboration AL ; DF  using the technique of a gadolinium loaded liquid scintillator tank. The efficiency of this type of detector is large for low energy neutrons and decreases rapidly above 20 MeV. Mean neutron multiplicities were obtained by NESSI at 1200 MeV on W and Pb targets with the same thickness as ours. Although the threshold of our detectors did not allow us to measure neutrons with energies lower than 2 MeV, it is tempting to compare the two results. For Pb at 1200 MeV, the NESSI collaboration measured for a 2cm thick, 15cm diameter target, 14.6 neutrons per reaction which after efficiency correction AL ; DF  amount to 20.3. In DF  it is indicated that the average multiplicities have not been corrected for the additional neutrons coming from secundary reactions in the liquid scintillator. This effect was however investigated and estimated to be 5% for the 2cm thick Pb target (see fig.8 of DF ). Taking this correction into account, the NESSI multiplicity is then 19.3. In our case, we can estimate the total number of neutrons by adding the experimental values measured between 2 and 1200 MeV to the \[0-2 MeV\] multiplicity given by the codes. If we take the average between INCL and BPQ values we find 16.9 neutrons per reaction. The error is estimated to be of the order of 15% taking into account the errors discussed above plus the extrapolation to low energies. While the thickness of the targets in both experiments is the same the diameter is somewhat larger in the NESSI case. We performed a simulation with LAHET to investigate this effect. It appeared that because of secondary reactions of particles emitted sidewards the number of neutrons is 5% larger than with a 3cm diameter target. This means that for a 15cm diameter target we would find 17.8 $`\pm `$ 2.7 neutrons to be compared to 19.3 $`\pm `$ 1.9 if we take an uncertainty of 10% for the NESSI result. The same can be done for the 1cm thick W target. In this case, the additional contribution due to the target diameter is 3%. Therefore, we obtain 16.9 $`\pm `$ 2.5 to be compared to 18.0 $`\pm `$ 1.0 for NESSI AL  after efficiency correction and if we assume, as in Pb, 5% additional neutron from secundary reactions in the liquid scintillator. As it can be seen in fig. 22, we can conclude that our results are compatible with those from NESSI within the error bars.

Also shown in tables 5, 6, 7 are the averaged kinetic energies carried away by the neutrons, extracted from the double-differential cross-sections multiplied by energy using the same procedure as for the multiplicities and compared with the calculations. For the 2-20 MeV bin, conclusions similar to what was stated for multiplicities can be drawn, reflecting the fact that Ex$`M_n`$ is governed by $`M_n`$ in so far as the same evaporation model is used in both calculations and thus gives an identical energy spectrum for the low energy neutrons. For the high energy bin, the compensating effect noticed for the multiplicities seems to be even stronger and, regarding our uncertainties, it is not possible to discriminate between the two models. It is interesting, nevertheless, to remark that these high energy neutrons carry out the major part (from 80% for Th to 98% for Al) of the emitted neutron energy and a large amount (about 30%) of the incident proton energy. In a thick target this will play an important role in the spatial distribution of the energy deposition and particle production.

## V Conclusion

In this paper, we have displayed double-differential cross-sections measured on a wide set of targets and at different energies. This has allowed a comprehensive comparison with some of the high energy models commonly used in high energy transport codes for applications. In particular, we have compared different models describing the first stage of the reaction (intranuclear cascade possibly followed by pre-equilibrium) keeping the same model for the de-excitation process. As we had shown in a previous paper that the Bertini intra-nuclear cascade was leading to too high excitation energies at the end of the first reaction stage, we have tried to use the option, available in LAHET and recommended by its authors, of adding a pre-equilibium stage. This largely improves the predictions of the code. However, discrepancies tend to appear and grow larger as the energy is increased above 800 MeV and as the target becomes lighter. The Isabel model was also tried at 800 MeV. It gave reasonable agreement with the lead data but less good one for iron. Unfortunately, since the use of Isabel is limited to 1 GeV in LAHET, it was not possible to test the energy dependence of this model. Finally, we have shown that the use of the Cugnon intranuclear cascade model, INCL, implemented in the TIERCE code, is able to fairly reproduce the whole bulk of our results. However, it should be recalled that this model still suffers from serious deficiencies mostly due to the fact this it does not treat correctly the diffuseness of the nuclear surface. This was the reason why we had to renormalize the calculation to the correct total reaction cross-section. Also, the sharp surface approximation makes it impossible to have a correct prediction of the most peripheral collisions: this was clearly seen, for instance, in the isotopic distributions of residual nuclei close to the projectile presented in ref. Wla . A new version of the Cugnon model is in progress INCL4  which is expected to solve this problem.

All the data presented will be given to the EXFOR data base or are available on request.

###### Acknowledgements.

This work is partly supported by the HINDAS EU-project, contract FIKW-CT-2000-00031. We would like to thank O.Bersillon for providing us with the TIERCE code and J.Cugnon and C.Volant with the INCL code. We are also grateful to A.Letourneau for fruitful discussions about the NESSI results.


format: text

text:
# Fluctuation-dissipation ratio of a spin glass in the aging regime

## Abstract

We present the first experimental determination of the time autocorrelation $`C(t^{},t)`$ of magnetization in the non-stationary regime of a spin glass. Quantitative comparison with the response, the magnetic susceptibility $`\chi (t^{},t)`$, is made using a new experimental setup allowing both measurements in the same conditions. Clearly, we observe a non-linear fluctuation-dissipation relation between $`C`$ and $`\chi `$, depending weakly on the waiting time $`t^{}`$. Following theoretical developments on mean-field models, and lately on short range models, it is predicted that in the limit of long times, the $`\chi (C)`$ relationship should become independent on $`t^{}`$. A scaling procedure allows us to extrapolate to the limit of long waiting times.

Almost half a century ago, derivation of the fluctuation dissipation theorem (FDT) Callen ; Kubo  which links the response function of a system to its time autocorrelation function, made it possible to work out dynamics from the knowledge of statistical properties at equilibrium. Nevertheless, this progress was limited by severe restrictions. FDT applies only to ergodic systems at equilibrium. Yet, such systems represent a very limited part of natural objects, and there is now a growing interest on non-ergodic systems and on the related challenging problem of the existence of fluctuation dissipation (FD) relations valid in off-equilibrium situations.

A way to extend equilibrium concepts to non-equilibrium situations is to consider systems in which single time dependent quantities (like the average energy) are near equilibrium values though quantities which depend on two times (like the response to a field) are not. Spin glasses Binder  are such systems. They remain strongly non-stationary even when their rate of energy decrease has reached undetectable values. In the absence of any external driving force, they slowly evolve towards equilibrium, but never reach it, even on geological time-scales. In these conditions, FDT is not expected to hold. A quite general FD relation can be written as CuKu2 ; CuKu3  $`R(t^{},t)=\beta X(t^{},t)C(t^{},t)/t^{}`$, where $`R(t^{},t)`$ is the impulse response of an observable to its conjugate field, $`C(t^{},t)`$ the autocorrelation function of the observable and $`\beta =1/k_BT`$. FDT corresponds to $`X=1`$. Determination of $`X`$, the fluctuation-dissipation ratio (FDR), or of an “effective temperature”, $`T_{eff}=T/X`$, is the aim of many recent theoretical studies which predicted a generalization of FDT CuKu1 ; CuKu2 ; CuKu3  in “weak ergodicity breaking” systems Bouchaud . In the asymptotic limit of large times, it is conjectured that the FDR should depend on time only through the correlation function: $`X(t^{},t)=X(C(t^{},t))`$ for $`t^{}`$ (and $`t>t^{}`$) $`\mathrm{}`$. The dependence of $`X`$ on $`C`$ would reflect the level of thermalization of different degrees of freedom within different timescales CuKu3 . Thus, the integrated forms of the FD relation would become $`\chi (t^{},t)=\beta _{C(t^{},t)}^{C(t,t)}X(C)𝑑C`$ (susceptibility function) and $`\sigma (t^{},t)=\beta _0^{C(t^{},t)}X(C)𝑑C`$ (relaxation function). They would depend on $`t`$ and $`t^{}`$ only through the value of $`C`$. The field cooled magnetization would read $`\chi _{FC}=\beta _0^{C(t,t)}X(C)𝑑C=\beta (1_0^1C(X)𝑑X)`$ (in the simplest Ising case with $`C(t,t)=1`$), formally equivalent to the Gibbs equilibrium susceptibility in the Parisi replica symmetry breaking solution for the Sherrington-Kirkpatrick model Mezard , with $`Cq`$ (overlap between pure states) and $`Xx`$ (repartition of overlap). Theoretical attempts, analytical Franz  (with the constraint of stochastic stability) and numerical Marinari  (with the problems of size effects), were made in order to confirm the above properties in short range models. Up to now, experimental investigations correspond only to the quasi-stationary regime Grigera  or are very indirect CuGre .

Here we report the result of an investigation of FD relation in the insulating spin glass CdCr1.7In0.3S4 Alba , an already very well known compound, with $`T_g=16.2K`$. Above $`T_g`$, the susceptibility follows a Curie-Weiss law $`\chi =𝒞/(T\mathrm{\Theta })`$ where $`𝒞`$ corresponds to ferromagnetic clusters of about 50 spins, and $`\mathrm{\Theta }9K`$ Vincent1 . The sample is a powder with grain sizes around $`10\mu m`$, embedded in silicon grease to insure good thermal contact between grains, and compacted into a coil foil cylindrical sample holder $`5mm`$ wide and $`40mm`$ long. The two times dependence of the magnetic relaxation (TRM) of this compound was extensively studied Vincent2 .

In principle, SQUID measurement of magnetic fluctuations is very simple Ocio ; Refregier1 . The difficulty lies in the extreme weakness of the thermodynamic fluctuations (of the order of the response to a field about $`10^7\mathrm{G}`$ in our case). Therefore, the setup is carefully screened against stray fields by superconducting shields, strict precautions are taken to suppress spurious drifts of the SQUID electronics, and the pick-up coil is a third order gradiometer. The result is that the proper noise power spectrum of the system without sample allows time analysis of the magnetic fluctuations signal over up to 2000 s of sample fluctuations with more than 20 dB of signal/noise ratio. Moreover, in the non-stationary regime, the time autocorrelation of magnetic fluctuations $`C(t^{},t)=\frac{1}{N}_i𝐦_i(t^{})𝐦_i(t)`$, where $`𝐦_i`$ is the elementary moment at site $`i`$, must be determined as an ensemble average over a large number of records of the fluctuation signal, each one initiated by a quench from above $`T_g`$ (“birth” of the system). And finally, we want to compare quantitatively correlation and relaxation data. The relaxation function $`\sigma (t^{},t)=\frac{1}{N}_i𝐦_i(t)/𝐇_i`$ is measured by cooling the sample at time zero from above $`T_g`$ to the working temperature in a small field, turning off the field at time $`t^{}`$ and recording the magnetization at further times $`t`$. Using a classical magnetometer with homogeneous field, quantitative comparison between $`C`$ and $`\sigma `$ is almost impossible due to the strong discrepancy between the coupling factors in both experiments. Therefore, we have developed a new bridge setup depicted in Fig.1a, allowing measurements of both fluctuations and response. The pick up (PU) coil of self inductance $`L_0`$ is connected to the input coil of a SQUID, of self inductance $`L_S`$. The whole circuit is superconducting. Relaxation measurements use a small coil $`l`$ inserted in the pick-up circuit, and coupled inductively with mutual inductance $`M`$ to an excitation winding. A current $`I_0`$ injected in the excitation results in a field induced by the PU coil itself ($`1\mathrm{m}\mathrm{G}`$ here, clearly in the linear regime though inhomogeneous), and the sample response is measured by the SQUID. To get rid of the term $`L_0`$, the sample branch is balanced by a similar one without sample, excited oppositely (see Fig.1a). The flux delivered to the PU by an elementary moment $`𝐦`$ at position $`𝐫`$ is given by $`\mathrm{\Phi }(t)=𝐦(𝐫,t)𝐡(𝐫)`$ where $`𝐡`$ is the magnetic field produced by a unit of current flowing in the PU. Flux conservation in the PU circuit results in a current $`I_S`$ flowing in the input coil of the SQUID whose output voltage is $`V_S=GI_S`$. Detailed analysis of the system will be published elsewere. The main features are as follows.

As the fluctuations of elementary moments in the sample are homogeneous and spatially uncorrelated at the scale of the PU, the SQUID output voltage autocorrelation is given by:

$$C_S(t^{},t)=V_S(t^{})V_S(t)=C(t^{},t)Q\frac{G^2}{(L_0+2L_S)^2}$$
(1)
$`Q=_i𝐡^2(𝐫_i)`$ where the index $`i`$ refers to a moment site, is the coupling factor to the PU, including demagnetizing field effects since $`𝐡`$ is the internal field.

The elementary moment response at site $`i`$ is $`R_i(t^{},t)=𝐦_i(t)/𝐡(𝐫_i,t^{})`$.Taking into account that the medium is homogeneous, the relaxation function of the SQUID output voltage is given by

$$\sigma _S(t^{},t)=\frac{V_S(t^{},t)}{I_0}=\sigma (t^{},t)Q\frac{MG}{L_0(L_0+2L_S)}.$$
(2)
Thus, the coupling factor $`Q`$ disappears in the relation between $`C_S`$ and $`\sigma _S`$, independently on the nature and shape of the sample. There remains only the inductance terms $`M`$, $`L_0`$ and $`L_S`$. These being difficult to determine with enough accuracy, absolute calibration was performed using a copper sample of high conductivity, by measuring $`\sigma _S(t^{},t)`$ and $`C_S(t^{},t)`$ — computed by standard fast Fourier transform algorithm — at $`4.2\mathrm{K}`$ (4He boiling temperature at normal pressure): with this ergodic material, the relation between both measured quantities is linear with slope $`A/T`$, where $`A`$ is the  sample independent calibration factor (see Fig.1b). From the knowledge of $`A`$, determined at $`4.2\mathrm{K}`$, the system is equivalent to a thermometer, i.e. the FDT slope is known exactly at any temperature.

In the spin glass sample, $`C_S(t^{},t)`$ and $`\sigma _S(t^{},t)`$ were measured at $`T=0.8T_g`$ after quench from a temperature $`T1.2T_g`$. To get a precise definition of the “birth” time, a minimum value of 100 s was chosen for $`t^{}`$. The autocorrelation was determined from an ensemble of 320 records of $`12000\mathrm{s}`$ of the fluctuation signal. The ensemble averages were computed in each record from the signal at $`t^{}`$, averaged over $`\delta t^{}t^{}/20`$, and the one at $`t`$, averaged over $`\delta t(tt^{})/10`$ — the best compromise allowing a good average convergence still being compatible with the non-stationarity — , and averaging over all records. As there is an arbitrary offset in the SQUID signal, the connected correlation was computed. Nevertheless, this was not enough to suppress the effect of spurious fluctuation modes of period much longer than 2000 s, giving a non-zero average offset on the correlation results. Thus, as a first step, we have plotted all correlation data, taking as the origin the value of $`V_S^2(t^{})`$. Due to the elementary measurement time constant this last term corresponds to an average over $`tt^{}`$ about $`10^2`$s, i.e. a range of $`(tt^{})/t^{}`$ corresponding to stationary regime. Thus, all $`C_S`$ data are shifted by a common offset $`C_0`$. The result is shown in Fig.2a (right sided scale), as a function of $`tt^{}`$ for values of $`t^{}`$ from 100 s to 10000 s. Residual oscillations —and large error bars— for $`t^{}=100s`$ reveal the limit of efficiency of our averaging procedure. Corresponding relaxation data are plotted on Fig.2b. In both results, one can see that the curves merge at low $`tt^{}`$ , meaning that they do not depend on $`t^{}`$ (stationary regime). At $`tt^{}t^{}`$, they strongly depend on $`t^{}`$, the slower decay corresponding to the longer $`t^{}`$.

The correlation offset must be determined. As zero of correlation is unreachable in experimental time, correction of the offset could be obtained from the knowledge of $`C(t,t)`$. Nevertheless, due to clustering, $`C(t,t)`$ depends on temperature and cannot be determined from the high temperature susceptibility. In canonical compounds like 1% Cu:Mn Binder , with negligible clustering, the field-cooled susceptibility is temperature independent in agreement with the Parisi-Toulouse hypothesis Marinari ; Parisi , yielding $`C(t,t)=T_g\chi _{FC}(T)`$. We used a generalization of this relation with the condition that a smooth dependence of $`C(t,t;T)/T`$ must result Parisi2 . This was obtained by using for $`T_g`$ a slightly different value, $`T_g^{}=17.2\mathrm{K}`$. Then, from the value of the calibration factor $`A`$, and writing $`C(t,t;0.8T_g)=17.2\chi _{FC}(0.8T_g)`$, $`C_S(t,t;0.8T_g)`$ can be determined, and suppression of the offset can be performed by using the $`\chi (C)`$ plot, first introduced by Cugliandolo and Kurchan CuKu2 . We plot the normalized susceptibility function $`\stackrel{~}{\chi }(t^{},t)=1\stackrel{~}{\sigma }(t^{},t)`$ where $`\stackrel{~}{\sigma }(t^{},t)=\sigma _S(t^{},t)/\sigma _S(t,t)`$ (note that $`\sigma (t,t)=\chi _{FC}`$) versus normalized autocorrelation $`\stackrel{~}{C}(t^{},t)\stackrel{~}{C}_0=(C_S(t^{},t)C_0)/C_S(t,t;T)`$ for all experimental values of $`t^{}`$. In this graph, the FDT line has slope $`T_g^{}/T`$ and crosses the $`\stackrel{~}{C}`$ axis at $`\stackrel{~}{C}=1`$. On the data, a clear linear range appears at large $`\stackrel{~}{C}`$ (small $`tt^{}`$), displaying the FDT slope with error $`<3\%`$ in the sector $`\stackrel{~}{C}0.47`$. This allows the suppression of the correlation offset by horizontal shift of the data. The result is shown in Fig.3. It is of course based on a rough ansatz on $`C(t^{},t;T)`$ which needs further justifications, but we stress that the induced uncertainty concerns only the position of the zero on the $`\stackrel{~}{C}`$ axis, and not the shape and slope of the curves. With decreasing $`\stackrel{~}{C}`$ (increasing $`tt^{}t^{}`$), the data points depart from the FDT line. Indeed, $`\stackrel{~}{C}(t^{},t\mathrm{})=0`$ and $`\stackrel{~}{\chi }_{FC}=\stackrel{~}{\chi }(t^{},t\mathrm{})=1`$. The mean slope of the off FDT data corresponds to a temperature of about $`30\mathrm{K}`$. This value is far above our annealing temperature, ruling out a simple interpretation in terms of a “fictive” temperature Jackle . Despite the scatter of the results, a tendency for the data at small $`t^{}`$ to depart the FDT line at larger values of $`\stackrel{~}{C}`$ is clear: it is experimentally impossible to fulfill the condition of timescales separation underlying the existing theories . Even if the long $`t^{}`$ limit for $`\stackrel{~}{\chi }(\stackrel{~}{C})`$ does exist, it is not reached in the plot of data in Fig.3 and a $`t^{}`$ dependence of the $`\stackrel{~}{\chi }(\stackrel{~}{C})`$ curves is expected.

The left sided scales in Fig.2a and b correspond to $`\stackrel{~}{C}(t^{},t)`$ and $`\stackrel{~}{\sigma }(t^{},t)`$ respectively. In former works, it was shown that the whole relaxation curves could be scaled as the sum of two contributions, one stationary and one non-stationary Vincent2

$$\stackrel{~}{\sigma }(t^{},t)=(1\mathrm{\Delta })(1+(tt^{})/t_0)^\alpha +\mathrm{\Delta }\phi (\zeta ),$$
(3)
where $`t_0`$ is an elementary time of order $`10^{11}`$ s, $`\phi `$ is a scaling function of an effective time parameter $`\zeta t^{1\mu }t^{1\mu }`$ depending on the sub-aging coefficient $`\mu <1`$ Vincent2 , and $`\alpha `$ can be determined with good precision from the stationary power spectrum of fluctuations $`S(\omega )\omega ^{\alpha 1}`$. The inset in Fig.2b displays the result of the scaling on the relaxation curves with $`\alpha =0.05`$, $`\mathrm{\Delta }=0.21`$ and $`\mu =0.87`$. As shown in the inset of Fig.2a, the scaling works rather well on the autocorrelation curves with the same exponents, but now, $`q_{EA}`$, the Edwards Anderson order parameter, replaces $`\mathrm{\Delta }`$. We get $`q_{EA}=0.37`$. These results show clearly that the stationary part of the dynamics is still important yet in the aging regime, i.e. that the limit of long $`t^{}`$ is not reached within the timescale of our experiments (in fact, timescale separation is realized if $`t^{}\tau `$ where $`\tau `$ is the observation time such that $`C_{stat}(\tau )q_{EA}`$).

If granted, the scaling gives the long time limit of the non-stationary part of the dynamics, allowing a plot of the long times asymptotic non-stationary part of the $`\stackrel{~}{\chi }(\stackrel{~}{C})`$ curve. Of course, here we verify it only over 2 decades of time, up to $`t^{}=10000\mathrm{s}`$, but it was proven to be relevant on TRM up to $`t^{}=100000\mathrm{s}`$ Alba . The dashed line in Fig.3 is obtained by plotting the smoothed curves of aging parts of $`\stackrel{~}{\chi }(\zeta )`$ versus $`\stackrel{~}{C}(\zeta )`$. According to theoretical conjectures, $`d\stackrel{~}{\chi }(\stackrel{~}{C})/d\stackrel{~}{C}`$ would represent the static quantity $`x(q)`$ Franz . One can see that the curve does not point exactly towards $`\stackrel{~}{\chi }=1`$ but about 5% below. Therefore i) either the ansatz used to determine $`C(t,t;T)`$ is not realistic enough, ii) or the time scaling is no longer valid at the very large $`t^{}`$ needed for timescales separation. For future progress, scaling developments outside the strict time range separations at the basis of the “adiabatic cooling” analysis Dotsenko  or the “weak memory” analysis CuKu3  are needed. It seems that such developments are out of the theoretical possibilities for the moment. Toy models Ocio2  presently under development could allow a phenomenological approach of the problem.

In conclusion, we have presented the first experimental determination of the non stationary time autocorrelation of magnetization in a spin glass, an archetype of a complex system. With the help of the time scaling properties of both the relaxation and the autocorrelation, we were able to propose a first experimental approach of a possible generalization of FDT to non-stationary systems. Results at several temperatures are now needed in order to get a complete description of the $`\stackrel{~}{\chi }(\stackrel{~}{C})`$ behavior in the whole temperature range.

We thank J. Hammann, E. Vincent, V. Dupuis, L. F. Cugliandolo, J. Kurchan, D. R. Grempel, M. V. Feigel’man, L. B. Ioffe and particularly G. Parisi for enlightening discussions and critical reading of the manuscript. We are indebted to P. Monod for providing the high conductivity copper sample.


text:
# The REFLEX Cluster Survey: Probing the Mass Distribution in the Universe

## 1. The REFLEX Survey

The REFLEX111The REFLEX Team includes: H. Böhringer (MPE), L. Guzzo (OAB), C.A. Collins (LJMU), P. Schuecker (MPE), G. Chincarini (OAB), R. Cruddace (NRAL), A.C. Edge (Durham), S. De Grandi (OAB), D.M. Neumann (Saclay), T. Reiprich (MPE), S. Schindler (LJMU), P.A. Shaver (ESO), W. Voges (MPE) (ROSAT-ESO Flux-Limited X-ray) cluster survey is currently the largest sample of X-ray selected clusters of galaxies from the ROSAT All-Sky Survey with 1) a statistically homogeneous and fairly well-understood selection function; 2) measured redshifts. The survey contains 452 clusters over the southern celestial hemisphere ($`\delta <2.5^{}`$), at galactic latitudes $`|b_{II}|>20^{}`$ and is more than 90% complete to a flux limit of $`3\times 10^{12}`$ erg s-1 cm-2 (in the ROSAT band, 0.1–2.4 keV). X-ray fluxes and source extensions were re-measured from the RASS with a dedicated algorithm, thus avoiding the limitations of the standard analysis software in the characterisation of extended sources. The details of the whole identification process and a critical discussion of potential sources of incompleteness are presented in Böhringer et al. (2001a). Redshifts for all but 3 REFLEX clusters have been measured during a long Key Programme (1992-2000) using ESO telescopes (e.g. Guzzo et al. 1999), that collected more than 3500 galaxy redshifts over almost 500 X-ray targets.

Fig. 1 shows the distribution on the sky of the 452 clusters in the REFLEX survey, centered on the South Galactic Cap region (the largest contiguous area covered), while Fig. 2 plots their X-ray luminosity against their redshift. Given the survey flux limit (defined by the lower boundary of the point distribution), at large redshifts we are allowed to detect only the very bright, massive clusters. On the other hand, the very large solid angle of the survey allows for an extremely large volume to be explored (4.24 steradians, corresponding to $`8.7\times 10^8`$ h-3 Mpc3 out to $`z=0.3`$ in an $`\mathrm{\Omega }_M=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ cosmology), such that these extremely rare objects on the tail of the cluster X-ray luminosity function have a significant probability to be detected. It is not by chance, in fact, that the most luminous cluster known is still that discovered by the REFLEX survey in 1994, i.e. RXCJ1347.4-1144 (Schindler et al. 1995). Another example is shown in Fig. 4. Such a large volume is also ideal for sampling the very large modes of density fluctuations in the Universe (Schuecker et al. 2001). From the distribution of Fig. 2 one can compute the mean density of clusters as a function of distance, which indicates a very good completeness (i.e. constant density) out to at least $`z=0.2`$, possibly above. It is also evident how one can use the REFLEX survey to select sub-samples of clusters within well-defined ranges in $`L_X`$ (i.e. mass) and redshift. A number of follow-up studies of this kind are indeed ongoing (see Böhringer et al. 2001b for an overview), and more will certainly follow when the catalogue is released in early 2002.

The filled circles in Fig. 2 represent 142 clusters which do not appear in the visually-selected Abell/ACO catalogue (Abell 1958; Abell et al. 1989). It is interesting to see that these clusters are distributed basically at any redshift. This is a demonstration of how the “richness” criterion (the reason why most of these objects did not enter the Abell selection) is a poor indicator of the cluster mass, and how incomplete in mass a purely visually-selected cluster survey could be. We have also performed a direct X-ray analysis of the RASS data at all Abell-ACO cluster positions, measuring their X-ray flux. This gave the rather encouraging result that the REFLEX survey in fact detects all Abell-ACO clusters within its flux and area boundaries222Even more, the REFLEX survey misses only 1 of the so-called Supplementary Abell clusters i.e. the extension to the main ACO catalogue, where those objects that did not meet all original criteria while still looking as bonafide clusters were listed..

## 2. The Cluster X-ray Luminosity Function

The X-ray luminosity function (XLF) of galaxy clusters is an excellent observable incarnation of the cluster mass function, given that in clusters of galaxies X-ray luminosity is well related to the cluster mass (e.g. Reiprich & Böhringer 2002). For this reason, estimating the XLF was one of the prime targets of the REFLEX survey. It was also one of the first results obtained by the project along the way, from the higher–flux RASS1 Bright Sample (De Grandi et al. 1999). The recent XLF from the whole survey (Böhringer et al. 2002) is compared to that early result in the left panel of Fig. 3. The improvement in the error bars of the latest REFLEX XLF (smaller than the size of the dots) gives an immediate idea of its quality. Such an accurate local measurement represents a standard reference of crucial importance for evolutionary studies, i.e. for comparison with XLF measured at high redshift. One such example is shown in the right panel of the same figure, where the REFLEX XLF is plotted together with those from the RDCS (Rosati et al. 1998) and EMSS (Henry et al. 1992) deep surveys. With the REFLEX data pivoting the local abundance, the evidence for mild evolution of the XLF bright end at high $`z`$ can be tested to a high accuracy (e.g. Borgani et al. 2001; Henry 2001).

## 3. The Clustering of REFLEX Clusters

Clusters have a long history as tracers of large-scale structure (see Nichol, this volume) and, in particular, surveys of X-ray clusters have their own specific advantages (see Lahav et al. 1989 and Romer et al. 1994 for early applications, and Borgani & Guzzo 2001, Guzzo 2001 and Henry 2001 for recent reviews). First of all, once the characteristics of the X-ray instrument used are specified, a clean selection function can be constructed. This is by all means not an obvious task for an optically-selected survey. Second, fairly precise predictions can be made from the models for both the mass function and the clustering of objects above a given mass (see e.g. Borgani, this volume and Moscardini et al. 2000). These can be reliably translated into “observer space” in terms of X-ray luminosities and fluxes. These advantages allow one to select samples uniformly by mass over large areas (avoiding also other effects that plague optically-selected samples as line-of-sight projections), and provide the baseline for the REFLEX survey clustering results.

The right panel of Fig. 4 plots the 3D distribution of REFLEX clusters within a radius of 600 h-1 Mpc. Despite the fading with distance due to the flux-limited selection function, a number of superstructures with sizes 100 h-1 Mpc are evident: clusters are clearly still strongly clustered on such scales. This can be quantified by the two-point correlation function, that we plot for the whole flux-limited survey in Fig. 5 (Collins et al. 2000).

The shape of $`\xi _{cc}`$ is tightly constrained by the REFLEX data, and closely reproduces an amplified version of the galaxy two-point correlation function over almost two decades of scales, a classic prediction of biasing theory (Kaiser 1984). In the same figure, we also plot the predictions of two low-density CDM models from Moscardini et al. (2000), which are able to correctly reproduce both the shape and amplitude of the observed REFLEX $`\xi _{cc}`$. We remark that in this case the amplitude is not a free parameter as for galaxy correlation functions, since the bias value for the specific REFLEX selection function can be computed using an appropriate theory for the clustering of massive haloes (e.g. Mo & White 1996 and subsequent refinements), given the fairly straightforward relation between X-ray luminosity and mass (see Borgani & Guzzo 2001 for more details and references). Analogous results have been obtained from the power spectrum, as discussed in the contribution by P. Schuecker to this volume (and in detail by Schuecker et al. 2001).

Another prediction which can be precisely made for a survey like REFLEX is the dependence of the correlation length on the survey flux limit. The expectation values for different CDM variants are compared to the observed correlation length from the full REFLEX in Fig. 6, again indicating that low-$`\mathrm{\Omega }_M`$ models provide the best consistency with observations (Collins et al. 2000).

Redshift-space anisotropies in the clustering pattern of clusters can be evidenced through the bidimensional correlation function $`\xi (r_p,\pi )`$ \[or analogously $`\xi (\sigma ,\pi )`$\]. Classical analyses of Abell samples show strong elongations along the line-of-sight $`\pi `$ direction (e.g. Miller et al. 2001), which are typically not seen in automatic cluster catalogues (e.g. Nichol et al. 1992). These arise by a combination of projection effects, redshift errors and true spatial anisotropies in small–size catalogues. The plot of $`\xi (r_p,\pi )`$ for REFLEX (Fig. 7), in fact shows no evidence for stretching along the line of sight, indicating how these spurious effects are negligible in this survey. On the other hand, the contours are significantly compressed at large separations, a typical signature of streaming motions (proportional to the quantity $`\beta =\mathrm{\Omega }_M^{0.6}/b`$) which is for the first time seen in a cluster sample (Guzzo et al. 2002; see also Padilla & Baugh 2001).

We would like to conclude by mentioning that the statistical quality of the results presented here will significantly improve in the near future, thanks to the extension of REFLEX to fainter fluxes (by a factor $`2`$). The REFLEX-2 sample is in fact being constructed and will bring the total number of clusters over the same area to $`800`$ (see e.g. Böhringer et al. 2001b).

## References

Abell, G.O., 1958, ApJS, 3, 211

Abell, G.O., Corwin, H.G. Jr., & Olowin, R.P., 1989, ApJS, 70, 1

Böhringer, H., et al. (the REFLEX Team), 2001a, A&A, 369, 826

Böhringer, H., et al. (the REFLEX Team), 2001b, The Messenger, 106, in press

Böhringer et al. (the REFLEX Team), 2002, ApJ, in press (astro-ph/0106243)

Borgani S., & Guzzo L., 2001, Nature, 409, 39

Borgani, S., et al., 2001, ApJ, in press (Nov 2001 issue, also astro-ph/0106428)

Collins, C.A., et al. (the REFLEX Team), 2000, MNRAS, 319, 939

De Grandi, S., et al. (the REFLEX Team), 1999, ApJ, 513, L17

Guzzo, L. 2001, in Where’s the Matter? Tracing Dark and Bright Matter with the New Generation of Large Scale Surveys, proc. of meeting held in Marseille (June 2001), M. Treyer & L. Tresse eds., Frontier Group (astro-ph/0111309)

Guzzo, L., et al. (the REFLEX Team), 1999, The Messenger, 95, 27

Guzzo, L., et al. (the REFLEX Team), 2002, in preparation

Henry, J.P., et al., 1992, ApJ, 386, 408

Henry, J.P., 2001, in “AMiBA 2001: High-z Clusters, Missing Baryons, and CMB Polarization”, in press (astro-ph/0109498)

Kaiser, N., 1984, ApJ, 284, L9

Lahav, O., et al., 1989, MNRAS, 238, 881

Miller, C.J., Nichol, R.C., & Batuski, D.J., 2001, ApJ, 555, 68

Mo, H.J. & White S.D.M., 1996, MNRAS, 282, 347

Moscardini, L., et al., 2000, MNRAS, 314, 647

Nichol, R.C., et al., 1992, MNRAS, 255, 21p

Padilla, N.D., Baugh, C.M., 2001, MNRAS, submitted (astro-ph/0104313)

Reiprich, T.H., & Böhringer, H., 2002, ApJ, in press (astro-ph/0111285)

Romer, A.K., et al., 1994, Nature, 372, 75

Rosati, P., et al., 1998, ApJ, 492, L21

Schindler, S.A., et al. (the REFLEX Team), 1995, A&A, 299, L9

Schuecker et al. (the REFLEX Team), 2001, A&A, 368, 86