text:
# Structural study of an amorphous NiZr2 alloy by anomalous wide angle X-ray scattering and Reverse Monte Carlo simulations

## I Introduction

Amorphous materials have a great potential for application in technological devices, but their uses are limited due to several factors. One of them is the difficulty to obtain information about their atomic structures. Since most chemical and physical properties of these materials are defined by the first neighbors, knowledge of their structures becomes necessary. Due to the lack of translational symmetry, the determination of their atomic structures involves the combination of different diffraction and spectroscopic techniques as well as simulations and modeling.

The structure of an amorphous binary alloy is described by three pair correlation functions $`G_{ij}(r)`$, which are the Fourier transformation of the three partial structure factors $`𝒮_{ij}(K)`$. The total structure factor $`𝒮(K,E)`$, which can be derived from scattering measurements, is a weighted sum of these three $`𝒮_{ij}(K)`$ factors Faber and Ziman (1965). Thus, to determine the three partial $`𝒮_{ij}(K)`$ at least three independent $`𝒮(K,E)`$ factors are needed. Usually, the following methods have been used to obtain these $`𝒮(K,E)`$ factors:

1. Isomorphous substitution: in this method one component of the amorphous alloy is substituted by a (expected) chemically similar one. For example, in amorphous Ni-Zr alloys the Zr atoms are replaced by Hf atoms or the Ni atoms are partially substituted by Co or Fe atoms Lee et al. (1985). This method is not very accurate since one cannot be sure of the chemical similarity of the two elements involved.
2. Isotope substitution: in this method one component of the amorphous alloy is substituted by an isotope or by an isotopic mixture and the scattering measurements are carried out using neutron sources. For example, in amorphous Ni-Zr alloys the Ni atoms are replaced by a Ni isotope Mizoguchi et al. (1985). The advantages of this method are limited by the availability of adequate isotopes.

With the development of the synchrotron radiation sources, anomalous wide angle x-ray scattering (AWAXS) and differential anomalous scattering (DAS) techniques became available for structural study of multicomponent disordered materials. AWAXS utilizes an incident radiation that is tuned close to an atomic absorption edge so that it interacts resonantly with the electrons of that particular atom. The atomic scattering factor of each chemical component can therefore be varied individually and the chemical environment about each component in the material can be investigated. Thus, in the case of an amorphous binary alloy, by using only one sample three independent $`𝒮(K,E)`$ factors can be obtained. However, the matrix formed by the weights of the three $`𝒮(K,E)`$ factors is ill-conditioned, becoming difficult the determination of the three $`𝒮_{ij}(K)`$ factors.

FuossFuoss (1980); Fuoss et al. (1981) tried to overcome this difficulty by implementing the DAS approach, which was proposed by Schevick Shevchik (1977a, b). The DAS approach consists in taking the difference between the scattering patterns measured at two incident photon energies just below the edge of a particular atom, so that all correlation not involving this atom subtract out since only the atomic scattering factor of this atom changes appreciably. Later, De Lima et al. de Lima et al. (1988), following a suggestion made by Munro Munro (1982), combined the differential scattering factors $`DSF_i(K,E_m,E_n)`$ and the $`𝒮(K,E)`$ factors. They observed that this combination reduces the conditioning number of the matrix formed by the weights of these factors, allowing more stable values of $`𝒮_{ij}(K)`$ to be obtained.

The reverse Monte Carlo (RMC) simulation technique McGreevy and Pusztai (1988); McGreevy (1995); RMC ; McGreevy (2001) has been successfully used for structural modeling of amorphous structures. $`𝒮(K,E)`$ factors or their Fourier transformations, the total reduced radial distribution functions $`G(r)`$, can be used as input data. Applications of this technique to polymeric Rosi-Schwartz and Mitchell (1994), crystalline Mellergard and McGreevy (1999) and amorphous Karlsson et al. (2000); Jóvári and Pusztai (2001); Keen and McGreevy (1990); Wang et al. (1997); Bionducci et al. (1999); Wicks and McGreevy (1995); Cicco et al. (2000); Machado et al. (2002); Tengroth et al. (2001); Iparraguirre et al. (1993); Ohkubo et al. (2001); Pusztai and Sváb (1993) materials are described in the literature.

In this paper, we report the determination of the three $`𝒮_{ij}(K)`$ factors obtained for amorphous NiZr2 alloy by two independent ways: by making a combination of AWAXS and DAS techniques and also by combining AWAXS and RMC simulations techniques.

## II Theoretical background

### II.1 Total and partial structure factors

According to Faber and Ziman Faber and Ziman (1965), the $`𝒮(K,E)`$ factor is obtained from the scattered intensity per atom $`I_a(K,E)`$ as follows:

$`𝒮(K,E)`$ $`=`$ $`{\displaystyle \frac{I_a(K,E)\left[f^2(K,E)f(K,E)^2\right]}{f(K,E)^2}},`$ (1)
$`=`$ $`{\displaystyle \underset{i=1}{\overset{n}{}}}{\displaystyle \underset{j=1}{\overset{n}{}}}W_{ij}(K,E)𝒮_{ij}(K),`$
where $`K`$ is the transferred momentum, $`E`$ is the incident photon energy, $`𝒮_{ij}(K)`$ are the partial structure factors and $`W_{ij}(K,E)`$ are given by

$$W_{ij}(K,E)=\frac{c_ic_jf_i(K,E)f_j(K,E)}{f(K,E)^2},$$
(2)
and

$`f^2(K,E)`$ $`=`$ $`{\displaystyle \underset{i}{}}c_if_i^2(K,E),`$
$`f(K,E)^2`$ $`=`$ $`\left[{\displaystyle \underset{i}{}}c_if_i(K,E)\right]^2.`$
Here, $`c_i`$ is the concentration and $`f_i(K,E)=f_0(K)+f^{}(E)+if^{\prime \prime }(E)`$ is the atomic scattering factor of the component of type $`i`$ and $`f^{}(E)`$ and $`f^{\prime \prime }(E)`$ are the anomalous dispersion terms. The total reduced distribution function $`G(r)`$ is related to the $`𝒮(K,E)`$ factor through the Fourier transformation

$$G(r)=\frac{2}{\pi }_0^{\mathrm{}}K\left[𝒮(K,E)1\right]\mathrm{sin}(Kr)𝑑K,$$

while the partial reduced distribution functions $`G_{ij}(r)`$ are related to the $`𝒮_{ij}(K)`$ by means of the Fourier transformation

$$G_{ij}(r)=\frac{2}{\pi }_0^{\mathrm{}}K\left[𝒮_{ij}(K)1\right]\mathrm{sin}(Kr)𝑑K.$$

From $`G_{ij}(r)`$ the partial radial distribution functions $`\text{RDF}_{ij}(r)`$ can be calculated by

$$\text{RDF}_{ij}(r)=4\pi \rho _0c_jr^2+rG_{ij}(r).$$

Interatomic distances are obtained from the maxima of the $`G_{ij}(r)`$ functions and coordination numbers are calculated by integrating the peaks of the $`\text{RDF}_{ij}(r)`$ functions.

### II.2 Differential structure factors and differential radial distribution functions

By using the Faber and Ziman formalism Faber and Ziman (1965), the $`DSF_i(K,E_m,E_n)`$ factor around the component $`i`$ is obtained as follows

$`DSF_i(K,E_m,E_n)`$ $`=`$ $`{\displaystyle \frac{\left[I_a(K,E_m)I_a(K,E_n)\right]\left[LS(K,E_m)LS(K,E_n)\right]}{f(K,E_m)^2f(K,E_n)^2}},`$ (3)
$`=`$ $`{\displaystyle \underset{j}{}}U_{ij}(K,E_m,E_n)𝒮_{ij}(K),`$
where

$$LS(K,E)=f^2(K,E)f(K,E)^2,$$

is the Laue scattering term and

$$U_{ij}(K,E_m,E_n)=\frac{c_ic_j\left[f_i(K,E_m)f_j(K,E_m)f_i(K,E_n)f_j(K,E_n)\right]}{f(K,E_m)^2f(K,E_n)^2}.$$
(4)
For $`ij`$ the weight $`U_{ij}(K,E_m,E_n)`$ must be multiplied by 2. The differential distribution function $`DDF_i(r)`$ can be calculated by

$$DDF_i(r)=4\pi \rho r^2+\frac{2r}{\pi }_0^{\mathrm{}}K\left[DSF_i(K,E_m,E_n)1\right]\mathrm{sin}(Kr)𝑑K.$$

For a binary amorphous alloy, as compared to the total radial distribution function $`\text{RDF}(r)`$, which is a sum of three $`\text{RDF}_{ij}(r)`$ functions, this function is the sum of only two $`\text{RDF}_{ij}(r)`$ functions and is sensitive only to the environment around the atoms of type $`i`$. In this way it is similar to EXAFS, though there are important differences due to the different sections of $`K`$-space measured Kortwright et al. (1982).

### II.3 The matrix inversion method

In order to determine the three $`𝒮_{ij}(K)`$ factors using the eight $`𝒮(K,E)`$ and four $`DSF_i(K,E_m,E_n)`$ factors, a least-squares method was used. First of all, we define the function

$$(K)=\underset{\mathrm{}=1}{\overset{M}{}}P_{\mathrm{}}\left[𝒮_{\mathrm{}}(K)\underset{i=1}{\overset{2}{}}\underset{j=1}{\overset{2}{}}h_{ij}^{\mathrm{}}(K)𝒮_{ij}(K)\right]^2.$$

For $`1\mathrm{}8`$, $`𝒮_{\mathrm{}}(K)`$ indicates the $`𝒮(K,E_{\mathrm{}})`$ factors and for $`9\mathrm{}12`$ it stands for the $`DSF_i^{\mathrm{}}(K,E_m,E_n)`$ factors, as defined by Eqs. 1 and 3. In the same way, for $`1\mathrm{}8`$, $`h_{ij}^{\mathrm{}}(K)`$ indicates the $`W_{ij}(K,E_{\mathrm{}})`$ weights and for $`9\mathrm{}12`$ it stands for the $`U_{ij}^{\mathrm{}}(K,E_m,E_n)`$ weights, as seen in Eqs. 2 and 4. $`P_{\mathrm{}}`$ is a weight which allows us to vary the contribution of any $`𝒮_{\mathrm{}}(K)`$ factor. Here, $`P_{\mathrm{}}=1`$, $`\mathrm{}`$, and $`M=12`$. By taking the derivatives of the $`(K)`$ function with respect to the $`𝒮_{ij}(K)`$ factors and putting them equal zero a linear system formed by three equations can be found. These three linear equations can be written as a matrix equation through

$$\stackrel{~}{𝒮}=\stackrel{~}{W}\stackrel{~}{𝒮}_{ij},$$

and by inverting this matrix equation the $`𝒮_{ij}(K)`$ can be found. There are several methods to invert a matrix, and we have adopted the Singular Value Decomposition (SVD) method Press et al. (1986).

### II.4 The RMC method

The basic idea and the algorithm of the standard RMC method are described elsewhereMcGreevy and Pusztai (1988); McGreevy (1995); RMC ; McGreevy (2001) and its application to different materials is reported in the literature Rosi-Schwartz and Mitchell (1994); Mellergard and McGreevy (1999); Karlsson et al. (2000); Jóvári and Pusztai (2001); Keen and McGreevy (1990); Wang et al. (1997); Bionducci et al. (1999); Wicks and McGreevy (1995); Cicco et al. (2000); Machado et al. (2002); Tengroth et al. (2001); Iparraguirre et al. (1993); Ohkubo et al. (2001); Pusztai and Sváb (1993). In the RMC procedure, a three-dimensional arrangement of atoms with the same density and chemical composition of the alloy is placed into a cell (usually cubic) with periodic boundary conditions and the $`G_{ij}^{\text{RMC}}(r)`$ functions corresponding to it are directly calculated through

$$G_{ij}^{\text{RMC}}(r)=\frac{n_{ij}^{\text{RMC}}(r)}{4\pi \rho _0r^2\mathrm{\Delta }r},$$

where $`n_{ij}^{\text{RMC}}(r)`$ is the number of atoms at a distance between $`r`$ and $`r+\mathrm{\Delta }r`$ from the central atom, averaged over all atoms. By allowing the atoms to move (one at each time) inside the cell, the $`G_{ij}^{\text{RMC}}(r)`$ functions can be changed and, as a consequence, the $`𝒮_{ij}(K)`$ and $`𝒮^{\text{RMC}}(K)`$ factors are changed. Thus, the $`𝒮^{\text{RMC}}(K,E)`$ factor is compared to the $`𝒮(K,E)`$ factor in order to minimize the differences between them. The function to be minimized is

$$\psi ^2(E)=\frac{1}{\delta }\underset{i=1}{\overset{m}{}}\left[𝒮(K_i,E)𝒮^{\text{RMC}}(K_i,E)\right]^2,$$

where the sum is over $`m`$ experimental points and $`\delta `$ is related to the experimental error in $`𝒮(K,E)`$. If the movement decreases $`\psi ^2`$, it is always accepted. If it increases $`\psi ^2`$, it is accepted with a probability given by $`\mathrm{exp}(\mathrm{\Delta }\psi ^2/2)`$; otherwise it is rejected. As this process is iterated $`\psi ^2`$ decreases until it reaches an equilibrium value. Thus, the atomic configuration corresponding to equilibrium should be consistent with the experimental total structure factor within the experimental error. By using the $`G_{ij}^{\text{RMC}}(r)`$ functions and the $`𝒮_{ij}^{\text{RMC}}(K)`$ factors the coordination numbers and interatomic distances can be calculated. In addition, the bond-angle distributions can also be determined.

## III Experimental Procedure

### III.1 Sample and holder

Amorphous NiZr2 thin films were prepared on a NaCl cooled substrate by the sputtering triode technique at the Argonne National Laboratory. In order to eliminate the substrate a metallic holder with a central rectangular hole was glued on the thin film, and the set was immersed in water to dissolve the NaCl substrate. The thickness and density of the thin films were 2 $`\mu `$m and 7.32 g/cm3, respectively.

### III.2 Apparatus and data collection

The AWAXS scattering experiments were performed some years ago on the D4 beamline at LURE (Orsay - France) using the DCI synchrotron radiation source. All measurements were performed using a two-circles diffractometer in the vertical plane described elsewhereLaridjani et al. (1985); M.Laridjani et al. (1987) equipped with a two crystal Si (220) or Ge (400) monochromator and a Si:Li energy sensitive detector which is sufficient to discriminate the large Kα resonant Raman/fluorescence signal when incident photon energy is tuned close to the Ni and Zr K edges but insufficient to distinguish the small elastic, Kβ and Compton signals. The Si (220) and Ge (400) crystal monochromators were used to collect data around the Ni and Zr edges, respectively. The sample was sealed, under vacuum, into a cell containing a large kapton window fixed around the diffractometer horizontal rotation axis to reduce the air scattering. In order to take into account the time-decrease of the beam, the incident intensity $`I_0`$ was measured by placing a NaI photomultiplier tube between the monochromator and the diffractometer. Absorption measurements were performed on the Ni and Zr K edges for energy calibration. For each incident photon energy listed in Table 1, the energy windows for Kα, resonant Raman/fluorescence, Kβ, Compton, and elastic signals were defined by measuring the scattered intensities at two different angular positions: the beginning of the windows was defined at the lowest scattering angle and their ends at the highest one. Scattering patterns were collected using the symmetrical reflection geometry. To improve the quality of the data, at least four scattering measurements were performed at the incident photon energies shown in Table 1.

### III.3 Determination of the anomalous scattering factors

In order to interpret the scattering data correctly, the real and the imaginary parts $`f^{}`$ and $`f^{\prime \prime }`$ of the atomic scattering factor were determined accurately following a procedure described by Dreier et al.Dreier et al. (1984) X-ray absorption coefficients were measured over a broad energy range near the Ni and Zr K edges on the sample and $`f^{\prime \prime }`$ was calculated using the optical theorem. The absorption measurements were performed on the XAS beamline at LURE, in transmission mode, and two standard ionization chambers were used as the detecting system. Outside the region of measurement, theoretical values mcm  of $`f^{\prime \prime }`$ were used to extend the experimental data set over a larger energy range and $`f^{}`$ was calculated using the Kramers-Kronig relation. For the measurements away from the K edges, the $`f^{}`$ and $`f^{\prime \prime }`$ values were taken from a table compiled by Sasaki Sasaki (1984). The obtained values for the incident photon energies 8330 eV and 17987 eV are listed in Table 1 together with Sasaki values. The atomic scattering factor away from the K edge, $`f_0(K)`$, of neutral Ni and Zr atoms were calculating according to the analytic function given by Cromer and Mann Cromer and Mann (1968).

### III.4 Data analysis

Each scattering pattern was treated separately. Due to the diffractometer characteristics, the polarization correction was disregarded. The measured scattered intensities away from the Ni and Zr K edges were corrected for detector non-linearities and reabsorption effects being then put on a per atom scale and the Compton scattering eliminated. For those measured close to the Ni and Zr edges, besides the corrections above, the contribution of the Kβ fluorescence to the elastic scattering was eliminated. This contribution was evaluated, after corrections for reabsorption effects, by using the Kβ/Kα ratio, which was measured at an incident photon energy above the edge, because they were completely separated. Two method were tested to put the scattered intensities on a per atom scale: that described by Krogh-Moe-Norman Krogh-Moe (1956); Norman (1957); Wagner (1972) and the one based on the high-angle procedure Wagner (1972); Gingrich (1943). For the elastic scattered intensity around the Zr edge, for which the $`K`$ data range is larger ($`K_{\text{max}}13`$ Å-1), both method gave the same result. The $`K`$ data range was reduced and the Krogh-Moe-Norman method was then applied. The normalized signal was identical to the one obtained previously. So, we have kept the Krogh-Moe-Norman method to put the scattered intensities on a per atom scale including those around the Ni edge, for which the $`K`$ data range is smaller ($`K_{\text{max}}8`$ Å-1). The good quality of the total structure factors shown in Fig. 1 and 2 corroborates the accuracy of the procedure used. The Compton scattering contribution was calculated according to the analytic approximation given by Pálinkas Pálinkas (1973).

## IV Results and Discussion

### IV.1 Total structure factors, differential structure factors and differential distribution functions

Figures 1 and 2 show the $`𝒮(K,E)`$ factors for the incident photon energies described in Table 1. The $`𝒮(K,E)`$ factors for the energies greater than 16101 eV (see Fig. 2) show a pre-peak located at about $`K=1.65`$ Å-1, which reaches the maximum intensity close to the Zr edge. It is absent in the $`𝒮(K,E)`$ factors for the smallest energies (see Fig. 1). Since $`𝒮(K,E)`$ is a weighted sum of $`𝒮_{ij}(K)`$, we calculated the weights $`W_{ij}(K,E)`$ for the $`𝒮_{\text{Ni-Ni}}(K)`$, $`𝒮_{\text{Ni-Zr}}(K)`$ and $`𝒮_{\text{Zr-Zr}}(K)`$ factors. Their contributions for the $`𝒮(K,E)`$ factor at 8330 eV are about 4%, 32% and 64% while for the $`𝒮(K,E)`$ factor at 17987 eV they are about 10%, 43% and 47%, respectively. Thus, the determination of Ni-Ni first neighbors is more difficult than for the Ni-Zr and Zr-Zr first neighbors.

We have used the scattered intensities on a per atom scale obtained at the incident photon energy pairs (8210, 8330 eV), (8330, 9659 eV), (16101, 17987 eV) and (17987, 19200 eV) to calculate four $`DSF_i(K,E_m,E_n)`$ factors around the Ni and Zr atoms. Figures 3 and 4 show them and their $`DDF_i(r)`$ functions. Figure 3 shows that the chemical environment around Ni and Zr atoms is very different. For instance, it is seen a broad shoulder on the right side of the first peak in the $`DSF_{\text{Ni}}(K,E_m,E_n)`$ factors which is absent in the $`DSF_{\text{Zr}}(K,E_m,E_n)`$ factors. The $`DDF_{\text{Zr}}(r)`$ functions show the first coordination shell split into two subshells, located at about $`r=2.64`$ Å and 3.23 Å, while the $`DDF_{\text{Ni}}(r)`$ functions display the first coordination shell located at a mean distance $`r=2.77`$ Å. The contribution of the Ni-Ni and Ni-Zr pairs to them is not resolved. This can be related to the small $`K`$ data range of the measured scattered intensities. The weight $`U_{\text{Ni-Ni}}(K)`$ to the $`DSF_{\text{Ni}}(K,E_m,E_n)`$ factor is about 25%, and it is almost four times larger than the weight $`W_{\text{Ni-Ni}}(K)`$ for the $`𝒮(K,E)`$ factors, allowing a more realistic determination of the number of Ni-Ni first neighbors.

By considering the atomic radius of Ni and Zr atoms, the subshells located at the distances $`r=2.64`$ Å and 3.23 Å in the $`DDF_{\text{Zr}}(r)`$ functions were directly attributed to Zr-Ni and Zr-Zr correlations. This fact allowed us to deconvolute their first shells and, consequently, also the first shells of the $`DDF_{\text{Ni}}(r)`$ functions, by assuming a Gaussian shape for them. To do this, we explicitly took into account the $`K`$ dependence by Fourier transforming into the $`K`$ space the Gaussian distributions multiplied by the weights $`U_{ij}(K,E_m,E_n)`$ and back Fourier transforming over the same $`K`$ range used to calculate the $`DDF_i(r)`$. The structural parameters extracted from the deconvolution are listed in Table 2.

### IV.2 Partial structure factors obtained from the matrix inversion and RMC methods

#### IV.2.1 Partial $`𝒮_{ij}(K)`$ factors obtained from the matrix inversion method

By considering the eight $`𝒮(K,E)`$ and the four $`DSF_i(K,E_m,E_n)`$ factors, a reduction of a factor of 15 (at $`K=1.65`$ Å-1) and 7 (at $`K=7.65`$ Å-1) in the conditioning number of the matrix $`\stackrel{~}{W}`$ was reached when compared to the use of the eight $`𝒮(K,E)`$ factors only. To obtain the three $`𝒮_{ij}(K)`$ factors, the $`𝒮(K,E)`$ and $`DSF_i(K,E_m,E_n)`$ factors were restricted to $`K_{\text{max}}=7.65`$ Å-1. The $`𝒮_{\text{Ni-Ni}}(K)`$, $`𝒮_{\text{Ni-Zr}}(K)`$ and $`𝒮_{\text{Zr-Zr}}(K)`$ factors (thin lines) obtained are shown in Fig. 5. From this figure, the $`𝒮_{\text{Ni-Zr}}(K)`$ and $`𝒮_{\text{Zr-Zr}}(K)`$ factors seem to be of good quality, specially between $`K=1.5`$ Å-1 and $`K=3.5`$ Å-1, where the conditioning number presents the largest changes. On the other hand, the $`𝒮_{\text{Ni-Ni}}(K)`$ factor displays a minimum apparently without physical meaning in this region. The $`𝒮_{\text{Ni-Zr}}(K)`$ factor shows a minimum at $`K=2.2`$ Å-1, which was already observed in the $`DSF_{\text{Ni}}(K,E_m,E_n)`$ factors, while the $`𝒮_{\text{Zr-Zr}}(K)`$ factor displays a shoulder at this $`K`$ value. Although the $`𝒮_{\text{Ni-Ni}}(K)`$ factor is very noisy, it is possible to see a weak pre-peak at about $`K=1.5`$ Å-1 which was already observed in the $`𝒮(K,E)`$ factors around the Zr edge and also some weak halos beyond $`K=3.5`$ Å-1.

Figure 6 shows the $`\text{RDF}_{ij}(r)`$ functions (thin lines) obtained by Fourier transformation of the $`𝒮_{ij}(K)`$ factors. From this picture one can see that the RDF$`{}_{\text{Ni-Ni}}{}^{}(r)`$ function shows a very weak first coordination shell, located at about $`r=2.64`$ Å, while the second and third ones, located at about $`r=4.5`$ Å and 6.7 Å, respectively, are very intense. The RDF$`{}_{\text{Ni-Zr}}{}^{}(r)`$ function displays a well isolated first coordination shell, located at about $`r=2.77`$ Å; a small shoulder at about $`r=3.83`$ Å, and well defined and intense second and third coordination shells, located at about $`r=5.4`$ Å and 7.9 Å. The RDF$`{}_{\text{Zr-Zr}}{}^{}(r)`$ function shows the first coordination shell well isolated, located at about $`r=3.3`$ Å; the second one splits into two subshells, located at about $`r=5.18`$ Å and 6.35 Å, which are very weak. The interatomic distances and coordination numbers for the first neighbors obtained from these functions are listed in Table 2.

#### IV.2.2 Partial $`𝒮_{ij}(K)`$ obtained from the RMC simulations

In order to perform the simulations we have considered a cubic cell with 1800 atoms (600 Ni and 1200 Zr), $`\delta =0.01`$, and a mean atomic number density $`\rho _0=0.054719`$ atoms/Å3. The minimum distances between atoms were fixed at the beginning of the simulations at $`r_{\text{Ni-Ni}}=2.2`$ Å, $`r_{\text{Ni-Zr}}=2.4`$ Å and $`r_{\text{Zr-Zr}}=2.6`$ Å. To make the simulations we used the RMC programs available at the Internet RMC . The $`𝒮(K,E)`$ factors obtained at the photon energies 8330 eV, 17398 eV and 19200 eV were used as input data. These $`𝒮(K,E)`$ factors were cut at the same $`K`$ value. The $`𝒮(K,E)`$ (solid lines) and $`𝒮^{\text{RMC}}(K,E)`$ factors (square lines) are shown in Fig. 7 and they show a good agreement.

The $`𝒮_{ij}^{\text{RMC}}(K)`$ factors (thick lines) are also shown in Fig. 5, together with those obtained from the matrix inversion method. The agreement among them is very good, mainly for the $`𝒮_{\text{Zr-Zr}}(K)`$ factor. For instance, the shoulder previously observed in $`𝒮_{\text{Zr-Zr}}(K)`$ factor and the peaks located around $`K=2.58`$ Å-1, 4.3 Å-1 and 6.2 Å-1 are also seen in the $`𝒮_{\text{Zr-Zr}}^{\text{RMC}}(K)`$ factor; the minimum at $`K=2.1`$ Å-1 and peaks around $`K=2.7`$ Å-1, 5.2 Å-1 and 7.6 Å-1 observed in the $`𝒮_{\text{Ni-Zr}}(K)`$ are also seen in the the $`𝒮_{\text{Ni-Zr}}^{\text{RMC}}(K)`$ factor; and it is interesting to note that the minimum previously observed in the $`𝒮_{\text{Ni-Ni}}(K)`$ factor is now replaced by a well defined maximum in the $`𝒮_{\text{Ni-Zr}}^{\text{RMC}}(K)`$ factor. Other peaks observed before in the $`𝒮_{\text{Ni-Ni}}(K)`$ factor are also seen in the $`𝒮_{\text{Ni-Ni}}^{\text{RMC}}(K)`$ factor. The RDF$`{}_{ij}{}^{\text{RMC}}(r)`$ (thick lines) functions are also shown in Fig. 6. The interatomic distances and coordination numbers for the first neighbors extracted from the $`G_{ij}^{\text{RMC}}(r)`$ and RDF$`{}_{ij}{}^{\text{RMC}}(r)`$ functions are also listed in Table 2.

It is interesting to note that although we have used eight independent $`𝒮(K,E)`$ and four $`DSF_i(K,E_m,E_n)`$ factors to obtain the three $`𝒮_{ij}(K)`$ factors by the matrix inversion method, the $`𝒮_{\text{Ni-Ni}}(K)`$ factor shows a minimum between $`K=1.5`$ Å-1 and 3.5 Å-1 which have no physical meaning. On the other hand, using only three $`𝒮(K,E)`$ factors as input data in the RMC method we have obtained excellent results, even for the $`𝒮_{\text{Ni-Ni}}(K)`$ factor.

In a previous study, Lee et al. Lee et al. (1985) and Mizoguchi et al. Mizoguchi et al. (1985) reported the $`𝒮_{ij}(K)`$ factors obtained for the amorphous Ni35Zr65 and Ni36Zr64 alloys, respectively. Lee et al used alloys containing mixtures of Co and Fe in small quantities while Mizoguchi et al used the isotopic substitution method. A comparison among the $`𝒮_{ij}^{\text{RMC}}(K)`$, $`G_{ij}^{\text{RMC}}(r)`$ and RDF$`{}_{ij}{}^{\text{RMC}}(r)`$ functions obtained in this work, by using only one sample, with those reported by them shows that they are in a good agreement.

By defining the partial bond-angle distribution functions $`\mathrm{\Theta }_{ijl}(\mathrm{cos}\theta )`$ where $`j`$ is the atom in the corner, we calculated the angular distribution of the bonds between first neighbor atoms. The six $`\mathrm{\Theta }_{ijl}(\mathrm{cos}\theta )`$ functions found are shown in Fig. 8. This kind of information cannot be obtained by the matrix inversion method because the partial distribution functions give only a one-dimensional description of the atomic arrangement.

Hausleitner and Hafner Hausleitner and Hafner (1992) investigated several amorphous alloys formed by transition metals using molecular dynamics simulations. They obtained the structure factors, coordination numbers, interatomic distances and also the bond-angle distribution functions. According to them, if the components of an alloy have a large difference in the number of $`d`$ electrons there is a pronounced non-additivity of the pair interactions and a strong interaction between unlike atoms. Therefore, in certain amorphous alloys the interatomic distance of heteropolar bonds can be shorter than the distance of homopolar ones Machado et al. (2002); Fukunaga et al. (1984). This also favors the formation of trigonal-prismatic units. The shortening effect was not observed in this study. The bond angle values of triangular and square faces of a distorted prism are found around $`\theta =60^{}`$, $`90^{}`$$`100^{}`$, $`109^{}`$ and $`147^{}`$. Thus, the bond angle values shown in Fig. 8 suggest that distorted polyhedral units with triangular faces are present in the atomic structure of the amorphous NiZr2 alloy and if there are square-faced polyhedral units, they are found in a very small quantity. The presence of this kind of units was already reported in amorphous Ni60Ti40 Machado et al. (2002), Ni33Ti67 Iparraguirre et al. (1993), Pd82Si18 Ohkubo et al. (2001) and Ni65B35 Pusztai and Sváb (1993) alloys.

Havinga et al. Havinga et al. (1972) reported crystallographic data for several compounds with CuAl2 (C16)-type structure. We have used their data for the NiZr2 compound in the CRYSTAL OFFICE98 softwareAto  to obtain the interatomic distances and coordination numbers for the first neighbors listed in Table 2. Its structure has angles around $`\theta =53`$–56, $`64^{}`$, $`96`$$`112^{}`$ and $`143`$$`154^{}`$. For the amorphous NiZr2 alloy the $`\mathrm{\Theta }_{i\text{Zr}l}(\mathrm{cos}\theta )`$ functions ($`i`$ and $`l`$ can be Ni or Zr atoms) found are very similar to those present in the NiZr2 compound. On the other hand, the $`\mathrm{\Theta }_{i\text{Ni}l}(\mathrm{cos}\theta )`$ functions show some important differences. The Ni–Ni–Ni sequence in the compound is linear ($`\theta =180^{}`$) while in amorphous NiZr2 alloy it can also be found at triangular angles ($`\theta =60^{}`$ and $`112^{}`$); the Zr–Ni–Ni and Zr–Ni–Zr sequences found in the amorphous NiZr2 alloy show bond angles around $`\theta =180^{}`$ that are not found in the compound. These data suggest that in the amorphous alloy the local structure around a Ni atom can be more distorted than around a Zr atom.

The chemical short-range order parameter is a quantitative measurement of the degree of chemical short-range order. There are several definitions of chemical short-range order parameters Steeb and Lamparter (1992). We have adopted the Warren chemical short-range order (CSRO) parameter $`\alpha _1`$ Maret et al. (1987) given below to compare the amorphous NiZr2 alloy and the NiZr2 compound.

$$\alpha _1=1\frac{N_{12}^1}{c_2(c_1N_2^1+c_2N_1^1)},$$
(5)
Here, $`N_{12}^1=N_{\text{Ni-Zr}}`$ and $`N_1^1=N_{\text{Ni-Ni}}+N_{\text{Ni-Zr}}`$ and $`N_2^1=N_{\text{Zr-Zr}}+N_{\text{Zr-Ni}}`$ are the total coordination numbers in the first shell. By using the coordination numbers given in Table 2, we found $`\alpha _1^{\text{RMC}}=0.029`$ and for the compound $`\alpha _1^{\text{NiZr}_2}=0.024`$. The chemical short-range order in the amorphous NiZr2 alloy is similar to that found in a solid solution, for which $`\alpha _1=0`$. The same behavior is observed in the NiZr2 compound.

## V Conclusion

The amorphous NiZr2 alloy was investigated by using AWAXS, DAS and RMC simulations techniques. The three $`𝒮_{ij}(K)`$ factors were determined by using two independent ways: a combination of AWAXS and DAS techniques and by combining AWAXS and RMC simulations techniques. The $`𝒮_{\text{Ni-Zr}}(K)`$ and $`𝒮_{\text{Zr-Zr}}(K)`$ factors obtained by both combinations showed an excellent agreement, but the simulations gave the best result for the $`𝒮_{\text{Ni-Ni}}(K)`$ factor. The coordination numbers and interatomic distances calculated from the RDF$`{}_{ij}{}^{}(r)`$ functions obtained by these two combinations show a good agreement among themselves and with those extracted from the deconvolution of the $`DDF_{\text{Ni}}(r)`$ and $`DDF_{\text{Zr}}(r)`$ functions.

The results achieved in this study suggest that the combination of the RMC simulations method and AWAXS can be used to obtain more reliable $`𝒮_{ij}(K)`$ factors than those obtained by the matrix inversion method, in particular for the Ni-Ni pairs, which have the smallest weighting factor. This was already observed in a neutron diffraction study performed in a CuBr alloy Pusztai and McGreevy (1998), which is also an ill-conditioned problem similar to ours. In addition, this combination gives us some evidence of the presence of distorted polyhedral units with triangular faces in the atomic structure of the amorphous NiZr2 alloy.

The calculated Warren CSRO parameter indicated that the chemical short-range order found in the amorphous NiZr2 alloy is similar to that found in a solid solution, and there is some resemblance with that found in the NiZr2 compound.

###### Acknowledgements.

Two of the authors, J. C. De Lima and K. D. Machado, would like to thank the Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq), Brazil, for financial support.


format: text

text:
# Brown dwarfs in the Pleiades cluster : clues to the substellar mass function Based on observations obtained at Canada-France-Hawaii Telescope

## 1 Introduction

The Initial Mass Function (IMF), i.e. the mass spectrum resulting from complex physical processes at work during star formation, is a formidable constraint for star formation models. Its determination at low stellar and substellar masses is therefore one of the main motivation for the rapidly expanding quest for brown dwarfs. Very low mass stars and brown dwarfs are also prime candidates to investigate the structure and dynamical evolution of large stellar systems, such as e.g. clusters. Searches for the lowest mass, isolated objects have been conducted in the galactic field, in young open clusters and in star forming regions, brown dwarfs being brighter when younger. The identification of substellar candidates is often based on color magnitude diagrams (hereafter CMD) built from deep wide-field imaging surveys. One of the major shortcomings of this selection method is the contamination by older and more massive late-type field dwarfs which may lie in the same region of the CMD.

The Pleiades cluster (RA=$`3^h46.6^m`$, Dec=$`+24^o04^{}`$) is an ideal hunting ground for substellar objects. It is relatively nearby with a distance of about 125 pc, with the apparent magnitude of massive brown dwarfs bright enough for being easily detected on intermediate-size telescopes. Its age of approximately 120 Myr (Stauffer et al. 1998, Martin et al. 1998) makes the lithium test particularly useful for identifying brown dwarfs : at this age, the lithium depletion boundary (in mass) coincides with the hydrogen burning mass limit (hereafter HBML). Also, the Pleiades is a rich cluster with about 1200 members and its galactic latitude is relatively high ($`b=23^o`$) which minimizes the possible confusion between members and red giants from the Galactic disk. Moreover, the cluster motion ($`\mu _\alpha =+19`$ mas/yr, $`\mu _\delta =43`$ mas/yr) is large compared to that of field stars so that cluster kinematic studies are a powerful way to recognize true members among the photometric candidates.

To date, several surveys have been conducted in this cluster. They are summarized in Table 1 in term of covered area, completeness limit and number of new brown dwarfs candidates. Only surveys with a magnitude limit larger or equal to $`I17.5`$ are listed here (the HBML corresponds to $`Ic=17.8`$, Stauffer et al. 1998). Starting in 1997111Prior to 1997, mass function estimates were based on strongly contaminated samples of brown dwarf candidates., several groups reported estimates for the Pleiades mass function in the upper part of the substellar domain, often represented by a single power law $`dN/dMM^\alpha `$ with $`0.5\alpha 1.0`$ (see Table 1).

While these various estimates agree within uncertaintites, the substellar samples on which they rely are still relatively small and, in some cases, not fully corrected for contamination by field dwarfs. Some surveys are also quite limited spatially, and the derived mass function may not be representative of the whole cluster. In order to put the determination of the Pleiades substellar mass function on firmer grounds, we performed a deep and large survey of the Pleiades cluster in the $`I`$ and $`Z`$-band, using the CFH12K camera at the Canada-France-Hawaii telescope. The survey covers 6.4 square degrees and reaches up to 3 degrees from the cluster center. It encompasses a magnitude range $`I13.524`$, which corresponds to masses from 0.025$`M_{}`$ to 0.45$`M_{}`$ at the distance and age of the Pleiades.

In section 2, we describe the observations and the data reduction. Results in the stellar and substellar domains are presented in Section 3. Forty brown dwarf candidates are identified, 29 of which are new discoveries. We discuss contamination of the photometric samples by field stars in the stellar and substellar domains, investigate the radial distribution of substellar objects, and derive the cluster substellar mass function. In Section 4, we discuss the overall shape of the cluster mass function over the stellar and substellar domains, the possible consequences of early cluster dynamical evolution on its present mass function and the implications for the brown dwarf formation process. A comparison is made between the Pleiades and the Galactic field mass functions, which reveals significant differences at the low mass end, a large part of which is probably attributable to unresolved cluster binaries.

## 2 Observations and data reduction

The $`I`$ and $`Z`$ observations were made using the CFH12K mosaic camera (Cuillandre et al. 2001) at the prime focus of the Canada-France-Hawaii telescope, from December 17 to 19, 2000. The camera is equipped with an array of twelve $`2048\times 4096`$ pixel CCDs with $`15\mu m`$ size pixels, yielding a scale of 0.206”/pixel and a field of view of 42’$`\times `$28’. The $`I`$ Mould and $`Z`$ filter profiles are shown in Fig 1. We obtained images of 17 fields in the Pleiades whose coordinates are given in Table 2 and the location on the sky relative to the brightest Pleiades members shown in Fig 2. Fields were selected in order to avoid the brightest stars which would produce reflection halos and scattered light on the detectors, i.e., a higher and non uniform background which is difficult to remove. We also avoided the southwestern region of the cluster where a small CO cloud causes relatively large extinction. The observing conditions were photometric with a $`I`$-band seeing of typically 0.5” FWHM, but occasionally up to 1.1” FWHM.

Short and long exposures were taken for each field in order to cover a large magnitude range. We obtained a set of four images in both $`I`$ and $`Z`$ filters : one 10 sec exposure plus three 300 sec exposures in the $`I`$-band and one of 10 sec exposure plus three of 360 sec exposures in the $`Z`$-band. The short exposures encompass the magnitude range $`13.5I21.5`$ whereas the long ones cover $`17I24`$. A comparison between the number of stellar like objects detected in long exposures of overlapping fields indicates that the survey is about 90% complete down to $`I22`$ (see Fig 3). Note that the completeness limit of the short exposures ($`I19.5`$) is much larger than the saturation limit of the long exposures ($`I17`$), so that the survey completely samples the continuous mass range from $`0.03M_{}`$ to $`0.45M_{}`$ in the Pleiades cluster.

### 2.1 Photometry

Each mosaic consists of 12 CCD images which were reduced and analysed separately. The images were overscan corrected, debiased and flat-fielded. The flats were normalized to a reference CCD to retain the appropriate relative scaling between chips. The same photometric zero-point can thus be used for all the CCDs of a mosaic. The images were also fringe corrected using patterns derived from a smoothed combination of more than 20 images in each band. Then, images of similar exposure time for a given field were stacked together using an optimal CCD clipping to remove the artefacts, and then taking the average of the remaining final stack for each pixel to preserve the photometry.

To detect all the sources in the frames, we averaged all long exposure $`I`$ and $`Z`$ images of a field previously shifted to the same location and used the automatic object-finding algorithm from the Sextractor package (Bertin & Arnouts 1996). Then, a PSF-fitting photometry on the $`I`$ and $`Z`$ images was performed for all the detected objects using the PSFex package. We discriminated between point-like and extended or corrupted sources using the FWHM distribution of all the detected objects on each field. In practice, we defined conservative lower and upper limits around the well defined stellar peak of this distribution and rejected all the sources located outside these boundaries.

Photometric standard Landolt fields SA98 and SA113 were observed and reduced in the same way as the science images. Three successive exposures of SA98 were obtained with an incremental offset of several arcminutes in RA. Thus, common sets of photometric standard stars were observed on every CCD of the mosaic and we checked that the photometric zero-point in $`I`$-band was the same for each CCD : we did find a scatter of only 0.03 mag. The photometry in the $`Z`$-band is sensitive to detector pass band and Landolt does not give $`Z`$ magnitudes for his standards. Thus we selected unreddened A0 standards and set their $`Z`$ magnitude equal to their $`I`$ magnitude assuming that their $`IZ`$ color was zero (by definition of an A0 star). Thanks to the several exposures of SA98 shifted of a few arcminutes, those A0 stars were observed on different chips so that we could derive $`Z`$ zero-points for some of the CCDs. We found a scatter of only 0.04 mag and we then used the mean value as a global zero-point, assuming that it could be used for all the CCDs as in the $`I`$-band. We verified this assumption a posteriori by measuring the difference between the magnitudes of the same objects detected in the common region of 2 overlapping pointings. We thus verified that there was no systematic error and estimated the photometric rms error up to the completeness limit ($`I22`$, see Fig. 4, top panel) which amounts to 0.07 mag or less at $`I22`$.

### 2.2 Astrometry

In order to obtain accurate coordinates for cluster brown dwarf candidates, we had to derive an astrometric solution for each CCD of the mosaic. We used the CFHT’s Elixir package (Magnier & Cuillandre 2002) to compute the astrometric solution for each image. The algorithm calculates the celestial coordinates for all the detected objects using the approximate solution given by the header, compares them with the USNO2 catalog and refines the solution. Most of the USNO2 stars were saturated in our long exposure images so that we had to use lists of refined coordinates derived from the short exposures as a reference catalog.

We estimated the astrometric error by comparing coordinates of stars present in overlapping regions and found an accuracy better than 0.5 arcsecond. The astrometric rms error is shown in Fig 4 (bottom panel) as a function of $`I`$ magnitude.

## 3 Results

The $`(I,IZ)`$ color magnitude diagrams of the point-like objects contained in the 17 Pleiades fields are shown in Fig 5 and 6. The short and long exposures have been analysed separately corresponding respectively to the stellar and the substellar domains. In both cases we present our photometric selection of candidates before dealing with the field star contamination. We examine the spatial distribution of cluster members and attempt to measure their core radius. We then use these estimates to derive the Pleiades mass function.

### 3.1 Stellar domain

The $`(I,IZ)`$ color magnitude diagram for the short exposure images of our survey is presented in Fig 5. The 120 Myr isochrone from the NEXTGEN models of Baraffe et al. (1998) shifted to the Pleiades distance ($`(mM)_o=5.53`$) is shown as a dashed line. On the basis of the location of this theoretical isochrone, we made a rather conservative photometric selection to include all possible stellar members between $`I=13.5`$ and $`I=17.5`$. This selection corresponds to the box drawn in Fig 5.

To remove all the contaminating field stars from our sample, we compared our list of candidates to the results of other large Pleiades surveys. Adams et al. (2001) performed a large search for Pleiades stellar members using the photometry from the Two Micron All Sky Survey (2MASS) and proper motions determined from Palomar Observatory Sky Survey (POSS) plates. This search extends to a radius of 10$`\mathrm{°}`$ around the cluster center, well beyond the tidal radius, which means that it covers the complete cluster area. The completeness limit of the POSS plates is $`I16.5`$, i.e. $`0.1M_{}`$. The authors analysed the proper motion of all the objects previously selected on the basis of their 2MASS JHK photometry and defined a membership probability $`p`$ (see Adams et al. 2001 for details). We cross-correlated our list of stellar candidates with the list of all the sources analysed by Adams et al. (2001) and we kept all the objects with $`p>0.1`$ so as to minimize the non-member contamination down to $`I=16.5`$ (Adams’ survey completeness limit). All those sources are shown as filled triangles in Fig 5.

For stars fainter than $`I>16.5`$ we compared our results with those from Hambly et al. (1999). They used photographics plates from the United Kingdom Schmidt Telescope to construct a $`6\mathrm{°}\times 6\mathrm{°}`$ proper motion survey centered on the Pleiades. To minimize the contamination, we chose all the objects out of our photometric candidates having proper motion within $`1\sigma `$ ($`20`$ mas/yr) of the known cluster motion ($`\mu _\alpha =+19`$ mas/yr, $`\mu _\delta =43`$ mas/yr, Robichon et al. 1999). They are indicated as open triangles in Fig 5.

We stopped our stellar selection at $`I=17.5`$ corresponding about to Hambly’s survey completeness limit but also to the HBML. A short list of those very probable low mass stellar members in our survey is presented in Table 3. We considered that the residual contamination of this sample is low enough to be neglected. The analysis of the fainter objects, i.e. substellar candidates, has been done from the long exposure images and is explained hereafter.

### 3.2 Substellar domain

#### 3.2.1 Photometric selection

The $`(I,IZ)`$ color magnitude of the point-like objects contained in the long exposures of the 17 Pleiades fields is shown in Fig 6. Candidates previously identified by Bouvier et al. (1998) and confirmed on the basis of spectroscopic data, infrared photometry (Martín et al. 2000) and proper motion (Moraux et al. 2001) are shown as open circles. These objects define the high mass part of the cluster substellar sequence from $`I17.8`$ down to about $`I19.5`$. We note that the location of two Pleiades members (CFHT-PL-12 and CFHT-PL-16) suggest that they are likely binaries as already suspected by Bouvier et al. (1998) and Martín et al. (2000).

We overplotted the 120 Myr isochrones from the NEXTGEN and DUSTY models from Baraffe et al. (1998) and Chabrier et al. (2000) respectively, assuming a distance modulus for the Pleiades cluster of $`(mM)_o=5.53`$, $`A_V=0.12`$ and a solar metallicity. At a Teff which corresponds to late-M and early L spectral types, dust grains begin to form, changing the opacity and resulting in objects having bluer $`IZ`$ colors than predicted by the NEXTGEN models. DUSTY models instead include a treatment of dust grains in cool atmospheres for T$`{}_{eff}{}^{}`$ 2300K. To build our sample of Pleiades brown dwarf candidates for $`17.8I19`$, we defined a line 0.1 mag bluer in $`IZ`$ than the NEXTGEN isochrone and we selected all the sources located on the right side of this line. For $`I19`$, we selected all the objects redward of the DUSTY isochrone and we stopped our selection around the completeness limit, i.e. $`M0.03M_{}`$. All the candidates are shown in Fig 6 as filled triangles. The photometry and coordinates of these objects are given in Table 4. Two sources are located $``$0.12 mag left on the NEXTGEN isochrone at about $`I=18.5`$ and have not been considered as brown dwarf candidates. The proper motion of the faintest of these two objects has been measured and indicates non-membership. The other source will be followed up but this will not change the mass function estimate.

Part of our survey overlaps with Bouvier et al.’s (1998) survey performed in 1996, so that we were able to derive proper motion for some of the objects identified in both surveys. The two epochs of observations are separated by approximately 4 years and the resulting proper motion uncertainty is typically $`1\sigma 7`$mas/yr. The details of the procedures used to derive proper motion are given in Moraux et al. (2001). Objects which have a proper motion less than $`2\sigma `$ from the cluster motion ($`\mu _\alpha =+19`$ mas/yr, $`\mu _\delta =43`$ mas/yr) are very likely Pleiades members. We have written those objects in bold characters in Table 4 and objects whose proper motion indicates non membership in parenthesis.
Two of our new candidates had already been identified as probable Pleiades members by Pinfield et al. (2000), CFHT-PLIZ-2=BPL 327 and CFHT-PLIZ-10=BPL 240, on the basis of their optical and infrared photometry and from their proper motion by Hambly et al. (1999). These objects are also written in bold characters in Table 4.

#### 3.2.2 Contamination

An ($`I,IZ`$) diagram alone cannot identify objects as certain Pleiades members as one expects some level of contamination by field stars. Due to the relatively high galactic latitude of the Pleiades cluster, heavily reddened distant objects should not contaminate our photometric sample of candidate members. However, the relative location in the CMD of the theoretical Pleiades isochrone and a zero age main sequence isochrone from DUSTY models indicates that some of the photometrically selected brown dwarfs candidates could in fact be field M-dwarfs at a distance about 30% closer than the Pleiades. Considering the whole selection range which extends from 0.5 mag below the cluster sequence to the binary cluster sequence, we find that contaminating field M-dwarfs can lie in a distance range from 60 to 125 pc. Then, taking into account the area of the survey, the volume occupied by contaminants is about 1150 pc3. The field star luminosity function for $`M_I=1214.5`$ can be approximated as a constant $`\varphi 0.003`$ stars/pc3 per unit $`M_I`$ as estimated from the DENIS survey (Delfosse 1997). We therefore expect to find about 7 field stars out of 21 candidates in the range $`I=17.819.8`$, i.e. a contamination level of about 33% as previously derived from proper motion measurements by Moraux et al. (2001).

At fainter magnitudes the contamination level cannot be derived from the field star luminosity function which is not very well known for $`M_I>14.5`$. However, we can use the number of stars identified in the DENIS survey down to $`I=18`$ in a given color (or temperature) range in order to estimate the number of contaminants in our sample for this color interval. For example, the brown dwarf candidates with $`I`$ between 20.2 and 21.7 have a temperature of $`2000`$K (Chabrier et al. 2000). We then consider the number of DENIS objects within a restricted temperature range around this value and $`I`$ between 16.5 and 18, and multiply this number by two factors : a) the ratio of our CFHT survey area to the DENIS survey area, and b) the ratio of the two volumes corresponding to the two magnitude ranges ($`I=`$20.2 to 21.7 and $`I=`$16.5 to 18) for the DENIS survey. We thus predict 5 or 6 field dwarfs to occupy the region of the Pleiades color magnitude diagram corresponding to $`0.04M_{}M0.03M_{}`$. This indicates again a contamination level of $`30\%`$.

For redder objects, the statistics of the DENIS survey is very low so that it is difficult to estimate the contamination. Moreover, our survey starts to be incomplete in this domain and that is why we limited our analysis to $`M=0.03M_{}`$.

### 3.3 Radial distribution

In order to deduce the total number of Pleiades members and derive the cluster mass function from a survey covering only a fraction of the cluster area, we need to investigate the spatial distribution of cluster members and its dependence on mass.

The spatial distribution of our Pleiades brown dwarf candidates is shown in Fig. 2. Overplotted are circles of radii 0.75 to 3.5 degrees centered on the cluster center. From this diagram we estimated the covered area within annulii of $`0.25\mathrm{°}`$ width and we counted the number of substellar objects found therein. We then obtained radial surface densities for brown dwarf candidates by dividing these numbers by the corresponding surveyed areas. We proceeded in the same way for low mass stars. The number of stellar and substellar objects per square degree as a function of the radial distance is shown in Fig 7.

The stellar distribution ($`13.5<I<17.5`$, i.e. $`0.48>M>0.08M_{}`$) shown on the left panel is well fitted by a King distribution (King 1962) :

$$f(x)=k\left[\frac{1}{\sqrt{1+x}}\frac{1}{\sqrt{1+x_t}}\right]^2$$
(1)
where $`k`$ is a normalisation constant, $`x=(r/r_c)^2`$ and $`x_t=(r_t/r_c)^2`$ with r the radius from the cluster center. The core radius $`r_c`$ increases as the stellar mass decreases and the tidal radius $`r_t`$ corresponds to the location where the gravitationnal potential of the galaxy equals the cluster potential. Using $`r_t=5.54\mathrm{°}`$ from Pinfield et al. (1998), we found $`k110`$ per square degrees and $`r_c2`$ degrees for a median mass of our stellar sample of $`M0.2M_{}`$. From this radial distribution, we obtain a total of 557 stars between 0.08 and 0.48$`M_{}`$. For a dynamically relaxed cluster, the core radius is expected to vary with stellar mass as $`M^{0.5}`$ and Jameson et al. (2002) derived the relationship $`r_c=0.733M^{0.5}`$ for the Pleiades stars. For $`M0.2M_{}`$, this yields $`r_c1.6\mathrm{°}`$, slightly smaller than our value.

The radial distribution of brown dwarf candidates is shown on the right side of the figure 7. The plain histogram corresponds to the whole list of candidates whereas the shaded histogram corresponds to objects having proper motion consistent with cluster membership (written in bold characters in Table 4). This histogram does not extend further than $`r=2.25^o`$ corresponding to the surveys of Pinfield et al. (2000) and Bouvier et al. (1998) from which proper motions have been derived. Those surveys were not as deep as ours so that only brown dwarf candidates brighter than $`I=20.9`$ were counted. A King profile fitted to this histogram yields $`r_c1.3`$ degrees as a lower limit to the cluster substellar core radius. The plain histogram also decreases in the first few radius bins but then increases further away from the cluster center. Note, however, that the uncertainties due to small number statistics are large and the plain histogram is not corrected for contamination for field stars, whose rate is expected to increase away from the cluster center. An illustrative King profile with $`r_c3.0`$ degrees is shown as a possible fit to the brown dwarf distribution, mainly based on the few first radial bins. The median mass for the brown dwarfs candidates is $`M0.05M_{}`$ and the value expected by Jameson et al.’s (2002) relationship $`r_c=0.733M^{0.5}`$ is 3.4 degrees.

With $`r_c3.0^o`$ and $`k=28.5`$ per square degrees, integration of the King distribution yields a total of $`130`$ brown dwarfs between $`17.8<I<21.7`$, i.e. between 0.07 and 0.03$`M_{}`$ in the whole cluster. From this distribution, we expect to find $`26`$ brown dwarfs Pleiades members in our survey out of the 40 selected candidates. This would correspond to a contamination level of 35%, quite consistent with our estimate above.

### 3.4 The Pleiades mass function

The Pleiades mass function can be estimated from our CFHT large survey over a continuous mass range from $`0.03M_{}`$ to $`0.45M_{}`$. For the stellar part (down to $`I=17.5`$) we use our sample of candidates derived from short exposure images and decontaminated as explained above. We derived masses from $`I`$-band magnitudes using the 120 Myr isochrone from Baraffe et al. (1998). Below the HBML ($`I17.5`$) we consider our selection of brown dwarf candidates (Table 4) and we apply a correction factor of 0.7, assuming a contamination level of 30%. We used the 120 Myr isochrone from the DUSTY models of Chabrier et al. (2000) to estimate masses.

In order to correctly estimate the mass function of the whole cluster from a survey which is spatially uncomplete, one has to take into account the different radial distribution of low mass stars and brown dwarfs. Our CFHT fields are located between 0.75 and 3.5 degrees from the cluster center. We now proceed to estimate the fraction of low mass objects located in this ring compared to the total number of such objects in the cluster. For a King-profile surface density distribution, the total number of stars seen in projection within a distance $`r`$ of the cluster center is obtained by integrating equation 1 :

$$n(x)=k\pi r_c^2\left[\mathrm{ln}(1+x)4\frac{\sqrt{1+x}1}{\sqrt{1+x_t}}+\frac{x}{1+x_t}\right]$$
(2)
We consider stellar core radii following Jameson et al. (2002) relationship $`r_c=0.733M^{0.5}`$ and assumed a substellar core radius $`r_c=3.0^o`$ (see previous section). We then deduce from equation 2 that 74% of the $`0.4M_{}`$ stars and 80% of the brown dwarfs are located within a distance from the cluster center between 0.75 and 3.5 degrees. This indicates that the relative number of brown dwarfs compared to low mass stars deduced from our survey is representative to their relative number over the whole cluster. In other words, the area covered by this survey is large enough so that any correction to the mass function for a mass-dependent radial distribution is negligeable.

The derived mass function is shown in Fig. 8 as the number of objects per unit mass. Within the uncertainties, it is reasonably well-fitted by a single power-law $`dN/dMM^\alpha `$ over the mass range from $`0.03M_{}`$ to $`0.45M_{}`$. A possibility to explain why the $`0.06M_{}`$ data point is low is discussed in Dobbie et al. (2002). A linear regression through the data points yields an index of $`\alpha =0.60\pm 0.11`$, where the uncertainty is the $`1\sigma `$ fit error222Using the lower limit of $`r_c=1.3\mathrm{°}`$ for the cluster core radius in the susbtellar domain would yield $`\alpha =0.63\pm 0.11`$ instead.. This result is consistent with previous estimates (cf. Table 1) and is affected by smaller uncertainties thanks to the combination of relatively large samples of low mass stars and brown dwarfs, proper correction for contamination by fields stars, and extended radial coverage of the cluster.

In Moraux et al. (2001) the mass function index was found to be $`\alpha =0.51\pm 0.15`$. Lacking a proper determination of the radial distribution of cluster members, the assumption was made that brown dwarfs and very low mass stars were similarly distributed. The index estimate was based on Bouvier et al.’s (1998) survey which covered fields spread between 0.75 and 1.75 degrees from the cluster center. From the radial distribution derived above, we find that this annulus contains 43% of the $`0.4M_{}`$ stars and 33% of the cluster brown dwarfs. Applying these correcting factors to the number of Pleiades members found in that survey, the mass function index becomes 0.63 instead of 0.51. This corrected value is in excellent agreement with our new estimate.

Extrapolating the power-law mass function down to $`M=0.01M_{}`$, we predict a total number of $`270`$ brown dwarfs in the Pleiades for a total mass of about $`10M_{}`$. Clearly, while brown dwarfs are relatively numerous, they do not contribute significantly to the cluster mass. Adams et al. (2001) derived a total mass of $`800M_{}`$ for the Pleiades which means that, even though brown dwarfs account for about 25% of the cluster members, they represent less than 1.5% of the cluster mass.

## 4 Discussion

One of the main motivations for the determination of the lower mass function is to constrain the star and brown dwarf formation processes. A pressing issue is then whether the Pleiades mass function (MF) observed at an age of $`120`$ Myr is representative of the initial mass function (IMF), i.e. the mass spectrum resulting from the formation process.

We compare our results to the mass function of other young open clusters and star forming regions in order to investigate the dynamical evolution of the Pleiades cluster from a few Myr to its present age of 120 Myr. We also compare the Pleiades MF to the Galactic disk mass function in order to constrain the brown dwarf formation process.

### 4.1 Comparison to other young clusters : clues to early dynamical evolution

We compare our results to those obtained recently for other young open clusters such as M35 ($`150`$ Myr, Barrado et al. in prep) and $`\alpha `$ Per ($`80`$ Myr, Barrado et al. 2002). The lower mass function of those clusters can be approximated by a single power-law with an index $`\alpha =0.58`$ between $`0.4`$ and $`0.1M_{}`$ for M35 and $`\alpha =0.56`$ between $`0.2`$ and $`0.06M_{}`$ for $`\alpha `$ Per. These values are very similar to the one determined here for the Pleiades ($`\alpha =0.60\pm 0.11`$). Results obtained for star forming regions such as $`\sigma `$-Orionis ($`5`$ Myr, Bejar et al. 2001) and IC348 ($`3`$ Myr, Tej et al. 2002) are also consistent with this value. The power-law indices are $`\alpha =0.8\pm 0.4`$ ($`0.2>M>0.013M_{}`$) for $`\sigma `$-Orionis and $`\alpha =0.7\pm 0.2`$ ($`0.5>M>0.035M_{}`$) for IC348. The mass functions of Pleiades-age clusters and star forming regions thus appear to have a similar shape across the stellar-substellar boundary.

Comparison over a larger mass domain is best achieved by plotting the number of stars per logarithmic mass units :

$$\xi (\mathrm{log}M)=\frac{dN}{d\mathrm{log}M}$$

(in this representation, $`dN/dMM^{1.0}`$ corresponds to $`dN/d\mathrm{log}M`$ constant). Fig. 9 shows the Pleiades mass function over the mass range from 0.03$`M_{}`$ to 10$`M_{}`$. The stellar portion of the mass function has been derived from the Prosser and Stauffer Open Cluster Database333Available at http://cfa-www.harvard.edu/stauffer/opencl/ constructed from several proper motion surveys. This catalog, however, becomes uncomplete below $`0.5M_{}`$. We therefore normalized our survey so that the total number of objects over the 0.08-0.48$`M_{}`$ mass range corresponds to the number computed from the King profile (see section 3.3). This normalization ensures the continuity of the Pleiades mass function shown in Fig. 9 from the brown dwarfs up to the most massive stars of the cluster.

The cluster mass function is fitted by a lognormal form over more than 2 decades in mass :

$$\xi (\mathrm{log}M)\mathrm{exp}\left[\frac{(\mathrm{log}M\mathrm{log}M)^2}{2\sigma _{logM}^2}\right]$$
(3)
with $`M0.25M_{}`$ and $`\sigma _{logM}=0.52`$. This mass function can be approximated by a single power-law above $`1.5M_{}`$, $`\xi (\mathrm{log}M)M^{1.7}`$ in logarithmic mass units (or $`dN/dMM^{2.7}`$ in linear mass units), and peaks at $`M0.25M_{}`$ before decreasing at lower mass. Luhman et al. (2000) found that $`\rho `$ Ophiuchus, IC348 and the Trapezium mass functions $`\xi (\mathrm{log}M)`$ also exhibit a maximum around $`0.25M_{}`$ and are quite similar to the Pleiades one in the stellar domain.

Within uncertainties, we thus find no evidence for significant differences in the mass function of the Pleiades, other Pleiades-age open clusters and star forming regions across the stellar-substellar boundary but also over the entire stellar domain. This suggests that the cluster population observed at an age of about 120 Myr is still quite similar to its population at only a few Myr. As a cluster evolves weak gravitational encounters occur leading eventually to the evaporation of the lowest mass members. The results above seem to indicate that such a dynamical evolution has not yet affected the mass function of the cluster. This is in agreement with dynamical models (e.g. de la Fuente Marcos & de la Fuente Marcos 2000, Adams et al. 2002) which predict that only about 10% of the total number of cluster members is lost after $`100`$ Myr, and that the fraction of brown dwarfs to stars remains nearly constant, so that the shape of the mass function across the stellar-substellar boundary is hardly affected.

The mass function of the Pleiades at an age of $`120`$ Myr thus seems to be representative of the cluster population at an early age of a few Myr. But does it also reflects the initial population of the cluster at the time it formed (IMF)? Rapid dynamical processes associated to cluster formation might conceivably lead to the prompt ejection of low mass stars and brown dwarfs during the very early stage of a cluster life, before a few Myr. If instrumental, such processes might result in a depletion of the low mass cluster population very early on and thus modify the shape of the mass function.

An example of conditions under which such processes may occur is when the gas is expelled from the cluster during its formation. The gravitational potential decreases drastically and all the low mass objects in the outer part of the cluster are ejected. However, current models (Kroupa 2001) predict that roughly as many massive stars as low mass stars are ejected in this process, unless mass segregation is present initially. Hence, the shape of the mass function should not be significantly affected.

Another violent dynamical process may be associated to the formation of brown dwarfs themselves. Reipurth & Clarke (2001) proposed that substellar objects are “stellar embryos” according to the following scenario : as molecular cloud cores fragment to form unstable protostellar multiple systems which decay dynamically, the lowest mass fragments are ejected from the core, and deprived of surrounding gas to accrete they remain substellar objects. In this scenario, the velocity dispersion may be expected to be larger for brown dwarfs than for stars and, for a fraction of substellar objects, the velocity may indeed exceed the escape velocity. Such a process could then be very efficient in quickly removing the lowest mass objects from the cluster, thus strongly modifying its initial mass function. A possibility to constrain this brown dwarf formation process would be to search for an observational signature of their primordial kinematics. Unfortunately, at an age of 120 Myr the Pleiades cluster is already largely relaxed and has lost memory of its initial evolution.

Furthermore, while several models have been developped to investigate the dynamical ejection of brown dwarfs, predictions regarding their initial velocity distribution differ. Sterzik & Durisen (1998) find that the velocity dispersion depends on both mass and binarity, being larger for low mass single objects than for massive binaries. Delgado-Donate et al. (2002) find that it depends only on binarity and that nearly all ejected objects are preferentially low mass and single ones. According to M. Bate (priv. comm.) the velocity dispersion depends neither on mass nor on binarity so that dynamical ejection of low mass objects should not affect the overall mass function of the cluster. Pending the resolution of these uncertainties, it is then difficult to assess whether the proposed ejection mechanism would lead to peculiar kinematical signatures in the substellar population and would thus possibly affect the cluster mass function. This issue is further discussed in the next section by comparing the Pleiades and the Galactic disk mass functions.

### 4.2 The Pleiades and Galactic disk mass functions : clues to the brown dwarf formation process

Below $`0.8M_{}`$ field stars have not had time to evolve off the main sequence, so that the observed Galactic disk mass function ought to be representative of the initial mass function. A log-normal fit to the Galactic disk mass function in the mass range 0.1-1.0$`M_{}`$ has been derived by Chabrier (2001) based on large scale photometric surveys of low mass stars in the solar neighbourhood. This log-normal approximation to the field mass function is shown as a dashed line in Fig. 9 where it has been normalized so as to have the same number of $`1M_{}`$ stars as in the Pleiades. Major differences clearly appear between the Pleiades and the field mass functions at low masses. The log normal mass function peaks at $`M0.1M_{}`$ for the disk population and at $`0.25M_{}`$ for Pleiades members. Furthermore, Chabrier (2002) estimates that the brown dwarf population in the disk is comparable in number to the stellar one, $`N_{BD}N_{}`$, and that the substellar mass contribution to the disk budget amounts to about $`10\%`$. In the Pleiades, we have instead $`N_{BD}N_{}/3`$ and a brown dwarf mass contribution to the cluster mass of only $`1.5\%`$. Indeed, the Pleiades mass function lies well below the disk’s one at the stellar-substellar boundary in Fig. 9.

Part of the observed difference may arise from the effect of binarity, which is not accounted for in the Pleiades mass function. While the Pleiades mass function derived above includes unresolved cluster binaries, the Galactic disk mass function has been derived for the single star population (i.e. all binaries are resolved, cf. Chabrier 2001). In order to correct the observed Pleiades mass function for unresolved binaries at low masses, we follow Luhman et al. (1998) who used a Monte Carlo technique to estimate the difference between single star and system mass functions. Assuming that the properties of field binaries derived by Duquennoy & Mayor (1991) apply to Pleiades systems (cf. Bouvier, Rigaut & Nadeau 1997), we form unresolved binaries by randomly pairing objects drawn from a segmented power law mass function which represents the single star mass distribution down to the stellar-substellar boundary444We do not attempt to correct the observed Pleiades mass function for unresolved binaries in the substellar domain because the brown dwarf binary statistics is still very poorly known, as is the statistics of brown dwarfs companions to low mass stars. For an assumed 50% binary fraction in the Pleiades cluster, we find no major differences between the star and system mass functions at $`M>0.6M_{}`$ since above this mass the presence of a low mass companion does not significantly affect the determination of the primary’s mass. Below $`0.6M_{}`$ down to $`0.07M_{}`$, however, we find $`\mathrm{\Delta }\alpha 0.5`$ between the power law exponents of the single star and system mass functions.

Then, approximating the log-normal form of the Pleiades mass function derived above by a three-segment power-law between $`0.07M_{}`$ and $`0.6M_{}`$, we estimate the binary corrected Pleiades mass function which is shown in Fig. 9 (dotted line). The binary corrected Pleiades mass function and the Galactic disk mass function now both peak around $`0.130.1M_{}`$ and their shapes are roughly identical in the stellar domain down to $`0.07M_{}`$. This suggests that the difference between the observed Pleiades and Galactic disk mass functions in the stellar range is merely the result of unresolved Pleiades binaries. When binarity is properly accounted for, the stellar mass function of the Pleiades is consistent with that of field stars of the solar neighbourhood. The comparison of the two mass functions in the substellar domain cannot be as detailed, due to the large uncertainties still affecting the derivation of the substellar mass function in the Galactic disk. Chabrier (2002) derives an upper limit of $`\alpha 1`$ for a power law approximation of the field substellar mass function (see also Reid et al. 1999). This upper limit is consistent with the Pleiades substellar mass function exponent $`\alpha 0.6`$ derived above.

The similarity of the Pleiades and field mass functions down to at least the substellar limit and possibly below suggests that low mass objects have remained in the cluster at the time of its formation. Hence, we do not find evidence for massive ejection of the lowest mass objects early in the life of the cluster, which suggests either that dynamical ejection from protostellar groups is not the dominant mode of brown dwarf formation, or that this process yields a velocity dispersion for brown dwarfs which is not different from that of stars.

## 5 Conclusion

We have conducted a deep wide field photometric survey of the Pleiades cluster to build a sample of probable cluster members with masses in the range $`0.03M_{}`$ to $`0.48M_{}`$. We have identified 40 brown dwarfs candidates, of which 29 are new discoveries. Taking into account the radial distribution of cluster members, we derive the cluster mass function accross the stellar-substellar boundary. We find that a single power-law $`dN/dMM^\alpha `$ with an index $`\alpha =0.60\pm 0.11`$ provides a good match to the cluster mass function in the $`0.030.48M_{}`$ range. This new estimate is based on a survey which combines a large radial coverage of the cluster and a realistic assessment of the contamination by field stars. Furthermore, the survey completely covers the $`0.030.48M_{}`$ mass range, so that the result does not rely on the combination of heterogeneous surveys, as has been the case before. We therefore believe this new estimate is reasonably robust. Small changes may be expected when our survey has been followed up with either infrared photometry and/or proper motions.

Over a larger mass domain, covering almost 3 decades in masses from $`0.03M_{}`$ to $`10M_{}`$, we find that the cluster mass function is better fitted by a log-normal distribution with $`M0.25M_{}`$ and $`\sigma _{logM}0.52`$. When unresolved Pleiades binaries are taken into account, the log-normal Pleiades mass function is not unlike the Galactic disk mass function. This suggests that the dynamical evolution of the cluster has had yet little effect on its mass content at an age of 120 Myr. It also suggests that the brown dwarf formation process does not lead to the dynamical evaporation of substellar objects at the time the cluster forms.

###### Acknowledgements.

We thank E. Bertin for allowing us access to PSFex before public release, E. Magnier for his help in the astrometric calibration of the frames, J. Adams and N. Hambly for providing us data in electronic form prior to publication, I. Baraffe for computing specific substellar isochrones for us and K. Luhman for providing his Monte-Carlo software program to estimate the effect of unresolved binaries on the mass function. We also gratefully aknowledge helpful discussions with X.Delfosse on estimating field star contamination from DENIS data, with M. Bate, C. Clarke, E. Delgado and M. Sterzik on the dynamical evolution of young brown dwarfs in clusters, and with G. Chabrier who brought our attention on the Galactic disk mass function and on the effect unresolved binary systems might have on the shape of the observed mass function.