E Khan

M Grasso

J Margueron

Constraining the nuclear pairing gap with pairing vibrations

Keywords: numbers: 21.60.Jz, 21.65.Cd, 25.40.Hs, 25.60.Je

Pairing interactions with various density dependencies (surface/volume mixing) are constrained with the two-neutron separation energy in the Tin isotopic chain. The response associated with pairing vibrations in very neutron-rich nuclei is sensitive to the density dependence of the pairing interaction. Using the same pairing interaction in nuclear matter and in Tin nuclei, the range of densities where the LDA is valid in the pairing channel is also studied.

I. INTRODUCTION

Studies on pairing effects in both nuclear matter and finite nuclei have known intensified interests in the recent years [1]. There are two main approaches for pairing, depending whether the mean field is based on Gogny finite range interaction or on Skyrme interaction. In the first approach, a similar functional is used in both the particle-hole channel and the pairing channel, although interactions are not exactly the same due to the density dependence of the pairing interaction. In the Skyrme approach, the functionals are meant to be different in the two channels, as witnessed for instance by their density dependence. The use of a different interaction in the particle-hole channel and in the pairing channel has been justified a decade ago [2]; this is for instance the case of employing the Skyrme interaction in the particle-hole channel and a zero range density dependent interaction in the pairing channel. We shall focus on the Skyrme approach: in this case the pairing density functional is difficult to constrain and it has not been possible to derive an universal pairing interaction during past decades, using for instance the odd-even mass staggering on finite nuclei (see e.g. [3,4]). This may indicate the need for another approach, using additional constrains: should the pairing density functional be extended, and are there additional relevant observables to constrain it ?

Nuclear matter could help in constraining the pairing functional. This requires however to bridge nuclei and nuclear matter through LDA in the pairing channel: its condition of validity should be more systematically analysed. It has been recently shown that the two paired neutrons are spatially localised in low density medium which corresponds to the surface of the nucleus [5]. The same conclusion is drawn by analysing the di-neutron configuration in the excited states [6,7], and also performing calculations in low density matter in [8,9], mainly renewing the possibility to link in some cases the nuclear matter and nuclei in the pairing channel through the LDA.

Concomitantly the pairing functional has been extended in order to study the condensation of the Cooper pairs (BEC-BCS crossover) in both symmetric and neutron matter. In nuclear matter the medium polarization increases the pairing gap at low densities in symmetric matter, whereas it reduces the gap in neutron matter, in-dicating an isospin dependence of the pairing functional [10]. The application to finite nuclei of extended pairing density functional have shown the relevance of the LDA in the pairing channel [11].

The pairing functional studies may thus enter in a new era, renewing the method to design the pairing interaction: i) using an isospin dependence of the pairing interaction ii) using eventually the nuclear matter as an additional constrain for the pairing interaction iii) looking for additional observables in nuclei than the odd-even mass staggering to constrain the pairing interaction. Point i) has been investigated in [10,11]. Point ii) requires the validity of the LDA in the pairing channel and shall be addressed in the present work.

In the case of point iii) an interesting observable is pairing vibrations, measured through two-particle transfer. It is well known that the transfer cross section crucially depends on the pairing interaction at work in the transferred paired [START_REF] Broglia | Advances in Nuclear Physics[END_REF]13]. However in the 70-80's the form factor of the transition has never been calculated fully microscopically. The first microscopic calculations has been performed only recently [14], allowing for a strong link between the pairing interaction and pairing vibrations. Several calculations followed [6,15], showing the renewed interest for such studies.

It is therefore meaningful to use pairing vibrations as a complementary observable to the masses, in order to constrain the pairing interaction, and study the implications to the nuclear matter. More generally, the isospin extension of the pairing interaction should be accompanied with additional constraints. One purpose of this work is to evaluate if pairing vibrations could play this role (Section III).

The method is to analyse the sensitivity of pairing vibrations to various pairing interactions which provide the same two-neutron separation energy in Tin isotopes, and evaluate the consequences on the pairing gap in symmetric and neutron matter. On this purpose it is necessary to determine the range of density where the LDA is valid in the pairing channel (Section II). After many years of study, there is still no unambiguous universal pairing functional ranging on the whole nuclear chart, and current efforts are aiming in that direction. The problem may be due to the method used to constrain it, namely comparing the pairing gap with odd-even mass differences, or evaluating the separation energies along a given isotopic chain. It therefore may be useful to consider a more general context: the evaluation of several pairing interactions constrained by odd-even mass difference, on nuclear matter on one side, and on additional observables on the other side, should shed a renewed light on the problem. To achieve this goal it is first necessary to determine the range of validity of the LDA in the pairing channel.

A. Method to determine the functional

The method is the following: we first consider surface and various mixed paring interactions. The parameters are determined so as to describe the two neutron separation energy. Then pairing vibrations are used in order to disentangle between the various pairing interactions (Section III). We choose 124 Sn and 136 Sn nuclei: these are spherical nuclei where pairing vibrations are likely to occur [13]. One is stable and the second has a large neutron excess.

The microscopic calculations for the ground state are based on the Hartree-Fock-Bogoliubov (HFB) model. The Skyrme interaction SLy4 [16] is chosen for the particle-hole channel of the HFB equations. The adopted pairing interaction is the usual zero-range densitydependent interaction

V pair = V 0 1 -η ρ(r) ρ 0 α δ (r 1 -r 2 ) ( 1
)
where η provides the surface/volume character of the interaction. We set α = 1 and ρ 0 = 0.16 fm -3 . The numerical cutoff for the microscopic calculations is given by E max = 60 MeV (in quasiparticle energies) and j max = 15/2. For each value of η, V 0 is chosen to fit the known experimental two-neutron separation energies for Sn isotopes. Surface and mixed interactions have been considered in this work and the used values of (η,V 0 ) are listed in Table I. As an illustration to visualize the features of the calculated pairing effects, we display in Fig. 1 the neutron pairing field for 124 Sn corresponding to the surface η = 1 and the mixed η = 0.35 interactions.

B. Pairing gap in uniform matter

The relation between the pairing gap in uniform matter at a given density and the pairing field at a given radius in nuclei has been explored in Ref. [11]. It has been found that in the case of mixed interactions, the LDA is in good agreement with the full microscopic HFB calculation (differences less than 15% on the pairing field). This might be related to the extension of the Cooper pair which is getting smaller at the surface of nuclei (about 2 fm) compared to that in the interior (about 5-6 fm) [5]. Close to the surface, pairing properties shall not be very different from that of a uniform piece of matter at the same density. It is then interesting to explore the low density properties of the different pairing interactions listed in Table I.

Fig. 2 displays the pairing gap in uniform matter for various pairing interactions. It is observed that the different interactions leads to very different pairing gap at low density while around saturation density, there is a density (ρ=0.11 fm -3 ) at which the pairing gap and pairing strength coincide for the three pairing interactions.

From Fig. 2, two conclusions can be drawn: i) the twoneutron separation energy used to adjust the parameters of the pairing interaction is an observable which provides a strong constrain on the pairing gap localized at the surface of the nuclei: the pairing gap in nuclear matter has been constrained for ρ=0.11 fm -3 , which corresponds to R ≃ 5 fm in Tin nuclei, as showed by the 124 Sn densities displayed on Fig. 1 ii) to better constrain the value of the parameter η, one shall find another observable sensitive to the pairing strength at low density (large radius, experimentally easier to probe). Indeed, in the very external part of the nuclei the pairing strength is very different from one interaction to another. The pure surface pairing interaction predict a pairing gap as high as 8 MeV at low density while the various mixed pairing interaction are grouped below 3 MeV (see Fig. 2). Therefore, one might expect that properties of collective modes sensitive to the external part of the nuclei could be changed by the properties of the pairing interaction at low density. Pair transferred reaction mechanisms like (p,t) or (α, 6 He) which is very surface peaked shall also help in extracting the value of the pairing gap in the external part of the nuclei or equivalently at low density.

C. Validity of the LDA

One aim of this work is to investigate the mapping of the pairing gap between nuclei and nuclear matter, implying the following relation:

∆ ∞ (ρ) = ∆ HF B (ρ(r) = ρ) (2)
Eq. ( 2) corresponds to the LDA in the pairing channel. However, it should be noted that in practice, the LDA might still be valid even if Eq. ( 2) is not verified for very low densities (ρ ≤ 0.1 fm -3 ): the major role is played by densities around the saturation density, where nucleus properties are the most important. asymmetry are taken as the one of 124 Sn. One can see that the LDA is not valid for a pure surface pairing, but becomes valid in the η=0.35 mixed case, for densities lower than 0.1 fm -3 . Therefore, in this case, the low density part of the nuclear matter gap could be constrained by nuclei calculations of the gap: a valid case of the LDA in the pairing channel is showed here quantitatively for the first time. This range of density corresponds to a radius larger than about 5.1 fm in 124 Sn and 5.2 fm in 136 Sn localised at the surface of the nuclei, as can be seen on Fig. 1. In the cases of the other pairing interactions, the LDA underestimates the pairing gap at low density and overestimate it at high density. However, it should be noted that the relative difference of the gaps between the various pairing interactions are respected.

As noted above, to assess in practice the validity of the LDA, the density range around the saturation density should be investigated, where nuclei bulk properties are at work. Therefore a large discrepancy at very low density (ρ ≤ 0.1 fm -3 ), as observed in the case of the surface pairing interaction in Fig . 3, does not preclude on the validity of the LDA. More generally it can be stated that the less important the surface contribution to the pairing interaction, the better the LDA.

III. PAIRING VIBRATIONS

As stated above, it may be useful to consider an additional observable than the separation energy, in order to constrain the pairing interaction, namely its density dependence. There are only few observables which could be relevant to constrain pairing effects. It has been shown that the first 2 + state in nuclei is sensitive to the pairing interaction [17]: both its position and strength depend of the pairing interaction. However this is mainly related to the pairing gap value, which is the same observable extracted from odd-even mass excess. It should be noted than none of these two observables (the first 2 + state and the odd-even mass staggering) can be directly linked to predictions. On one side there is the difficulty to modelize excited states. On the other side, the difficulty is to describe odd nuclei.

Pairing vibrations may be more adequate observable. They can be probed for instance with two neutron transfer in nuclei close to shell closure. We refer to [START_REF] Broglia | Advances in Nuclear Physics[END_REF]13] for details on pairing vibrations. Basically, these modes corresponds to the (collective) filling of subshells, in transition from an A to A+2 nuclei.

With pairing vibrations, pairing effects are probed by 3 ways. The first one is the magnitude of the pairing gap ∆ (average of the pairing field): a large pairing gap implies strength at larger energies, following the formula E 2 ≃(ǫλ) 2 + ∆ 2 . This component is also present in the first 2 + state in the ph response as well as in the odd-even mass staggering. But in the case of the pairing vibrations, there are two additional contributions: first, the transition densities generating the strength are the pairing one, which means that the unperturbed response as well as the perturbed response are sensitive to the impact of the pairing on the wave functions. Finally, the residual interaction, generating the QRPA response, is also sensitive to pairing.

Pairing vibrations are therefore expected to be very sensitive to the pairing interaction. On the other hand, it may also be difficult to disentangle between the three above mentioned However, the first one can be evaluated using the energies of the unperturbed response, the second one by studying the pairing transition densities, and the last one by comparing the unperturbed and the QRPA responses. It should be noted that a related study will also be performed in [START_REF] Matsuo | Proc. of COMEX 3 Conference[END_REF].

A. Method: QRPA in the pp channel

The method is described in [14,17]. Namely the QRPA equations are solved in coordinate space, using the green function formalism. The variation of the generalized density R' is expressed in term of 3 quantities, namely ρ ′ , κ ′ and κ′ , which are written as a column vector:

ρ ′ =    ρ ′ κ ′ κ′    , (3)
where

ρ ′ ij = 0|c † j c i | ′ is the variation of the parti- cle density, κ ′ ij = 0|c j c i | ′ and κ′ ij = 0|c † j c † i | ′
are the fluctuations of the pairing tensor associated to the pairing vibrations and | ′ denotes the change of the ground state wavefunction |0 > due to the external field. In contrast with the RPA where one needs to know only the change of the ph density (ρ ′ ), the variation of the three quantities (3) have to be calculated in the QRPA. In the three dimensional space introduced in Eq. ( 3), the first dimension represents the particle-hole (ph) subspace, the second the particle-particle (pp) one, and the third the hole-hole (hh) one. The response matrix has therefore 9 coupled elements in QRPA, compared to one in the RPA formalism.

The variation of the HFB Hamiltonian is given by:

H ′ = Vρ ′ , ( 4
)
where V is the matrix of the residual interaction expressed in terms of the second derivatives of the HFB energy functional, namely:

V αβ (rσ, r ′ σ ′ ) = ∂ 2 E ∂ρ β (r ′ σ ′ )∂ρ ᾱ(rσ)
, α, β = 1, 2, 3. (5)

In the above equation the notation ᾱ means that whenever α is 2 or 3 then ᾱ is 3 or 2.

The QRPA Green's function G can be used for calculating the strength function associated with the twoparticle transfer from the ground state of a nucleus with A nucleons to the excited states of a nucleus with A+2 nucleons. This strength function is :

S(ω) = - 1 π Im F * (r)G 22 (r, r ′ ; ω)F (r ′ )dr dr ′ (6)
where G denotes the (pp,pp) component of the Green's function and F in the external perturbating field associated with the addition of two particles.

In the QRPA calculations the full HFB quasiparticle spectrum up to 60 MeV is included. These states are used to construct the unperturbed Green's function G 0 . The residual interaction is derived from the two-body force used in HFB according to Eq. ( 5). The contribution given by the velocity-dependent terms of the Skyrme force to the residual interaction is calculated in the Landau-Migdal approximation, which is shown to be accurate [START_REF] Mizuyama | [END_REF]. The strength function for the two-neutron transfer is calculated using Eq. ( 6). The unperturbed Green's function is calculated with an averaging interval equal to 0.15 MeV. All details can be found in Ref. [14].

The response function is calculated for the pp channel. All the calculations are performed in a box of size 22.5 fm. It should be noted that exact continuum treatment is much heavier, especially for nuclei such as Sn isotopes. Moreover the aim is not to study the impact of the continuum treatment (see [14] for such a study). Finally the Sn isotopes under study are far from the drip line, and continuum effects are expected to play a negligible role.

B. Unperturbed response results

The HFB solutions are used in the quasiparticle random-phase approximation (QRPA) scheme to analyze self-consistently the excitations modes associated to the pair transfer reactions. Since we study here two-neutron transfers, we focus on the neutron HFB quasiparticle states that are used to construct the elementary configurations of the excited modes. We work with positiveenergy quasiparticle states. Once calculated the quasiparticle spectrum, it is possible to deduce some properties of the unperturbed response function.

The quasiparticle states with energy less than 6 MeV and an occupation probability ≤ 80 % are presented in Tables II and III for 124 Sn and 136 Sn, respectively. Let us discuss the two cases η = 0.35 and η = 1 (for η = 0.65, results are similar to those obtained with η = 0.35). For 124 Sn, in the case of a mixed pairing interaction, η = 0.35, all the quasiparticle states with energy lower than 5 MeV are totally occupied with the exception of a h 11/2 state at 1.5 MeV which is 42% occupied. This is the only lowenergy state that can contribute to some extent to the excitation mode. The states that are completely empty and can thus contribute more to the excitation are located at higher energies. The first is an f 7/2 state at 5.8 MeV. The others have larger energies (at least 1 MeV more). One can thus expect that the unperturbed response profile starts with a peak at twice 5.8 MeV, i.e. at ∼ 11.6 MeV (with some small contribution at 3 MeV). In the case of a surface interaction, η = 1, again, all the states between 0 and 5 MeV are occupied with the exception of a h 11/2 state at 2.2 MeV (42% of occupation). This time there are several unoccupied states just above 5 MeV, the lowest energy being at 5.4 MeV (p 3/2 state). Hence, the unperturbed response is expected to have some structure starting from ∼ 10.8 MeV with a small contribution at ∼ 4.4 MeV.

For the nucleus 136 Sn the situation is different: there are several low-lying unoccupied states. For η = 0.35 the lowest energy for a completely unoccupied state is 1.9 MeV (p 3/2 state). At 0.8 MeV one also finds a f 7/2 state with 45% of occupation. In the case η = 1 the lowest energy for a totally unoccupied state is 1.7 MeV (p 3/2 state) and a f 7/2 state is found at 1.6 MeV with 32% of occupation. The unperturbed response is expected to start at ∼ 3.8 and 3.2 MeV for η = 0.35 and 1, respectively. In the former case a small contribution at ∼ 1.6 MeV is also expected.

In order to disentangle between the various pairing effects, the unperturbed response in the two neutrons addition mode is first shown on Fig. 4 for 124 Sn. The unperturbed response is built on the HFB single quasiparticle (QP) spectrum, for the three pairing interactions. It should be noted that the spectrum is showed above 10 MeV, because there is only the h 11/2 subshell which can welcome two neutrons to make a low energy state: all the other configurations belong to the next major shell (see Table II), explaining this high energy feature of the spectrum, as stated above. For all the mixed pairing interaction, the unperturbed spectrum is similar, showing that both the single quasiparticle energy and wave functions are close to each other in that case. However, in the procedure. The larger magnitude comes from the wave functions, and will be studied with the QPRA response.

It can be already stated that the QRPA response will also have more strength at lower energy, due to this peculiar feature of the unperturbed spectrum for the pure surface pairing force. Fig. 5 shows the unperturbed response for the two neutron addition mode in 136 Sn. In this case, at the beginning of an open neutron shell several low energies configurations can welcome the two neutrons (see Table III). As in the case of 124 Sn, the response exhibits larger strength at low energy in the specific case of the pure surface pairing interaction, compared to others pairing interaction. This is related to the pairing field profile as shown on Fig. 1. It should be noted that in order to clearly see the effect due to the surface pairing, not only the first 0 + state, but also the energy area of a few MeV above should be explored since the results are different from 0 to 4 MeV on Fig. 5.

C. Perturbed response results

Fig. 6 shows the QRPA response for 124 Sn, with a pure surface and the two mixed interactions. As expected the residual interaction plays a similar role in all the cases, gathering strength and shifting it to lower energy. In the case of 124 Sn, a peak around 9 MeV is the strongest for the surface pairing interaction, to be compared with the one around 10 MeV for the other interactions. Hence it is expected that the pairing vibration transition strength should be larger in the case of a pure surface force. However it is known that it is difficult to describe accurately the magnitude of these transitions, especially for absolute cross section calculations [20]: one-step or sequential two-step process, triton wave function, zero-range or finite-range DWBA have to be considered. It is therefore necessary to rely on the angular distribution, calculated from the form factor, related itself to the pairing transition density [13].

The pairing transition density is defined as:

κ ν (r, σ) = 0|c (r, σ) c (r, σ) |ν (7)
where c † (r, σ)= -2σc † (r, -σ) is its time reversed counterpart.

It allows to calculate the form factor in the zero-range DWBA approximation. The pairing transition densities of Fig. 7 show, in the case of 124 Sn, a difference, going from surface to other modes: the transition density decreases at the surface. However the difference is not dramatic and may be overruled by the experimental uncertainties. The larger strength of the 9 MeV peak in the pure surface pairing interaction is due to a larger density at the surface.

For the 136 Sn neutron-rich nucleus, the low energy spectrum displayed on Fig. 8 is dramatically changed from using surface to other interactions, on a several MeV area. A three peaks structure appears in the surface case, compared to the 2 peak structure of the other cases. The integrated strength is also larger in the surface case. Fig. 9 and 10 show the corresponding transition densities. They exhibit very different shapes, comparing results with the pure surface pairing interaction and the mixed pairing interaction. Hence 136 Sn is a good test case to probe the pairing interaction through pairing vibrations. For instance in the case of the most intense peak, the central part is dominant in the transition density for the mixed case, whereas the surface part of the transition density dominates in the pure surface interaction. Hence a measurement of the angular distributions associated with the pairing vibration strength in very neutron rich-nuclei such as 136 Sn seems more decisive to disentangle between the pairing interactions than with 136 Sn. This may be due to the larger neutron skin in 136 Sn, consisting of low density neutron-rich matter.

It has been shown in a previous work how the pairing transition densities allow to calculate the two neutron form factor in order to predict angular distributions for pairing vibrations [14]. Work along these lines should be undertaken in order to bring additional constrains on the pairing interaction.

IV. CONCLUSIONS

The impact of various pairing interactions on pairing vibrations predictions has been analysed for the first time using a HFB+QRPA approach. They should provide a good sensitivity from a pure surface interaction compared to mixed interactions, especially in the case of very neutron-rich nuclei such as 136 Sn. Moreover nuclear matter gap calculations show that the LDA in the pairing channel is valid in the surface of the nuclei for the η = 0.35 mixed surface/volume interaction, but not rigorously valid in other cases. In the case of exotic nuclei, pairing vibrations are found more sensitive to the surface/volume type of the pairing interaction, than in the case of stable nuclei. This may be due to the larger extension of the neutron density in very neutron-rich nuclei.

The same study using an isospin dependent pairing interaction will be undertaken. The hope is to come one step closer to a more global pairing interaction, using odd even mass staggering, pairing vibrations, and nuclear matter as constraints. Experimentally, the pairing transition densities can be tested through the form factor used to calculate the two neutrons transfer cross section. This implies to use a adequate reaction model. Work along these lines will be undertaken in an near future.

FIG. 2 :

 2 FIG. 2: Pairing gap versus the density for uniform matter for different pairing interactions.

Fig. 3 ηFIG. 3 :

 33 FIG. 3: Pairing gap versus the density for nuclear matter and for 124 Sn (see text)

1 )FIG. 5 :

 15 FIG.5: Unperturbed response function for 136 Sn in the two neutrons 0 + addition mode. The pure surface mode is in solid line, the η=0.65 mode is in dotted line, and the η=0.35 mode in dashed-dotted lines.

1 )FIG. 6 : 3 )FIG. 7 : 1 )FIG. 8 :

 163718 FIG.6: QRPA response function for 124 Sn in the two neutrons 0 + addition mode. The pure surface mode is in solid line, the η=0.65 mode is in dotted line, and the η=0.35 mode in dashed-dotted lines.

FIG. 9 :

 9 FIG.9: Neutron transition density in the two neutrons addition mode for 136 Sn for the first two peaks of the strength, in the case of the mixed η=0.65 interaction

FIG. 10 :

 10 FIG. 10: Neutron transition density in the two neutrons addition mode for 136 Sn for the first 3 peaks of the strength, in the case of the pure surface interaction

TABLE I :

 I Values of η and V0 of the pairing interaction.

	η	V0 (MeV fm -1 )
	0.35	-285
	0.65	-390
	1	-670
	II. VALIDITY OF THE LDA IN THE PAIRING
		CHANNEL

TABLE II :

 II Neutron quasiparticle states with E ≤ 6 MeV and occupation less than 80%. The nucleus is 124 Sn.

	η	State	E (MeV)	occ
	0.35	h11/2	1.5	0.42
		f7/2	5.8	0.01
	0.65	h11/2	1.7	0.42
		f7/2	5.7	0.01
	1	h11/2	2.2	0.42
		p3/2	5.4	0.003
		f7/2	5.5	0.02
		p1/2	5.6	0.002
		s1/2	5.7	0.002

TABLE III :

 III Same as in Table II but for 136 Sn. Unperturbed response function for 124 Sn in the two neutrons 0 + addition mode. The pure surface mode is in solid line, the η=0.65 mode is in dotted line, and the η=0.35 mode in dashed-dotted lines

		14
		12
	)
	-1
	* ) (MeV	8 10
	S(E	6
		4
		2
		0 10 11 12 13 14 15 16 17 18 19 20 E * (MeV)	η 0.35	State f7/2	E (MeV) 0.8	occ 0.45
				p3/2	1.9	0.01
				p1/2	2.4	0.006
	FIG. 4:		f5/2	2.9	0.01
				s1/2	3.3	0.0005
				d5/2	4.0	0.0002
				d3/2	4.1	0.0005
				g9/2	5.6	0.0001
				g7/2	5.6	0.0001
			0.65	f7/2	0.9	0.43
				p3/2	1.9	0.02
				p1/2	2.4	0.008
				f5/2	2.9	0.02
				s1/2	3.2	0.0004
				d5/2	3.9	0.0001
				d3/2	3.9	0.0003
				g9/2	5.4	0.0001
				g7/2	5.5	0.00003
			1	f7/2	1.6	0.32
				p3/2	1.7	0.02
				p1/2	1.9	0.01
				s1/2	1.9	0.0004
				d5/2	2.6	0.0003
				d3/2	2.6	0.0002
				f5/2	3.0	0.01
				g9/2	4.1	0.0002
				g7/2	4.1	0.0001
			case of the pure surface pairing, the spectrum is changed.
			The energies are shifted to lower values, and the overall
			strength is increased. The lower energy shift can be un-
			derstood by more single QP states located at low energy.
			This can be explained by a lower pairing gap and a dif-
			ferent energy spectrum found in the HFB self-consistent

Acknowledgements The authors thank D. Beaumel, M. Matsuo and A. Vitturi for fruitful discussions


Delphine Bosch
email: bosch@gm.univ-montp2.fr

Janne Blichert-Toft

Frédéric Moynier

Bruce K Nelson

Philippe Telouk

Pierre-Yves Y Gillot

Francis Albarède

Pb, Hf and Nd isotope compositions of the two Reunion volcanoes (Indian Ocean): A tale of two small-scale mantle "blobs"?

Keywords: Réunion island, plume, Pb, Hf, Nd isotopes, blob, Piton de la Fournaise: Piton des Neiges, Indian ocean, C component, EM-1 component

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Studies of geochemical and isotopic variations through time of plume-lavas erupted from a single volcano or from neighbouring volcanoes are expected to shed light on the composition and evolution of mantle regions involved in a hotspot environment and to document the processes taking place during magma ascent.

Réunion Island is one of the best known hot-spot volcanoes and has often been compared to the Hawaiian and Kerguelen archipelagos (e.g., Blichert-Toft and Albarède, 1999;[START_REF] Eiler | Oxygen isotope constraints on the sources of Hawaiian volcanism[END_REF][START_REF] Eisele | The role of sediment recycling in EM-1 inferred from Os, Pb, Hf, Nd, Sr isotope and trace element systematics of the Pitcairn hotspot[END_REF][START_REF] Hofmann | Source characteristics derived from very incompatible trace elements in Mauna Loa and Mauna Kea basalts, Hawaii Scientific Drilling Project[END_REF][START_REF] Ingle | Hf isotope constraints on mantle sources and shallow-level contaminants during Kerguelen hot spot activity since ~120 Ma[END_REF][START_REF] Lassiter | Isotopic evolution of the Mauna Loa volcano: Results from the initial phase of the Hawaii Scientific Drilling Project[END_REF][START_REF] Neal | Mantle sources and the highly variable role of continental lithosphere in basalt petrogenesis of the Kerguelen plateau and Broken Ridge LIP: Results from ODP Leg 183[END_REF]. This island is considered as one of the best places worldwide to document the petrological history of relatively simple volcanoes during the different cycles of volcanic eruptions (e.g. [START_REF] Fisk | Geochemical and experimental study of the genesis of Réunion island, Indian Ocean[END_REF][START_REF] Graham | A 360,000 year helium isotope record from Piton de la Fournaise, Réunion Island[END_REF][START_REF] Mcdougall | Strontium isotopic composition and potassium-rubidium ratios in some rocks from Réunion and Rodriguez, Indian Ocean[END_REF]. A detailed comparison of the geochemistry between the two volcanoes of the island, Piton des Neiges (PDN) and Piton de la Fournaise (PDF) (Fig. 1), however, is not available. Published Sr, Nd, Pb and He isotopic data from the Réunion volcanoes define a narrow range of variation [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF][START_REF] Albarède | Modelling the recent geochemical evolution of the Piton de la Fournaise volcano, Réunion Island[END_REF][START_REF] Deniel | Le magmatisme du Piton des Neiges. Monographie[END_REF][START_REF] Fisk | Geochemical and experimental study of the genesis of Réunion island, Indian Ocean[END_REF][START_REF] Graham | A 360,000 year helium isotope record from Piton de la Fournaise, Réunion Island[END_REF][START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF][START_REF] Oversby | Genetic relations among volcanic rocks of Réunion: chemical and lead isotopic evidences[END_REF][START_REF] Staudacher | Noble gas systematics of Réunion island, Indian Ocean[END_REF]. Recently, subtle patterns of isotopic variability within this range have been identified in Nd and Pb isotope studies from Piton de la Fournaise [START_REF] Hanyu | Noble gas study of the Réunion hotspot: evidence for distinct less-degassed mantle sources[END_REF][START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF][START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF].

In this study we report Pb, Hf, Nd and Sr isotope analyses of a large sample set of lavas from Réunion Island. This constitutes the first extensive isotopic study of lavas from the two volcanoes at the scale of the island, and in particular the first report of Hf isotope data. In order to detect changes of the isotopic composition through time between the two volcanoes but also between the different sequences of each volcano (i.e., from 2 Ma to present) samples were selected to cover most of the sub-aerial history of the volcanoes. The new data document the spatial and temporal evolution of Réunion plume geochemistry and allow for a comparison with other hotspot environments.

Geological setting

Réunion Island lies 750 km east of Madagascar and is considered to be related to the hot-spot track that trends northeast through Mauritius Island, along the Mascarene Plateau and ends in the Deccan Traps of western India (e.g., [START_REF] Duncan | Réunion hot spot activity through Tertiary time: initial results of the Ocean Drilling Program, Leg 115[END_REF][START_REF] Morgan | Hotspots tracks and the opening of the Atlantic and Indian oceans[END_REF]. [START_REF] Burke | The African plate[END_REF] challenged this model and proposed that Réunion Island results from the activity of a younger hot-spot unrelated to the track of the Deccan hot-spot. In 1997, Albarède and co-workers [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF] proposed that each volcano results from impingement at the base of the lithosphere of a blob representing a solitary wave of hot-spot material with a radius of approximately 100 km.

The island consists of two volcanoes, 30 km apart, that were active simultaneously over at least 500 000 years, from ca. 530 to 30 ka [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF]. They show strong petrological similarities with a predominance of basalts and hawaiites [START_REF] Albarède | Modelling the recent geochemical evolution of the Piton de la Fournaise volcano, Réunion Island[END_REF][START_REF] Fisk | Geochemical and experimental study of the genesis of Réunion island, Indian Ocean[END_REF][START_REF] Ludden | Magmatic evolution of the basaltic shield volcanoes of the Réunion island[END_REF][START_REF] Nativel | Etude pétrographique et géochimique des volcans de la Réunion[END_REF][START_REF] Sobolev | Petrology of long-lived mantle plume magmatism: Hawaii, Pacific, and Réunion island, Indian Ocean[END_REF].

Piton des Neiges is the oldest volcano of the island with the earliest eruptive cycle dated at ca. 2 Ma and the youngest at around 30 ka [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF][START_REF] Mcdougall | The geochronology and evolution of the young volcanic island of Réunion, Indian Ocean[END_REF][START_REF] Rocher | Structural evolution of Piton des Neiges volcano[END_REF]. This volcano is now extinct and occupies the northwestern two thirds of the island. The oldest series (previously called the Oceanite Series) belongs to the area known as La Montagne [START_REF] Mcdougall | The geochronology and evolution of the young volcanic island of Réunion, Indian Ocean[END_REF] and consists of transitional to silica-saturated alkaline lavas. This series, erupted during the shield-building volcanic phase between 2.1 Ma and c. 440 ka, can be divided into three major groups: a first group mainly comprised of olivine basalts with ages between 2.1 and 1.9 Ma, an intermediate phase with eruption of basalts and, a last younger phase, corresponding to the end of shield-building between 580 and 440 ka [START_REF] Deniel | Le magmatisme du Piton des Neiges. Monographie[END_REF]. The late-stage capping lavas, corresponding to the construction of a strato-volcano within the eroded shield volcano, are made of ignimbrites and alkaline under-saturated lavas with ages of 180 to 68 ka [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF][START_REF] Rocher | Structural evolution of Piton des Neiges volcano[END_REF]. The younger Piton des Neiges basalts are silica over-saturated alkaline lavas and are contemporaneous with the oldest lavas of Piton de la Fournaise (from 570 to 300 ka) [START_REF] Rançon | The Grand Brûlé exploration drilling: new data on the deep framework of the Piton de la Fournaise volcano. Part I: Lithostratigraphic unit and volcanostructural interpretation[END_REF]. Fifteen lavas from Piton des Neiges, selected to cover most of its lifetime (from 2150 to 374 ka), were selected for this study.

The large Piton de la Fournaise volcano erupted successive series of lavas that evolved from mildly alkalic very early on to tholeiitic for most of its volcanic history [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF].

On the basis of K-Ar data [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF], three distinct prehistoric eruptive cycles have been recognized and dated at 527-290 ka (Rivière des Remparts, RP), 180-105 ka (Nez de Boeuf, NB) and 70-40 ka (Langevin, LGV). We also studied modern lava flows erupted between 1700 and 1956 (Figure 1). Picrites, basalts, plagioclase-rich basalts and differentiated samples have been recognized within the four sequences. The oldest series (RP) consists of basalts, picrites, differentiated lavas and plagioclase-cumulates with scarce gabbroic fragments [START_REF] Bachélery | Evolution volcano-structurale du Piton de la Fournaise, depuis 0,53 Ma[END_REF]. A drill hole in the Grand Brulé (Enclos Fouqué) encountered gabbros, wehrlites and dunites with traces of hydrothermal alteration [START_REF] Augé | The Grand Brûlé exploration drilling: new data on the deep framework of the Piton de la Fournaise volcano. Part 3: Mineral chemistry of the cumulate rocks[END_REF][START_REF] Lerebour | The Grand Brûlé exploration drilling: new data on the deep framework of the Piton de la Fournaise volcano. Part 2: Secondary minerals[END_REF]. The hiatus between the RP series and the three other sequences (no samples dated between 290 ka and 219 ka) coincides with a major discontinuity between the older series and the most recent shield, possibly reflecting low volcanic activity, intense erosion or formation of a new caldera [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF]. The NB series consists of olivine-bearing basalts and is thought to represent the volcanic infilling of an ancient canyon. This sequence is made of lavas erupted after the major discontinuity observed between the RP series (520-290 ka) and the more recent shield (built at approximately 150-130 ka). The NB series incorporates the occurrence of major extended discontinuities between the various lava flows [START_REF] Bachélery | Evolution volcano-structurale du Piton de la Fournaise, depuis 0,53 Ma[END_REF]. The LGV sequence is made up of basalts and picrites with few aphyric samples and one ankaramite (LGV14). The historic lavas (M) correspond to samples taken from the east flank of the Dolomieu volcano. They consist of basalts and picrites related to major volcanic pulses disrupting pre-existing olivine-rich cumulates [START_REF] Albarède | Modelling the recent geochemical evolution of the Piton de la Fournaise volcano, Réunion Island[END_REF]. Eighty-seven lavas that encompass the total range of lavas erupted during the entire lifetime of Piton de la Fournaise have been analysed.

Analytical methods

Powdered samples were weighed to obtain approximately 100 to 200 ng of Sr, Nd, Pb and Hf. A leaching step with 6N HCl during 30 minutes at 65°C was done before acid digestion. Samples were dissolved during 36-48 hours on a hot plate with a mixture of HF and HNO 3 . After evaporation to dryness, 1ml of HNO 3 was added to the residue and kept at about 90°C for 12-24 h. For Pb separation, after a complete evaporation, 0.5 ml of 8N HBr was added to the sample and kept at 70°C for 2-3 h before complete evaporation. The chemical separation of Pb was done using 50l of anion exchange resin (AG1X8, 200-400 mesh) and samples were loaded and washed in 0.5N HBr. Lead was then eluted in 6N HCl. Lead blanks were less than 40pg and are considered negligible for the present analyses. Hafnium chemical purification followed the method outlined by [START_REF] Blichert-Toft | Separation of Hf and Lu for high-precision isotope analysis of rock samples by magnetic sector multiple collector ICP-MS[END_REF]. Strontium isotopes were separated using Sr Eichrom resin [START_REF] Pin | Concomitant separation of strontium and samarium-neodymium for isotopic analysis in silicate samples, based on specific extraction chromatography[END_REF]. Neodymium isotopes were separated following [START_REF] Richard | 143 Nd/ 146 Nd, a natural tracer: an application to oceanic basalts[END_REF] including a first step of REE separation (using AG50WX12 cation exchange resin) followed by two steps of Nd purification using the HDEHP technique. Strontium, Nd and Hf total procedural blanks were less than 50 pg, 15 pg and 20 pg, respectively. Strontium isotope ratios were measured at the University of Toulouse. The other Sr isotopic data, labelled as # in Table 1, are from a previous study [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF]. Neodymium isotopes were measured either by thermal ionization mass spectrometry at the University of Montpellier II or, for Nd ratios labelled as § in Table 1, by multi-collector inductively-coupled plasma mass spectrometry (MC-ICP-MS; VG Plasma 54) at the Ecole Normale Supérieure in Lyon (ENSL). Details about analytical parameters, including reproducibility, accuracy and standards, are available online in Supplemental-File-1. Lead, Hf and some Nd (samples labelled § in Table 1) isotope compositions were determined by MC-ICP-MS using the VG Plasma 54 at ENSL [START_REF] White | High-precision analysis of Pb isotope ratios by multi-collector ICP-MS[END_REF][START_REF] Blichert-Toft | Separation of Hf and Lu for high-precision isotope analysis of rock samples by magnetic sector multiple collector ICP-MS[END_REF][START_REF] Luais | Precise and accurate neodymium isotopic measurements by plasma-source mass spectrometry[END_REF][START_REF] Albarède | Precise and accurate isotopic measurements using multiple-collector ICPMS[END_REF].

Results

The Pb, Sr, Hf and Nd isotope data are given in Table 1. 206 3a) characterized by a slope of 1.51 and a correlation coefficient of 0.91 (n=13).

No obvious relationship exists between Pb isotope ratios and age of the lavas. The seven oldest lavas encompass the whole range of variation of the Piton des Neiges samples, whereas the younger lavas yield less radiogenic Pb (Fig. 3). With respect to the Northern Hemisphere Reference Line [START_REF] Hart | Heterogeneous mantle domains: signatures, genesis and mixing chronologies[END_REF] 3a). Most samples from the four Piton de la Fournaise series plot along this line and exhibit a broad correlation between the age of the series and the Pb isotope compositions, the oldest sequence being characterized by the most radiogenic Pb isotope compositions. A few samples have significantly less radiogenic Pb and define an additional trend (referred to as Trend 2' in Fig. 3a), which clearly points toward the ankaramite LGV14 from the Langevin sequence. However, these samples yield  Hf (+9.0 and +9.2) and  Nd (+4.2 and +4.3) values indistinguishable from those of the other Piton de la Fournaise samples. The Trend 2' comprises at least two historical lavas (T1602 and M1701) but also some Langevin and Nez de Boeuf samples and two of the three youngest lavas from Piton des Neiges. The Sr isotope compositions measured for the Piton de la Fournaise lavas vary from 0.704034 to 0.704227 (this study and [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF][START_REF] Albarède | Modelling the recent geochemical evolution of the Piton de la Fournaise volcano, Réunion Island[END_REF] and are within the range of values reported earlier [START_REF] Fisk | Geochemical and experimental study of the genesis of Réunion island, Indian Ocean[END_REF][START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF][START_REF] Ludden | Magmatic evolution of the basaltic shield volcanoes of the Réunion island[END_REF]. Previous Nd isotope data from Piton de la Fournaise [START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF] showed a temporal increase of the Nd isotope ratios (Fig. 2b). Although our study is based on fewer measurements such an increase of 143 Nd/ 144 Nd through time is not observed. Instead, we find a nearly constant Nd isotopic composition. Our

Remparts Nd isotope data are 1  Nd higher than those reported by [START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF], whereas historic samples are 1  Nd lower. The temporal increase of Nd isotopes at nearly constant Sr ratios is a distinctive feature of these data and is not shown by other OIB. Moreover, very recent samples from Piton de la Fournaise [START_REF] Salters | Hf isotope constraints on mantle evolution[END_REF][START_REF] Debaille | Geochemical component relationships in MORB from the Mid-Atlantic Ridge, 22-35 degrees[END_REF]I. Vlastelic, unpub. data) as well as very old submarine samples (taken from seamounts covered by thick Fe-Mn deposits; [START_REF] Fretzdorff | Geochemistry and petrology of lavas from the submarine flanks of réunon Island (western Indian Ocean): implications for magma genesis and the mantle source[END_REF] display the same narrow Nd isotope range as our Piton de la Fournaise Nd isotope variation (Fig. 2b). Taken together, these observations do not support a 143 Nd/ 144 Nd increase through time in Piton de la Fournaise lavas. The trend reported by [START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF] may be an analytical artefact from not separating Nd from the other REEs prior to MC-ICP-MS analysis.

In Pb-Pb, Hf-Nd, Hf-206 Pb/ 204 Pb and Nd-206 Pb/ 204 Pb diagrams (Fig. 3, 4 and 5), the RP lavas display more restricted and homogeneous compositions in comparison with all other sequences (including the Piton des Neiges samples). One sample containing disaggregated gabbroic fragments (RP31) yields slightly less radiogenic Nd, Hf and Pb isotope compositions compared to other RP samples (Table 1). The ankaramite (LGV 14), which is the most magnesian of all of the analyzed PDF lavas (MgO=26.67wt.%; [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF], has the lowest 208 Pb/ 204 Pb and 206 Pb/ 204 Pb isotope compositions and one of the highest 87 Sr/ 86 Sr ratios of the present data set.

Nevertheless, its  Nd and  Hf values are not distinct from the rest of the PDF lavas (Table 1).

LGV-14 has high CaO/Al 2 O 3 and Sc/Yb ratios combined with very low FeO/MgO (Albarède et al., 1997). In Figure 3a, it is located at the lower extremity of Trend 2'. Two historical samples (T1602, M1701) display significantly lower 206 Pb/ 204 Pb and 208 Pb/ 204 Pb values than samples from other PDF sequences (Fig. 3a). These two samples partly overlap the 1998-2002 historical samples analysed by [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF] (Fig. 3a).

5.Discussion

Identification of various mantle components in the Réunion plume

Plume models suggest that plumes are compositionally composite [START_REF] Campbell | Implications of mantle plume structure for the evolution of flood basalts[END_REF][START_REF] Farnetani | Petrological models of magma evolution and deep crustal structure beneath hotspots and flood basalt provinces[END_REF]. If plume sources contain a variety of recycled material (igneous and sedimentary sections of the oceanic crust, oceanic lithospheric mantle, plume heads, etc…) and, in addition, if this material is recycled at different times during Earth's history, plume material is not expected to be isotopically homogeneous.

The RP lavas may represent a dominant plume component at Piton de la Fournaise. The RP lavas represent a long time interval (520-290 ka) and are conspicuously homogeneous. As shown in studies of other hot-spots, the contribution of a plume component is often positively correlated with lava production rates [START_REF] Hanan | The dynamic evolution of the Iceland mantle plume: the lead isotope perspective[END_REF][START_REF] Schilling | Fluxes and excess temperatures of mantle plumes inferred from their interaction with migrating ridge[END_REF]. The RP volcanic sequence erupted during the shield-building stage of Piton de la Fournaise activity and corresponds to high volcanic production rates (more than 1000 m of lavas now exposed at the surface It remains to be determined whether the isotopic variability of Réunion lavas reflects a heterogeneous source, as in Hawaii and Iceland (e.g. [START_REF] Hofmann | Mantle geochemistry: The message from oceanic volcanism[END_REF], or whether shallow level entrainment of upper mantle and contamination by lithosphere has modified an originally homogeneous isotopic signal [START_REF] Eiler | Oxygen isotope constraints on the sources of Hawaiian volcanism[END_REF][START_REF] Hauri | The Nature of Entrainment in Mantle Plumes: A boundary layer model including the effects of temperature and stress-dependent rheologies and depth-dependent physical properties[END_REF]. Participation of a small amount of MORB-type mantle has been proposed for the composition of the plume-magmas supplying recent shield lavas in Mauritius [START_REF] Peng | Drillhole lavas from the northwestern Deccan Traps, and the evolution of Réunion hotpsot mantle[END_REF][START_REF] Paul | Geochemistry of Mauritius and the origin of rejuvenescent volcanism on oceanic island volcanoes[END_REF][START_REF] Nohda | Systematic Variation of Sr-, Nd-and PbIsotopes with Time in Lavas of Mauritius, Réunion Hotspot[END_REF]. Entrainment of the ambient mantle as envisioned by [START_REF] Hauri | The Nature of Entrainment in Mantle Plumes: A boundary layer model including the effects of temperature and stress-dependent rheologies and depth-dependent physical properties[END_REF] is dynamically

improbable because temperature-dependent viscosity efficiently decouples the plume conduit from surrounding mantle [START_REF] Farnetani | Numerical investigation of the mantle plume initiation model for flood basalt events[END_REF][START_REF] Olson | Joint pattern development: Effects of subcritical crack growth and mechanical crack interaction[END_REF]. Shallow contamination of upwelling hot-spot material ponding beneath the lithosphere by the upper mantle and the oceanic crust would be more likely. In the case of Hawaii, however, Blichert-Toft et al. (1999) observed a strongly curved mixing array between  Hf and 206 Pb/ 204 Pb, that clearly leaves asthenospheric and lithospheric material off the trend. [START_REF] Huang | Trace element abundances of Mauna Kea basalt from phase 2 of the Hawaii Scientific Drilling Project: Petrogenetic implications of correlations with major element content and isotopic ratios[END_REF] pointed out that more than two components are necessary to account for the isotope geochemistry of Hawaiian lavas, but this additional complexity does not alter the simple geometrical relationship that excludes the depleted mantle from the mixing array. For the Réunion volcanoes, unfortunately, the range of isotopic variation is insufficient to obtain a well-defined mixing array.

We first consider the possibility that the Nd, Hf, and Pb isotope variability of the modern Réunion plume results from minor addition of depleted mantle to the pure plume-source material.

Such depleted MORB-source mantle would have a composition close to that of the basalts from site 706-701 of ODP Leg 115 [START_REF] White | Isotope geochemistry of Leg 115 basalts and inferences on the history of the Réunion mantle plume[END_REF], which have isotopic characteristics similar to that of typical Indian Ocean basalts [START_REF] Xu | Geochemical and Nd-Pb isotopic characteristics of the Tethyan asthenosphere: implications for the origin of the Indian ocean mantle domain[END_REF][START_REF] Zhang | Evidence for a widespread Tethyan upper mantle with Indian ocean-type isotopic characteristics[END_REF]. The contribution of this depleted component, however, must be limited compared to Hawaii and some other areas of the western Indian Ocean [START_REF] White | Isotope geochemistry of Leg 115 basalts and inferences on the history of the Réunion mantle plume[END_REF]. For example, in the 143 , 2006). The first and second principal components can be identified with, respectively, the C component of [START_REF] Hanan | Lead and helium isotope evidences from oceanic basalts from a common deep source of mantle plume[END_REF] and a minor EM-1 component, the hallmark of the Dupal anomaly, which is interpreted as either recycled pelagic sediments [START_REF] Chauvel | Himu Em -the French-Polynesian Connection[END_REF] or recycled plume heads [START_REF] Gasperini | Evidence from Sardinian basalt geochemistry for recycling of plume heads into the Earth's mantle[END_REF].

Small-scale geochemical variations in relation to shallow-level processes

Previous isotopic analyses of lavas from Piton de la Fournaise have revealed that compositional small-scale isotopic variations occur at very different time-scales, from thousands of years to a few tens of years, but also during the course of a single eruption [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF]. Such small-scale compositional changes potentially may be linked to specific processes such as rate of eruption, volume of lavas, or style of the volcanic eruption.

The Rempart sequence has, on average, the most radiogenic 208 Pb/ 204 Pb and 207 Pb/ 204 Pb ratios of all the sequences studied. This series, erupted early in the formation of the old Piton de la Fournaise shield-volcano, is more than 1000 m thick [START_REF] Bachélery | Evolution volcano-structurale du Piton de la Fournaise, depuis 0,53 Ma[END_REF]. The rate of volcanic eruption of this thick lava pile is difficult to evaluate given the lack of geochronological data and the numerous gaps throughout the series [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF]. Nevertheless, the decrease of the 208 Pb/ 206 Pb ratio from the oldest Rempart series to 1330 m with a subsequent increase in the younger lavas is particularly striking (Supplemental-File-2). In principle, such a Pb isotope distribution can be explained by variation of the plume composition on a very short timescale. However, this disagrees both with the PCA model, which suggests that 98.2% of the Piton de la Fournaise Pb isotope composition is explained by mixing between two major components, and with the steady decrease of 208 Pb/ 206 Pb observed from bottom to top of the section.

Alternatively, the cause may be limited shallow-level contamination/assimilation processes either with (1) the old Paleocene oceanic crust surrounding the magma chamber or with ( 2 surrounding old oceanic crust [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF]. Such assimilation of country rocks is typically related to periods of volcanic quiescence, for example during refilling of the magma reservoir between two eruption events [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF] or, more likely in the present case, to the onset of activity at a new volcanic edifice when early magmas more efficiently assimilate country-rocks. The correlation of the decrease in 208 Pb/ 206 Pb with sample location in the volcanic pile (i.e., age of the flow) identifies a temporal decrease of crustal assimilation processes (i.e., plagioclase-rich basalts and differentiated lavas with aphyric textures comprising samples below ~1000 m). Basaltic and picritic samples at the top of this series display more homogeneous Pb isotope ratios and slightly higher 87 Sr/ 86 Sr ratios, with an average of 0.70411 (samples > 900 m) versus 0.70408 (samples < 900 m). Thus, the upper samples are interpreted as having the isotopic composition closest to that of the plume feeding Réunion volcanism. Nevertheless, it is noteworthy that throughout the Rempart lavas, erupted during approximately 237 ky, the variation of the 208 Pb/ 206 Pb ratio is lower ( 208 Pb/ 206 Pb = 0.0058) than that registered by historic lavas emitted by Piton de la Fournaise between 1998 and 2002 ( 208 Pb/ 206 Pb = 0.0085; [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF]). This highlights the generally homogeneous character of this prehistoric volcanic sequence. Moreover, the Nd and Hf isotope compositions of the Remparts basalts, excluding plagioclase-rich sample RP31, are the most homogeneous of the studied lavas from Piton de la Fournaise. A long period of little to no activity, approximately 70 ky (between 290 and 219 ka), separates the Rempart from the Nez de Boeuf series [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF]. The Nez de Boeuf sequence is over 200 m thick and spans 85 ky [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF] with extensive discontinuities.

Field observations show that the Nez de Boeuf lavas cross-cut pre-existing Rempart dykes. The Pb isotope distribution of these lavas is complex. It displays a relatively large range of Pb isotope ratios compared to the three other sequences from Piton de la Fournaise and partly overlaps compositions of both older and younger lava sequences. This feature can be related to the rapid variation of magma composition over a short time-scale but appears difficult to reconcile with the conclusions reached for the thicker and longer Rempart series. It can be consistent with a model of lava flows reacting with older solidified magmas present as a dispersed dyke system located at shallow level, though this does not explain the less radiogenic samples of the Nez de Boeuf series (i.e. samples with Pb isotope ratios out of the range of older Remparts series). Thus, the observed Pb isotope shift is more likely related to a major modification of the volcanic system between the Remparts and Nez de Boeuf eruptions, an assumption which is strengthened by the hiatus separating the two series.

The Langevin sequence, which comprises a 485 m thick sequence of lavas erupted over 30 ky, was erupted after a quiescent period of almost 40 ky [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF]. Compared to the Rempart sequence, it shows, on average, a clear tholeiitic affinity and significantly higher 208 Pb/ 206 Pb and lower 206 Pb/ 204 Pb and 208 Pb/ 204 Pb ratios, as well as more scattered  Hf and  Nd (Table 1). No clear decrease or increase of Pb isotope ratios versus height is evident in this sequence. In Fig. 3, all samples plot along Trend 2 except for a few samples (including the ankaramite) that define Trend 2' together with two historic lavas. On a very short time-scale ( 4 years), [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF] observed an inverse correlation between volcanic eruption rate and the 206 Pb/ 204 Pb ratio for historic lavas. On a thousand year time-scale, we observe a general decrease in the 206 Pb/ 204 Pb and 208 Pb/ 204 Pb ratios broadly positively correlated with the proportions of picritic samples erupted from Piton de la Fournaise (i.e., 1/30 for Rempart, 4/23 for Langevin and 5/21 for the modern lavas). This general decrease in Pb isotope ratios is also roughly correlated with a decrease in the alkalinity index. In a hot-spot environment, most of the primary melts generated at depth tend towards a picritic composition [START_REF] Farnetani | Numerical investigation of the mantle plume initiation model for flood basalt events[END_REF][START_REF] Mckenzie | The volume and composition of melt generated by extension of the lithosphere[END_REF]. In the case of the Réunion volcanoes, only a few picrites and MgO-rich samples have been found. Nevertheless, the high seismic velocities observed at depths greater than 15 km [START_REF] Driad | Crustal and upper mantle structure under the Réunion hotspot from seismic wide angle data[END_REF] could be related to the presence of abundant clinopyroxene-rich cumulates. A model including a deeper origin, rapid ascent to the surface and contamination with hydrothermal products lying on top of the oceanic basement has been proposed to explain the Pb isotope composition exhibited by picrites erupted from Piton de la Fournaise volcano in January 2002 [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF]. In this model, Trend 2' is regarded as a mixing line between deeper cumulates, such as the ankaramite LGV14, and more "common" steady-state basalts. Following Vlastelic and co-authors ( 2005) such a "mixing process" is thought to occur preferentially during inter-eruption periods but is difficult to evaluate due to the lack of precise volume and rate data for prehistoric lavas. The increase of picritic basalts suggests a possible relationship with emptying of the deep magma reservoir.

The ankaramite LGV14

Two main modes of fractionation appear to be present in Piton de la Fournaise basalts: (1) olivine and Nd isotope compositions. This suggests that this particular ankaramite either formed by remelting clinopyroxene-rich cumulates hydrothermally altered by seawater or at least interacted with such rocks. An observation consistent with this interpretation is the evidence for deep hydrothermal interaction between seawater derived fluids and plutonic materials through magmatic conduits under the southwestern flank of the volcano [START_REF] Augé | The Grand Brûlé exploration drilling: new data on the deep framework of the Piton de la Fournaise volcano. Part 3: Mineral chemistry of the cumulate rocks[END_REF][START_REF] Lénat | Geoelectrical structure of the central zone of Piton de la Fournaise volcano (Réunion)[END_REF][START_REF] Lerebour | The Grand Brûlé exploration drilling: new data on the deep framework of the Piton de la Fournaise volcano. Part 2: Secondary minerals[END_REF][START_REF] Rançon | The Grand Brûlé exploration drilling: new data on the deep framework of the Piton de la Fournaise volcano. Part I: Lithostratigraphic unit and volcanostructural interpretation[END_REF]. A similar model of deeper origin and rapid ascent to the surface [START_REF] Albarède | Modelling the recent geochemical evolution of the Piton de la Fournaise volcano, Réunion Island[END_REF] has been proposed to explain the Pb isotope composition of picrites erupted from Piton de la Fournaise in January 2002 [START_REF] Vlastelic | Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: evidence for variable crustal contamination[END_REF]. An alternative explanation is to consider the ankaramite as a foreign 'nugget', such as a fragment of oceanic lithosphere, within the mantle plume. However, this would require that anhydrous melting of fertile or refractory lherzolite can potentially yield magmas with particularly high CaO/Al 2 O 3 .

Does a model of two small-scale mantle blobs explain the Réunion volcanoes?

One interesting pattern in our new data set is that the oldest lavas from Piton des Neiges (7 samples studied) have the same range of Pb isotope variations as nearly all Piton de la Fournaise lavas (80 samples). As previously mentioned, Piton des Neiges samples display significantly less radiogenic 207 Pb/ 204 Pb and 208 Pb/ 204 Pb isotope ratios at constant 206 Pb/ 204 Pb than samples from Piton de la Fournaise, with differences significantly larger than analytical errors (Table 1; Fig. 2a). Moreover, another striking feature is the greater scatter observed for Piton des Neiges samples compared to those from Piton de la Fournaise for all isotopic systems (Fig. 3,4,5).

The lavas erupted during the shield-building phase of each volcano are very voluminous and bear distinct isotopic signatures (i.e., La Montagne sequence and Les Remparts series). Early shallowlevel contamination/assimilation processes alone therefore cannot account for the range of isotopic differences observed between Piton de la Fournaise and Piton des Neiges lavas.

Likewise, a progressive temporal change of the plume composition between Piton des Neiges and Piton de la Fournaise volcanic activity is unlikely, since both volcanoes contemporaneously erupted lavas with markedly different isotopic ratios (i.e., younger lavas from the PDN and older lavas from the PDF). The differences in isotopic composition between lavas from Piton de la Fournaise and Piton des Neiges volcanoes therefore suggest that each volcano tapped different parts of a single large heterogeneous plume comprised of two distinct small-scale blobs with different geochemical composition. The shift in location of the volcanic center between Piton des Neiges and Piton de la Fournaise activity supports the existence of two small-scale blobs. The contemporaneous eruption of lavas from both volcanoes took place during a short time span (approximately 400 ky) and indicates that the magma production from the blob feeding Piton des

Neiges waned as the next blob started to feed Piton de la Fournaise. The distance between the two volcanic edifices (around 35 km) suggests that, if magma rises vertically, each blob had a relatively small radius (< 10 km), which is consistent with estimates of c. 14 km for waves in conduits fed by a pulsating or blob-like plume [START_REF] Whitehead | Magmas waves and diapiric dynamics[END_REF]). In the case of the Hawaiian plume, where the shield volcanoes are distributed into two curvilinear parallel trends with distinct geochemical characteristics, the radius of the axial part of the plume has been estimated at around 20 km [START_REF] Bryce | Geochemical structure of the Hawaiian plume: Sr, Nd, and Os isotopes in the 2.8 km HSDP-2 section of Mauna Kea volcano[END_REF]. The occurrence of these two distinct "trends" is interpreted in different ways. It could be related to melting of a plume concentrically zoned with a radial systematic variation of isotopic ratios [START_REF] Hauri | The Nature of Entrainment in Mantle Plumes: A boundary layer model including the effects of temperature and stress-dependent rheologies and depth-dependent physical properties[END_REF][START_REF] Lassiter | Isotopic evolution of the Mauna Loa volcano: Results from the initial phase of the Hawaii Scientific Drilling Project[END_REF]DePaolo et al., 2001) or it could be ascribed to a left-right asymmetry of the plume composition [START_REF] Abouchami | Lead isotopes reveal bilateral asymmetry and vertical continuity in the Hawaiian mantle plume[END_REF]. The position where the volcano trend passes over the plume (more or less closer to the plume axis) partly controls the systematic variation of the geochemical composition of the Hawaiian lavas. But a "blob" of anomalous plume material along the plume axis has to be present, in particular to take into account Mauna Loa isotopic compositions [START_REF] Bryce | Geochemical structure of the Hawaiian plume: Sr, Nd, and Os isotopes in the 2.8 km HSDP-2 section of Mauna Kea volcano[END_REF].

In the case of Réunion island, the successive ascent of such small-scale blobs (as small diapirs) could explain the discrepancy between the isotopic composition of the two volcanoes as well as the distance between the centers of the two volcanoes.

In the model presented for Réunion volcanoes in Fig. 9, we assume the existence of two distinct waves from a single pulsating plume under the Réunion lithosphere. The oldest Piton des Neiges volcano corresponds to the impingement of the first blob. The oldest lavas from this volcano were erupted during the more alkalic phase of shield-building when the lithosphere was disrupted and the plume magmas percolated through it. In some instances, parts of this lithosphere were assimilated, which is consistent with the presence of numerous dunitic and, more rarely, wehrlitic xenoliths in Piton des Neiges primitive shield basalts [START_REF] Upton | Peridotitic and gabbroic rocks associated with the shieldforming lavas of Réunion[END_REF]. Magmas erupted during disruption of the lithosphere are contaminated by shallow-level interactions, but the high volume of melt tends to reduce the influence of the contaminant end-member (Fig. 9a).

In a second stage (Fig. 9b), the eruption rate is still high, but the lithosphere has been disrupted and depleted by the lavas erupted during shield growth. The isotopic characteristics of lavas emitted during this phase (intermediate series) are the most homogeneous of those erupted by Piton des Neiges. Later melting of the remaining blob material, associated with lower eruption rates, tends to shift the isotopic composition of the erupted lavas (younger series) away from the C-component composition and towards the depleted mantle component. Although lavas from the Piton de la Fournaise volcano have isotopic compositions that fit this scheme, the oldest lavas (RP) are surprisingly the most homogeneous and the closest to the C-component composition.

This can be explained by postulating that each blob is slightly different in composition (see principal component analysis) with a more important contribution of a depleted mantle component in the blob feeding the Piton des Neiges volcano. If the second blob impinged on the lithosphere very close to the first blob, this would allow the earliest magmas to percolate through an already depleted lithospheric column (Fig. 9c). This is supported by gravity anomalies suggesting a higher anomaly underlying Piton des Neiges than Piton de la Fournaise (Malengreau et al., 1999). At last, during the final stage (Fig. 9d), Piton des Neiges is extinct and Piton de la Fournaise erupted lavas with geochemical compositions shifted towards a more depleted mantle component.

Is the Réunion plume really more homogeneous than other OIB provinces?

As discussed by [START_REF] Fisk | Geochemical and experimental study of the genesis of Réunion island, Indian Ocean[END_REF] and [START_REF] Graham | A 360,000 year helium isotope record from Piton de la Fournaise, Réunion Island[END_REF], the Réunion island lavas show a restricted range of isotope compositions and can be considered as relatively homogeneous with respect to most other hot-spots (Fig. 6, 7, 8 and Supplemental-File-3). In each isotopic space, the restricted field defined by Réunion lavas is distinct and does not overlap the large domains defined by both the Kerguelen and Hawaiian hot-spots, for example. This is relevant since in both cases the plume material and the processes that affected it are distinct. Indeed, lavas from Kerguelen result from mixing between the plume-source material and several other components, in particular old continental crust [START_REF] Mahoney | Geochemical Characteristics of Lavas from Broken Ridge, the Naturaliste Plateau and Southernmost Kerguelen Plateau -Cretaceous Plateau Volcanism in the Southeast Indian-Ocean[END_REF][START_REF] Weis | Kerguelen Archipelago Revisited -Geochemical and Isotopic Study of the Southeast Province Lavas[END_REF] and lithospheric mantle modified by complex, multi-stage metasomatic processes [START_REF] Mattielli | Evolution of heterogeneous lithospheric mantle in a plume environment beneath the Kerguelen archipelago[END_REF]. The geochemical characteristics of this hot-spot clearly identify temporal changes in the various mantle sources [START_REF] Ingle | Hf isotope constraints on mantle sources and shallow-level contaminants during Kerguelen hot spot activity since ~120 Ma[END_REF]. These mixing processes and interactions, as well as the longevity of this plume, account for its isotopic heterogeneity and its affinity with EM-1 and EM-2 components (Fig. 6, 7, 8 and Supplemental-File-3). The Hawaiian volcanoes have isotopic features clearly distinct from those of the Kerguelen lavas. Previous studies have demonstrated that the Hawaiian plume material mainly consists of three distinct components -Loihi, Koolau and Kea - [START_REF] Abouchami | High precision lead isotope systematics of lavas from the Hawaiian Scientific Drilling Project[END_REF]Blichert-Toft and Albarède, 1999;[START_REF] Eiler | Oxygen isotope constraints on the sources of Hawaiian volcanism[END_REF][START_REF] Hofmann | Source characteristics derived from very incompatible trace elements in Mauna Loa and Mauna Kea basalts, Hawaii Scientific Drilling Project[END_REF][START_REF] Lassiter | Isotopic evolution of the Mauna Loa volcano: Results from the initial phase of the Hawaii Scientific Drilling Project[END_REF]. In the various isotopic spaces, Hawaiian lavas Réunion lavas plot inside a mixing triangle of the DMM, EM-1 and HIMU components and almost perfectly in the center of the OIB field (Fig. 6, 7, 8 and Supplemental-File-3). The average composition of the modern Réunion plume, as defined in this study, has isotopic characteristics corresponding to the average of ocean island basalts worldwide. The modern Réunion plume presents strong isotopic similarities with the C-component, which is thought to be an end-member for MORB and plume basalt arrays and an internal component for oceanic island basalts [START_REF] Hanan | Lead and helium isotope evidences from oceanic basalts from a common deep source of mantle plume[END_REF][START_REF] Hanan | Depleted Iceland mantle plume geochemical signature: artifact of multi-component mixing?[END_REF]. Compared to lavas produced by the Kerguelen and Hawaiian hot-spot activities, the modern Réunion plume generates lavas characterized by a more homogeneous and 'unmodified' plume composition. The original isotopic signal of the plume has not been significantly altered during shallow level processes.

Conclusion

The average Pb, Hf, Nd and Sr isotope compositions determined for the modern Réunion plume support a relatively "pure" plume mantle-source composition similar to the C-component.

A small amount of depleted mantle component, closely resembling that currently measured in MORB from the Central Indian Ocean Ridge, also has been identified using principal component analysis. The imprint of this component appears to be stronger in the Piton des Neiges volcano than in the Piton de la Fournaise volcano. Through time, and independently for each of the two volcanoes, the isotopic composition of the last-stage lavas becomes more depleted than those of the shield-building stage. This coincides with a lower eruption rate and can be related either to increased shallow-level contamination processes or, more probably, to melting of heterogeneities present in the blob feeding the volcano. The isotopic features and available geophysical data from the two Réunion island volcanoes are best explained by a two-blob model ascending through an ambient mantle and generating one volcano each. On a global scale, the modern Réunion plume appears more homogeneous than other hot-spots and, in particular, more homogeneous than the Hawaii or Kerguelen hot-spots. For the Piton de la Fournaise lavas, the altitude is broadly correlated with the ages of the lava flows: with increasing height, Les Remparts (527-290 ka), Nez de Boeuf (180-105 ka), Langevin (70-40 ka). K-Ar ages are from [START_REF] Gillot | Evolution volcano-structurale du Piton de la Fournaise, depuis 0.53 Ma[END_REF]. For the historical lava flows, the elevation at which the samples were taken are not reported and samples have been arbitrarily positioned up to the Langevin series in the order of dates of eruption. Symbols are the same as those in Fig. 1.

Figure Captions

(b)  Nd and  Hf variations. Neodymium isotope results from previous studies are reported for comparison: grey diamond: Les Remparts [START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF], grey square: Langevin [START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF], grey small dot: historical lavas [START_REF] Luais | Temporal changes in Nd isotopic composition of Piton de la Fournaise magmatism (Réunion Island, Indian Ocean)[END_REF], grey triangle: historical lavas (Vlastelic, unpublished data), cross: lavas dredged from the submarine flanks of Réunion Island [START_REF] Fretzdorff | Geochemistry and petrology of lavas from the submarine flanks of réunon Island (western Indian Ocean): implications for magma genesis and the mantle source[END_REF]. The grey domain encompasses the totality of our Piton de la Fournaise Nd isotope data and includes data from both Vlastelic (unpublished data) and [START_REF] Fretzdorff | Geochemistry and petrology of lavas from the submarine flanks of réunon Island (western Indian Ocean): implications for magma genesis and the mantle source[END_REF]. DMM, HIMU, EM-1, EM-2 and C are from [START_REF] Zindler | Chemical geodynamics[END_REF], [START_REF] Salters | Hf isotope constraints on mantle evolution[END_REF] and [START_REF] Hanan | Depleted Iceland mantle plume geochemical signature: artifact of multi-component mixing?[END_REF]. Symbols as in Fig. 4. LGV 14

206 Pb/ 204 Pb

FIGURE 2a

LGV

  Nd/ 144 Nd versus206 Pb/ 204 Pb plot (Fig.6), Réunion plots much farther from the depleted mantle component than does Hawaii (also for Nd-Hf, Pb/Pb and Hf-206 Pb/ 204 Pb plots; Fig.7, 8 and Supplemental-File-3).An alternative interpretation is heterogeneity within the plume itself. In this case, Trend 2 of Piton de la Fournaise in Pb isotope space requires at least two components, while one additional component is needed to account for Trend 2' toward the ankaramite and the other samples plotting below the main trend. A principal component analysis (PCA) carried out on the Piton de la Fournaise Pb isotope data indicates that the first principal component (a mixing line of two geochemical components) accounts for 98.2 percent of the total variance, while the second principal component (a third geochemical component) accounts for 1.4 percent, with the remaining variance essentially representing noise. Adding Sr, Nd, and Hf isotope data to the PCA increases the apparent number of components, but it has been argued that this is an artefact of the non-linear character of mixing involving the isotope compositions of different elements[START_REF] Blichert-Toft | Geochemical segmentation of the Mid-Atlantic Ridge north of Iceland and ridge-hot spot interaction in the North Atlantic[END_REF][START_REF] Debaille | Geochemical component relationships in MORB from the Mid-Atlantic Ridge, 22-35 degrees[END_REF]

  The third component accounts for the ankaramite Trend 2' and points towards a depleted mantle component. Adding Piton des Neiges to the PCA analysis drastically changes the situation: the first principal component now accounts for 90.8 percent of the total variance and the second principal component for 8.8 percent. The three-geochemical component representation still holds but Piton des Neiges lavas clearly have a stronger imprint of a depleted mantle component and a weaker Dupal flavor.

  ) the plumbing network constituted by dykes from previous flows and cross-cut by new magmas. The top and bottom of the Rempart series are made of different types of lavas. Plagioclase-rich basalts contaminated by disaggregated gabbroic fragments have been recovered exclusively in this sequence and are located mainly at the bottom part between 700 and 900 m. These samples (RP31, RP25...) are characterized by higher 208 Pb/ 206 Pb and lower 206 Pb/ 204 Pb and 208 Pb/ 204 Pb ratios than the upper Rempart basalts and have been linked to contamination processes by the

  define a large domain characterized by high  Hf and  Nd coupled with low Pb isotope ratios with few samples trending towards the EM-1 end-member. As observed for the modern Réunionplume lavas, the most enriched lavas have been produced during the shield-building phase, whereas the most depleted lavas correspond to the late stage of volcanic activity. In Pb isotope space, Réunion lavas are characterized by a significantly 206 Pb-richer plume component compared to Kerguelen and Hawaiian lavas. The Hf and Nd isotope signatures of Réunion lavas are almost mid-way between those defining the Kerguelen field on the one hand, and the Hawaiian field on the other hand.

Figure 1 :

 1 Figure 1: Map of Réunion island. The island is composed of two large volcanic edifices, Piton des Neiges (PDN), now extinct and located in the northwest, and Piton de la Fournaise (PDF), currently active and situated in the southeast.

Figure 2 :

 2 Figure 2: (a) -206 Pb/ 204 Pb and 208 Pb/ 204 Pb versus age (in ka) for the Piton des Neiges volcano and versus elevation (m) for the Piton de la Fournaise volcano. The symbols encompass the error margins. For Piton des Neiges, the age range for the various lava series are the following: La Montagne (2 Ma), Intermediate series (from 2 Ma to 600 ka), Younger lavas (less than 577 ka).

Figure 3 :

 3 Figure 3: (a) 206 Pb/ 204 Pb versus 207 Pb/ 204 Pb and (b) 206 Pb/ 204 Pb versus 208 Pb/ 204 Pb for Piton de la Fournaise and Piton des Neiges lava flows. Trend 1 corresponds to an apparent mixing line for samples from Piton des Neiges. Trend 2 and 2' are apparent pseudo-binary mixing lines defined by samples from Piton de la Fournaise. See text for discussion.

Figure 4 :

 4 Figure 4:  Hf versus  Nd for Réunion lavas. Piton de la Fournaise lavas display a very restricted field (plain line) compared to the Piton des Neiges domain (dashed-line). The grey domain is considered as closest to the Réunion-plume composition as defined in this study.

Figure 5 :

 5 Figure 5: (a) 206 Pb/ 204 Pb versus  Hf and (b) 206 Pb/ 204 Pb versus  Nd for Piton des Neiges and Piton de la Fournaise lavas. Symbols as in Fig. 4.

Figure 6 :

 6 Figure 6: 206 Pb/ 204 Pb versus 143 Nd/ 144 Nd for Réunion lavas. The OIB, Hawaii islands and Kerguelen fields are from the GEOROC database [http://georoca.mpch-mainz.gwdg.de/georoc/].

  Figure2aClick here to download Figure: Figure2a-Bosch.pdf

  Figure3Click here to download Figure: Figure3-Bosch.pdf

  Figure9Click here to download Figure: Figure9-Bosch.pdf

  Pb/ 204 Pb, 208 Pb/ 204 Pb, 143 Nd/ 144 Nd and 176 Hf/ 177 Hf isotope ratios are reported as a function of age (Piton des Neiges lavas) or elevation Pb/ 204 Pb. In the 208 Pb/ 204 Pb vs 206 Pb/ 204 Pb diagram, the remaining samples define an alignment (Trend 1 in Fig.

	(Piton de la Fournaise samples) (Figs. 2a and 2b). [ 208 Pb/ 206 Pb and 207 Pb/ 206 Pb data are available
	in supplemental-file-2]
	Piton des Neiges lavas. Taken as a whole, the fifteen analyzed Piton des Neiges lavas yield a
	limited yet significant range of Pb isotope ratios with 206 Pb/ 204 Pb = 18.723-19.075, 207 Pb/ 204 Pb =
	15.541-15.589 and 208 Pb/ 204 Pb = 38.630-38.962. Omitting two samples (PN69T1, PN60A) with
	the least radiogenic Pb isotope ratios significantly reduces the scatter observed for 206 Pb/ 204 Pb
	(18.918-19.075) and 208 Pb/ 204 Pb (38.722-38.962). The range remains essentially unchanged for

207

  , all samples yield high 208 Pb/ 204 Pb and 207 Pb/ 204 Pb ratios for a given

	206 Pb/ 204 Pb = 18.838-19.012, 207 Pb/ 204 Pb = 15.574-15.607, and 208 Pb/ 204 Pb = 38.866-39.075,
	respectively. In the 208 Pb/ 204 Pb versus 206 Pb/ 204 Pb diagram, these samples define a correlation-line
	(referred to as Trend 2) nearly parallel to that defined by the Piton des Neiges samples
	(slope=0.98, r 2 =0.96, n=80) (Fig.
	206 Pb/ 204 Pb value, which is typical of OIB and MORB from the Indian Ocean (see Indian ocean
	MORB	and	OIB	compilation	in	GEOROC	database	[http://georoc.mpch-
	mainz.gwdg.de/georoc/]). The five samples from the intermediate series display the most
	homogeneous Pb isotope compositions of all analyzed samples. The  Hf and  Nd values of all
	Piton des Neiges lavas display a significant range of variation from +8.6 to +11.0 and +3.1 to
	+5.0, respectively. Lavas from the younger series have the most scattered and the highest  Hf
	(from +9.4 to +11.0).
	Piton de la Fournaise samples. The 87 lavas analysed from the Piton de la Fournaise volcano
	display a limited but complex Pb isotope distribution within the different isotope spaces. As for
	the Piton des Neiges samples, lavas erupted by Piton de la Fournaise display similarities with
	Indian Ocean OIBs and MORBs. Hafnium and Nd isotopes define a small domain of variation
	with  Hf = +8.5 to +10.1 and  Nd = +3.1 to +4.5. Most analyses define a restricted field defined by

  ). It may therefore be regarded as little contaminated at shallow level and as representing an unmodified plume composition. Representative isotopic values for this component are  Nd = +4.4,  Hf = +9.1,  Sr/ 86 Sr = 0.70411, 206 Pb/ 204 Pb =18.97, 207 Pb/ 204 Pb = 15.59, and 208 Pb/ 204 Pb = 39.03.

  [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF] described ankaramite as an exception, not just because it is the only sample on the island to exhibit very large clinopyroxene phenocrysts, but also because its chemical and Sr isotope compositions make it difficult to assign this lava to any of the two fractionation trends. Its high CaO/Al 2 O 3 and Sc/Yb ratios, combined with very low FeO/MgO, indicate the presence of pyroxene in the cumulus. These authors suggested that ankaramite cumulates had a deeper origin than other rock types but our new isotopic evidence

(± clinopyroxene) fractionation dominates the transitional lavas of the shield stage; (2) gabbroic (olivine+clinopyroxene+plagioclase) fractionation dominates the early, mildly alkalic lavas. From the remarkable buffering of Mg and Ni concentrations, Albarède and co-authors (1997) concluded that fractionation takes place not in a magma chamber but by percolation of melts through thick piles of cumulates. presents a slightly different perspective. With respect to the rest of Piton de la Fournaise lavas, this sample has particularly radiogenic Sr, unradiogenic 208 Pb/ 204 Pb, and essentially normal Hf

here to download Figure: Figure1-Bosch.pdf

	Figure1
	1400 Click 1000 1800 2200 Elevation (m)						Les Remparts Nez de Boeuf Morne Langevin Historical lavas	P F O U R N A I I T O N D E L A
								S
								E
	ISLAND REUNION Age (ka)	5 km 1300 800	Saint Denis La Montagne	Piton des Neiges	19°95	Intermediate lavas Younger lavas	E D O N P I T
								S
								N
								E
		2300 1800	Cirque PDN Salazie Cirque Mafate	Saint Benoit	Piton de la Fournaise	La Montagne	I G E S
			18.75 18.80 18.85 18.90 18.95 19.00 19.05	38.70	38.80	38.90	39.00
	Ocean Indian		Cilaos Cirque	Plaine des Sables				208 Pb/ 204 Pb
					Enclos
	La Montagne				PDF
	Langevin Nez de Boeuf Remparts Younger lavas Intermediate lavas	PDF PDN	Saint Pierre Remparts Rivière des	Langevin Rivière de	21°20
	Historical lavas
	Summit		55°30	55°40		55°50
			FIGURE 1

Figure 7: 143 Nd/ 144 Nd versus 176 Hf/ 177 Hf for Réunion Island lavas. The OIB, Hawaii islands and Kerguelen fields are from the GEOROC database [http://georoca.mpch-mainz.gwdg.de/georoc/].

Acknowledgements

We thank B. Galland for maintenance of the chemical laboratory of Géosciences Montpellier and P. Brunet (LMTG, Toulouse) for running Sr isotope analyses. We acknowledge Richard Carlson for the editorial handling and two anonymous reviewers for constructive comments. We are grateful to Ivan Vlastelic for providing unpublished Nd isotope data from historical Réunion flows and for fruitful discussions.

DMM, HIMU, EM-1, EM-2 and C are from [START_REF] Zindler | Chemical geodynamics[END_REF], [START_REF] Salters | Hf isotope constraints on mantle evolution[END_REF] and [START_REF] Hanan | Depleted Iceland mantle plume geochemical signature: artifact of multi-component mixing?[END_REF]. Symbols as in Fig. 4. [http://georoca.mpch-mainz.gwdg.de/georoc/]. NHRL: Northern Hemisphere Reference Line [START_REF] Hart | The DUPAL anomaly: a large-scale isotope anomaly in the southern hemisphere mantle[END_REF]. DMM, HIMU, EM-1, EM-2 and C are from [START_REF] Zindler | Chemical geodynamics[END_REF], [START_REF] Salters | Hf isotope constraints on mantle evolution[END_REF] and [START_REF] Hanan | Depleted Iceland mantle plume geochemical signature: artifact of multi-component mixing?[END_REF]. Symbols as in Fig. 4.

Table Caption

Table 1: Pb, Sr, Nd and Hf isotope ratios for lavas from Piton des Neiges and Piton de la Fournaise volcanoes. Uncertainties reported on Sr, Nd and Hf measured isotope ratios are 2n 1/2 analytical errors in last decimal place, where n is the number of measured isotopic ratios.

Analytical errors on Pb isotope ratios are better than 300 ppm. The 87 Sr/ 86 Sr isotope ratios labelled # are from [START_REF] Albarède | The geochemical regimes of Piton de la Fournaise volcano (Réunion) during the last 530,000 years[END_REF]. The 143 Nd/ 144 Nd ratios labelled § were analysed by MC-ICP-MS.  Sr ,  Nd and  Hf were calculated using ( 87 Sr/ 86 Sr) CHUR(0) = 0.7045, ( 143 Nd/ 144 Nd) CHUR(0) = 0.512638 and ( 176 Hf/ 177 Hf) CHUR(0) = 0.282772, respectively. Italic ratios are duplicates. -: not analysed.

Appendix A. Supplementary data are available online at the following address... Supplemental-File-1: Summary of the different parameters used during isotope analyses. *: references giving details about parameters used during isotope measurements with MC-ICP-MS (i.e., mass-fractionation, standard bracketing, precision, accuracy, correction of instrumental bias, standard,…). $ : internal mass bias correction using exponential correction to 205 Tl/ 203 Tl ratios. ¥ : example of NBS981 standard [START_REF] Todt | Evaluation of a 202 Pb-205 Pb double spike for high-precision lead isotope analysis[END_REF]