Maria S Contreras

Theodorus De Bruin

Pascal Mougin

Hervé Toulhoat

Thermochemistry of 1-methylnaphthalene hydroconversion : comparison of group contribution and ab-initio models

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

To satisfy the increasing global energy demand, it has become mandatory to make better use of heavy crude oils, vacuum residues, bitumen, and asphalts, by cracking these oils and/or cuts into lighter fractions such as gasoline, diesel, and kerosene. [START_REF]From Molecular Theory to Industrial Application[END_REF] Hydrocracking combines the cracking process together with the removal of sulfur, nitrogen impurities, and the saturation of aromatic and unsaturated bonds. To gain more insight in this process, numerous experimental and theoretical studies have been conducted. 2 However, all these studies are hampered by the large variety of complex chemical structures. Despite a vast body of work on characterization techniques and experimental information on molecular structures of heavy petroleum fractions, the current analytical techniques do not allow characterization of individual molecules or isomer groups in heavy oils. The widely used SARA method of characterization classifies the components of the oil into four big families: saturates, aromatics, resins, and asphaltenes. Aromatics, resins, and asphaltene fractions contain one or several nuclei, or cores, composed of aromatic and/or naphthenic cycles, onto which aliphatic chains are grafted. Since these molecules are very complex, a detailed molecular description is generally not available. Molecular structures of reconstructed aromatic, asphaltene and resin molecules have therefore been proposed. 3 The hydroconversion process is generally carried out at high temperatures (320-380 °C) and H 2 pressures up to 150 bar. These reaction conditions are the result of a compromise between the thermodynamics and kinetics: the thermodynamic equilibrium shifts at lower temperatures to the right in favor of the hydrogenated products, while the reactions become intrinsically faster at higher temperatures. Consequently, to find the optimum reaction conditions, a detailed knowledge of the thermodynamics and kinetics of the elementary reactions for the hydrogenation is mandatory. Once these properties have been determined, microkinetic models can be developed that predict the final product composition. The identification of all the essential elementary reactions for asphaltene or resin type molecules (or representative models thereof) hydrogenation is an almost impossible task. Instead, we have selected a welldefined model compound, 1-methylnaphtalene (1MN), for which Scheme 1 depicts the most essential hydrogenation reactions. In this scheme, also two cracking reactions are presented: the demethylation of 1MN and the ring-opening reactions of the aliphatic cycle.

In order to model the system as shown in Scheme 1, reaction rates for all direct and reverse reactions must be represented as functions of temperature and H 2 pressure. Since the ratios of direct to reverse rate constants are equal to the equilibrium constants, a useful constraint will be provided by the knowledge of the latter. To that end, it is mandatory to calculate the Gibbs energies of reactants and products. Since experimental data are available only for a limited number of molecules (Table 1) in the literature, we turn to computational methods to complete this data set.

A large body of work has been performed over the past four decades to estimate thermochemical properties of organic compounds. For example, the group additivity (GA) method developed by Benson et al. 4 allows estimating thermochemical properties like the heat of formation, heat capacity, and entropy, but their accuracy is not guaranteed for any possible structure. Benson's method does not incorporate stabilizing nonadditive effects of delocalization or resonance, and Herndon et al. 5 concluded that including these parameters improves the accuracy of Δ f H°predictions for polycyclic aromatic hydrocarbons (PAH) with three or more aromatic rings. A method for estimating ideal gas heat capacities has been proposed recently, 6 which makes use of a single parameter, the ratio between the number of atoms in the molecule and its molar mass: since in the present work we deal with isomers of compounds having very similar molar masses and compositions, we have preferred to use the older Benson method with structural contributions up to the third order, as a reference calculation to be compared with quantum chemical models.

Quantum chemical models coupled to statistical mechanics can approach experimental accuracies for molecular properties. However, polycyclic aromatics exhibit strongly correlated and extended π electrons systems, which are challenging for theory. Post Hartree-Fock methods are very accurate but computationally expensive.

Density Functional Theory (DFT) based methods on the other hand are computationally cheaper, but some functionals (like the popular B3LYP 7 ) give a serious overestimation of the enthalpies of formation with respect to the experimental values. 8 A systematic increase in the deviation with an increasing number of carbon atoms in a molecule has been reported by Van Speybroeck et al. at the B3LYP6-31G(d)/ level. 9 Composite methods (e.g., Gn-methods, G3B3, CBS-QB3) combine the results of several calculations, allowing more accurate estimation of heats of formation than DFT based methods.

In particular, the G4 theory has proven to be one of the most advanced methods to obtain accurate absolute energies. 10 However, we are interested in reaction energies and not absolute energies. In order to compare the precision of ab initio methods with the group-based method of Benson, we have calculated the heats of formation of all compounds and compared them to experimental values when possible. We will present different methods to compute the heat of formation, using atomization energies of iso-and homodesmotic reactions. Isodesmic reactions conserve the number of types of bonds, while in homodesmotic reactions the hybridization of atoms in bonds is also conserved. [START_REF] Mcnaught | Compendium of Chemical Terminology[END_REF] This type of reaction scheme aids in rescinding some of the systematic errors that arise in quantum mechanics, due to an incomplete capture of electronic correlation. [START_REF] Yu | [END_REF]13 However, since different isodesmic/homodesmotic reactions can be set up for the molecule of interest, the calculated thermodynamic property may be reaction dependent, as we will see later on.

In the next section, we detail the applied computational methods. We then discuss the obtained equilibrium constants from the computations and compare them with available experimental data where possible. To test the viability of our Scheme 1. Simplified Reaction Pathway for 1-Methylnaphthalene Hydroconversion a a Hydrocracking reactions are omitted except for ring opening. Naphthalene is also produced by demethylation of 1-methylnaphthalene. Reaction numbers are in parentheses; compounds are numbered in bold. model, i.e., the completeness of the reaction scheme presented in Scheme 1 and the accuracy of the computed data, we use the calculated thermodynamic data from Flash calculations in order to predict the final product composition under experimental reaction conditions and compare them with experiments by Bouchy et al., who studied the hydrogenation and hydrocracking of a model light cycle oil feed. 2

METHODS

The experimental data used in this study were obtained from the DIPPR database [START_REF] Rowley | DIPPR Data Compilation of Pure Compound Properties; Design Institute for Physical Property Data[END_REF] and NIST database. [START_REF] Afeefy | Phase Thermochemical Data in NIST Chemistry WebBook[END_REF] Benson's GA method was used to estimate the heats of formation (Δ f H°), entropy (S°), and heat capacity (C p °) of all the molecules presented in Scheme 1.

For DFT calculations, the Jaguar package [START_REF] Linstrom | National Institute of Standards and Technology[END_REF] from Schrodinger was used, and the input structures were generated using Maestro. [START_REF] Maestro | version 9.3[END_REF] Geometry optimization calculations were performed at the B3LYP/6-31G(d,p) level with frequency calculations that use the harmonic oscillator approach. The resulting geometries at the B3LYP level were used as input structures for the Gaussian-4 (G4) [START_REF] Curtiss | [END_REF] composite method calculations with Gaussian 09. [START_REF] Frisch | Revision B.01[END_REF] An atomization approach (eq 1) can be used to calculate the heat of formation and entropy of C n H m from DFT calculations, as shown in the eqs 2 and 3. Equation 2 is visualized in Figure 1.

↔+ nm CH C H nm (1) Δ °= Δ °+Δ°- °+ °- °Hn H m H nH mH H (C H ) (C ) (H ) [( C ) ( H ) (C H )] nm nm f ,298K f exp 298K f exp 298K calc 298K calc 298K calc ,298K
(2)

°= °+ °- °+ °- °S nS mS nS mS S (C H ) (C ) (H ) [( C ) ( H ) (C H )] nm nm ,298K exp 298K exp 298K calc 298K calc 298K calc ,298K (3)
Important systematic errors in B3LYP and G4 calculations have been reported that result in deviations observed with respect to experimental enthalpies and entropies of carbon and hydrogen atoms. An overestimation of 30 kJ•mol -1 with respect to the experimental values of the heats of formation was observed by Saeys et al. 8 at the B3LYP/6-31G(d) level. The underestimation of atomization energies and the overestimation of enthalpies of formation using the B3LYP functional are due to the underestimation of the bond energies. 9 Therefore, a different approach was used in parallel. In this approach, we establish the "real" reactions for the enthalpy of formation by taking experimental values, and eq 1-3 were modified as follows:

↔+ n m CH C 2 H nm 2 (4)
Δ °= Δ °+Δ°- °+ °- °⎡ ⎣ ⎢ ⎤ ⎦ ⎥ Hn H m H nH m H H (C H ) (C ) 2 (H ) (C ) 2 (H ) (C H ) nm nm f ,298K f exp 298K f exp 2,298K calc 298K calc 2,298K calc ,298K (5)
°= °+ °- °+ °- °⎡ ⎣ ⎢ ⎤ ⎦ ⎥ S nS m S nS m S S (C H ) (C ) 2 (H ) (C ) 2 (H ) (C H ) nm nm ,298K exp 298K exp 2,298K calc 298K calc 2,298K calc ,298K (6)
For this approach, the following experimental data were used for carbon at 298.15 K: Δ f H°= 716.7 kJ•mol -1 and S°= 158.1 J•mol -1 •K -1 . [START_REF] Afeefy | Phase Thermochemical Data in NIST Chemistry WebBook[END_REF] The experimental data for H 2 are shown in Table 1. For comparison, the heat of formation and entropy were also calculated using the atomization approach but values are not shown here; parity plots for this approach are provided in the Supporting Information (SI), Figures S1 andS2.

Flash calculations were carried out taking into account liquid-vapor equilibrium, which means that it is necessary to Table 1. Experimental Data Available for S°, Δ f H°, and calculate the equilibrium constants for each reaction. These constants can be calculated as follows:

Δ f G°a t 298 K a molecule b S°(J•mol -1 •K -1 ) Δ f H°(kJ•mol -1 ) Δ f G°(kJ•mol -1 )
∏ °= -Δ °= ̂υ ⎛ ⎝ ⎜ ⎜ ⎞ ⎠ ⎟ ⎟ KT GT RT a () e x p () () j j i i p, R i j , (7)
where K p,j °is the equilibrium constant, R is the universal ideal gas constant, T is the temperature, aî is the activity for the component i, ν is the stoichiometry number, j represents the reaction number, and Δ R G j °(T) is the change in the Gibbs energy for a reaction.

In order to calculate Δ R G°(T), we used the thermodynamic data calculated by the methods presented above. In the case of Benson's GA method, an extrapolation was made in order to calculate the thermodynamic properties at higher temperatures (above 298 K), using eqs 8 and 9.

∫ Δ °= Δ °+ °HT H C T T () ( 2 9 8 K ) () d T 298K p (8) ∫ °= °+ °S TS CT T T ( ) (298 K) () d T 298K p ( 9
)
where C p °(T) is the polynomial form of constant pressure heat capacity. Δ R G°(T) was calculated using the enthalpy and entropy as follows:

Δ °= Δ °-°GT HT T ST () () () R (10)
In the case of ab initio models, DFT calculations were carried out directly at the temperature desired.

Activities at equilibrium in the liquid and gas phase were calculated using the Peng-Robinson equation of state. If component i is present in the liquid phase, the activity is represented as

ϕ ̂= °a x PT P x P (, , ) ii i l i (11)
But, if component i is in the gas phase, the activity is calculated as ϕ ̂= °a y PT P y P (, ,)

i i i v i ( 12
)
where P°is the pressure under standard conditions, i.e., 1 bar, and ϕ i l and ϕ i v are the fugacity coefficients in the liquid and gas phase, respectively.

The parameters of the Peng-Robinson equation of state are expressed in terms of the critical properties and the acentric factor; this information is reported in the 2011 version of the DIPPR 801 database (BYU DIPPR Lab). [START_REF] Rowley | DIPPR Data Compilation of Pure Compound Properties; Design Institute for Physical Property Data[END_REF] The critical properties not reported in this database can be obtained by the group contribution model proposed by Marrero and Gani. [START_REF] Marrero | [END_REF] The acentric factor (ω) was estimated by using the method proposed by Constantinou et al. 21 The interaction parameters k ij for hydrogen and hydrocarbon binaries were calculated by using the method proposed by Moysan and coworkers. 22 Fugacity coefficients depend on the final composition of the system, which is initially unknown, and a trial-and-error iterative method is required to calculate fugacity coefficients as well as the final composition. The criteria of mass balance by elements must be fulfilled. The SIMPLEX method 23 is used to minimize the error between the thermodynamic equilibrium constants calculated using the Gibbs energy and those calculated in terms of activities (calculated using the equation of state). In order to achieve the calculation coupling phase and chemical equilibria, a first free energy minimization loop under mass balances constraints assuming a single gas phase is implemented; then a flash calculation is achieved and confirmed by a stability analysis in order to get a first guess of the amount of liquid phase if any; a second free energy minimization loop under constraints taking the gas an liquid phases into account allows one to predict the composition of both phases at equilibrium. Details of this algorithm have been published by Mougin et al. 24

RESULTS AND DISCUSSION

The gas phase standard state entropies and heats of formation (S°and Δ f H°) have been calculated using Benson's GA method and QM methods at the B3LYP and G4 levels for all the compounds present in Table 1, and a comparison of the calculated values and experimental data is represented in Figure 2 (for Δ f H°) and Figure 3 (for S°). See also Tables S1 andS2 in the Supporting Information, which also give examples of calculations using Benson's GA method.

In the case of the standard heat of formation, the linear regression for these models (Benson, B3LYP, and G4) shows an overestimation of Δ f H°by B3LYP, while a good (R 2 = 0.990) to excellent (R 2 = 0.999) agreement is found for, respectively, G4 and Benson.

Benson's method also calculates the entropy more accurately than B3LYP or G4 (Figure 3). For the QM calculations, the entropy contributions are calculated from the partition functions that use frequencies that were evaluated from the relaxed structure obtained at the B3LYP/6-31G(d,p) level and in the case of the G4 calculations at the B3LYP/6-31G(2df,p) level. The small difference between the two methods explains why the results for the entropy are that similar. The harmonic oscillator approximation is used in the assessment of the vibrational frequencies. This approximation might no longer be sufficiently precise for molecules with a large degree of flexibility and a high number of rotatable groups. 9 The shape of the rotational energy curve can give rise to various stable conformers for the investigated molecules. Furthermore, the presence of rotatable groups in the molecule can cause the wrong assignment of the symmetry number of the molecule, thereby introducing serious errors in the evaluation of the rotational entropy.

A remarkable improvement to calculate the Gibbs energy of formation and absolute entropies can be achieved by the utilization of isodesmic and homodesmotic reactions. However, the final results will depend on the choice of the homodesmotic reaction. This is, for example, illustrated in Table 2, where five different isodesmic reactions are given (of which four are homodesmotic) to calculate the thermodynamic properties of tetralin. Here, the total energy of tetralin is calculated at the G4 level, while experimental heats (or Gibbs energy) of formation and absolute entropies are used for the small molecules, e.g., methane, ethane, and cyclohexane. We refer to hybridimproved-G4-iso (HI-G4-iso) if the absolute energies come from G4 calculations and hybrid-improved-B3LYP-iso (HI-B3LYP-iso) if they are calculated at B3LYP/6-31G (d,p).

One way to describe isodesmic reactions in a systematic way, is the "bond separation reaction" in which all formal bonds between non-hydrogen atoms are separated into the simplest parent molecules containing these same kinds of bonds. All other bonds to the desired atom pair are replaced by hydrogen atoms, and the chemical equation is balanced by including the required number of single heavy atom hydrides, i.e., methane for hydrocarbons. For example, tetralin would be described as reaction a in Table 2. Yet, this method suffers from the fact that a relatively small error in the heat of formation of methane can produce a significant error in the heat of formation of 1methyltetraline, due to a large number of methane molecules in equation a. 25 It thus turns out that this reaction relatively poorly describes the thermodynamical properties of tetralin.

The thermodynamic properties are better reproduced with homodesmotic reactions, e.g., reactions b-e in Table 2. Generally, the use of homodesmotic reactions improves the accuracy of the calculated properties, as observed as well by Sivaramakrishnan et al. 26 It was also noticed that the conservation of naphthenic rings ameliorates the precision of the calculated thermodynamic properties. We also noticed that the results for entropy were further improved if one olefin bond was introduced in the naphthenic cycle to introduce some kind of rigidity, if this ring is fused to an aromatic ring in the reactant, as for example reaction e in Table 2.

Figures 4 and5 show parity plots for respectively Δ f H°and S°where the calculated values (Benson, HI-B3LYP-iso, and HI-G4-iso) are compared to experimental data. The use of homodesmotic reactions has significantly improved the accuracy of the calculated properties. This is also illustrated in Tables 3 and4, where the absolute average relative deviation (AARE) and the root-mean-square error (RMSE) have been reported for Δ f H°and S°. From these tables, it follows that Benson's GA empirical model has the smallest deviation for Δ f H°and S°, while among quantum chemical models, HI-G4iso reproduces better both Δ f H°and S°experimental values. Additionally, HI-G4-iso is completely universal, while Benson's method is not.

The equilibrium constant for a given process or reaction can be expressed in terms of Gibbs energy (see eq 7). This property (Δ R G 673K °) was calculated by using Benson's model and the QM models (G4, B3LYP, HI-B3LYP-iso, and HI-G4-iso) for each reaction present in Scheme 1. Figure 6 shows a comparison between those models (Benson, HI-G4-iso, and HI-B3LYP-iso models) and experimental data available, and Table 5 shows the deviation of calculated values from experimental data.

A smaller dispersion for Δ R G 673K °was observed for Benson's model (confirmed by a slope close to 1 for Benson's model), but errors of the same magnitude were obtained by using the HI-G4-iso model. Significant improvement was obtained by using DFT models with isodesmic reactions, especially at the G4 level, where the RMSE is reduced from 9.4 to 6.4.

The smallest deviations for Δ f H°and S°are obtained for Benson's model. However, remember that Benson'sG A method cannot always distinguish between position isomers observed in the hydroconversion process. Nonetheless, HI-G4iso could allow filling this gap, to estimate the equilibrium constant for isomerization reactions if experimental data are not available for position isomers.

Product Composition Prediction. In order to validate the methods studied, Benson's model and the HI-G4-iso method were chosen to estimate the equilibrium constant for the reactions represented in Scheme 1. Flash calculations were made to determine the equilibrium state for 1-methylnaphthalene hydroconversion at 40 bar and 400 °C.

Experimental results obtained by Bouchy et al. 27 under these conditions (40 bar and 400 °C) suggest a thermodynamic equilibrium for the (direct) hydrogenation reactions of 1MN. The discrepancy between the experimental 2MN/1MN ratio of about 0.17 and the theoretical value of 1.5 at thermodynamic equilibrium clearly demonstrate that isomerization reaction ( 5) is kinetically limited. Since the catalytic system does not contain any (strong) acid sites to catalyze this reaction, this reaction becomes kinetically prohibited, thereby explaining the much lower ratio experimentally observed.

However, in the flash calculations, which assume thermodynamic equilibrium for all reactions, the products related to this isomerization reaction (as well as other isomerization reactions and cracking reaction that are acid-catalyzed) thus will be overestimated with respect to the experimental values. This is shown in Table 6, in which hydrogenation and conversion were calculated as follows:

=--XX conversion (%) 100( 1)

26 (14) To remedy this discrepancy, we have simplified the reaction pathway (Scheme 1) by omitting the isomerization and cracking reactions, focusing exclusively on hydrogenation of the first and second rings of 1MN (reactions 1-4, as shown in Figure 7). By taking into account only the hydrogenation reactions, it follows that Benson's GA method better reproduces the experimental data of Bouchy than HI-G4-iso. The poorer agreement of the HI-G4-iso method can be explained by the differences observed for equilibrium constants values. Table 8 shows the equilibrium constants calculated at 673 K for Benson's and HI-G4-iso methods. Bouchy et al. 27 reported an equilibrium constant of 9.4 × 10 -4 for reaction 1 and a value of 8.2 × 10 -6 for reaction 2. This overestimation (or underestimation) of equilibrium constants is also observed by the overestimation of methyldecalin composition calculated using the HI-G4-iso method.

The integration of kinetic data in a microkinetic model that would use the rate constants for all essential elementary reactions would be a more elegant way to describe these experiments. These kinetic data can be obtained from experiments or by quantum mechanical calculations. The ratios of direct and inverse kinetic constants should then be equal to equilibrium constants in order to satisfy the microreversibility principle. Therefore, the present work is a necessary first step allowing one to introduce consistently such constraints in the microkinetic model. We have started to explore if the rate constants from full DFT calculations and those that are incorporated in a microkinetic model are sufficiently precise to describe these experiments. Once we have identified the reactions that are most rate-limiting, we can make suggestions to fine-tune this reaction to the desired products. This will be the object of a future publication.

CONCLUSIONS

We report an assessment of various DFT models for the calculation of molecular thermochemical properties in comparison with Benson's group contribution method for reactants, intermediates, and products involved in the hydrogenation of 1methylnaphatalene. Parity plots are compared with reference to available experimental data.

For this particular set of structures, Benson's GA method affords consistently the smallest deviations with respect to experimental data, reflecting the inherent correlative and interpolative nature of this method. Its prediction upon extrapolation to ill documented classes of structures will however remain questionable.

Among QM models, we have found HI-G4-iso, which combines evaluations of enthalpies and entropies of formation through improved homodesmotic schemes, as the most reliable, combining minimal RSME and bias.

Although small deviations are observed for Δ f H°and S°for Benson's and HI-G4-iso methods, discrepancies were observed between experimental product composition and those obtained by a simple thermodynamic model. As mentioned above, depending on the catalyst properties, the observed distribution of products may be far from full thermodynamic equilibrium. It is well-known for instance that the completion of hydrocracking, hydroisomerization, and hydrodecyclization reactions requires strong acid sites, together with sites able to activate hydrogenation and dehydrogenation reactions. Since acid sites were lacking in the experiments to which we referred, data reconciliation was possible assuming a scheme strictly limited to the latter reactions. For the interpretation of experimental results in the general case, a comprehensive microkinetic model will be needed, including consistently the thermodynamics constraints. Our work was a first step toward the goal of establishing such a model from first principles.

ASSOCIATED CONTENT

Tables showing a compilation of available experimental data and estimated heats of formation (Table S1) and entropies (Table S2) for the set of molecules used in this study; a compilation of the estimated values using the atomization approach of standard heat of formation (Figure S1) and standard entropy (Figure S2); tables showing a compilation of experimental data available, isodesmic reactions, and estimated values at the HI-G4-iso level for the set molecules used in this study (Table S3); and estimated values at 673 K of Gibbs energy (Table S4) and equilibrium constants (Table S5) of reactions shown in Scheme 1. This material is available free of charge via the Internet at http://pubs.acs.org.

Figure 1 .

 1 Figure 1. Illustration of the calculation of standard heats of formation from ab initio standard heats of atomization.

Figure 2 .

 2 Figure 2. Parity plot for Δ f H°:( ⧫ ) Benson, ( ■ )G 4 ,( •) B3LYP. Linear regression (y = mx + b): Benson: m = 1.004, R 2 = 0.999. B3LYP: m = 0.892, b = 74.962, R 2 = 0.962. G4: m = 0.939, b = 11.743, R 2 = 0.993.

Figure 3 .

 3 Figure 3. Parity plot for standard entropy (S°): ( ⧫ ) Benson, ( ■ )G 4 , (•) B3LYP. Linear regression (y = mx + b): Benson: m = 0.994, R 2 = 0.998. B3LYP: m = 1.148, b = 22.182, R 2 = 0.954. G4: m = 1.174, b = 18.562, R 2 = 0.956.

Table 2 .

 2 S°, Δ f H°, and Δ f G°Calculated by HI-G4-iso Method Using Different Iso/Homodesmotic Reactions for Tetralin as Compared to the Experimental Values of S°= 369.6 J•mol -1 •K -1 , Δ f H°= 26.61 kJ•mol -1 , and Δ f G°= 167.1 kJ•mol -1

Figure 4 .

 4 Figure 4. Parity plot for Δ f H°:( ⧫ ) Benson, ( □ ) HI-G4-iso, (○)HI-B3LYP-iso. Linear regression (y = mx + b): Benson: m = 1.004, R 2 = 0.999. HI-G4-iso: m = 1.008, R 2 = 0.997. HI-B3LYP-iso: m = 0.982, b = 7.326, R 2 = 0.999.

Figure 5 .

 5 Figure 5. Parity plot for standard entropy (S°): ( ⧫ ) Benson, ( □ ) HI-G4-iso, (○) HI-B3LYP-iso. Linear regression (y = mx + b): Benson: m = 0.994, R 2 = 0.998. HI-G4-iso: m = 1.024, b = -12.522, R 2 = 0.992. HI-B3LYP-iso: m = 1.048, b = -25.204, R 2 = 0.961.

Figure 6 .

 6 Figure 6. Parity plot for Gibbs energy estimation at 673K (Δ R G 673K °): ( ⧫ ) Benson, ( □ ) HI-G4-iso, (○) HI-B3LYP-iso. Linear regression (y = mx): Benson: m = 1.007, R 2 = 0.980. HI-G4-iso: m = 0.92, R 2 = 0.935. B3LYP-iso: m = 1.005, R 2 = 0.735.

Figure 7 .

 7 Figure 7. Hydrogenation of the first and second aromatic rings of 1MN.

Table 3 .

 3 Deviation of the Δ f H°Predicted by Different Methods to the Experimental Data Available for Δ f H°in kJ• mol -1

				Δ f H°B
		enson	B3LYP	HI-B3LYP-iso	G4	HI-G4-iso
	AARE (%)	6.2	237.0	28.9	38.2	10.4
	RMSE	2.9	78.9	8.8	17.4	6.3

Table 4 .

 4 Deviation of the S°Predicted by Different Methods to the Experimental Data Available for S°in J•mol -1 •K -1

				S°B
		enson B3LYP HI-B3LYP-iso	G4	HI-G4-iso
	AARE (%)	0.8	21.3	1.7	23.1	1.8
	RMSE	4.3	85.3	7.7	92.2	7.8

Table 5 .

 5 Deviation of the Δ R G 673K °Predicted by Different Methods from the Experimental Data Available forΔ R G 673K °in kJ•mol -1 Δ R G 673K °(kJ•mol -1 )Results from this calculation are shown in Table7. Calculated values using Benson's GA method are in good agreement with experimental results obtained by Bouchy et al. (an error around 5%).

		Benson	B3LYP	HI-B3LYP-iso	G4	HI-G4 -iso
	RMSE	3.2	10.1	14.0	9.4	6.4
	bias	-2.0	4.9	-0.1	-4.7	3.1

Table 6 .

 6 Flash Calculation Results at 40 bar and 400 °C Compared to Empirical Data

	liquid composition (% mol)
	molecule	Benson	HI-G4-iso	exptl. a
	2	1.55	0.08	37.90
	6	1.75	0.11	6.30
	3	0.95	0.05	22.20
	4	0.10	0.01	9.80
	7	0.11	0.04	1.60
	9	0.94	0.26	3.70
	20	0.60	0.01
	5 + 8	0.34	0.05	5.20
	21	0.04	0.08
	16	0.09	0.81	1.50
	18	0.12	0.09	0.50
	14	1.58	0.04
	17 +19	5.90	1.86	2.7
	15	20.38	6.23
	22 + 24	52.23	78.91	0.5
	23 + 25	9.67	11.04	0.1
	12	1.05	0.02	1.70
	13	0.46	0.04
	10	0.34	0.01	1.20
	11	1.81	0.27	1.10
	hydrogenation 1MN (%)	94.6	99.5	49.5
	conversion (%)	96.7	99.8	55.8

a The experimental values were taken from Bouchy et al.

27

Table 7 .

 7 Compositions Obtained from Flash Calculations at 40 bar and 400 °C for 1MN Hydrogenation Compared to Empirical Data The experimental values were taken from Bouchy et al.27 X 1MN , X xMT , and X 1MD refer to the molar fractions of 1MN, 1-methyl and 5methyltetralin, and trans-1-MD, respectively, with X 1MN + X xMT + X 1MD =1 .

			calculated
		exptl. a	Benson	HI-G4-iso
	X 1MN	0.52	0.54	0.44
	X xMT	0.40	0.39	0.33
	X 1MD	0.07	0.07	0.23

a

Table 8 .

 8 Equilibrium Constants Calculated at 673 K × 10 -04 6.70 × 10 -04 2 7.73 × 10 -06 3.66 × 10 -05 3 7.37 × 10 -05 1.94 × 10 -04 4 7.70 × 10 -05 1.27 × 10 -04

			K p °reaction
		Benson	HI-G4-iso
	1	7.34

AUTHOR INFORMATION

Corresponding Author *E-mail: herve.toulhoat@ifpen.fr.

Notes

The authors declare no competing financial interest.

ABBREVIATIONS


A Suffix

R Suffix

A Raji

J Rosato
email: joel.rosato@univ-amu.fr

R Stamm

Y Marandet

New analysis of Balmer line shapes in magnetic white dwarf atmospheres

We present an analysis of hydrogen Balmer line shapes observed in a white dwarf spectrum taken from the SDSS database. The spectrum presents lines with a triplet structure, featuring the presence of a strong magnetic field. Using a spectrum fitting model, we estimate the effective temperature of the star and we estimate the magnetic field present in the atmosphere. We also estimate the electron density from Stark broadening calculations.

Dear Author,

Here are the proofs of your article.

• You can submit your corrections online, via e-mail or by fax.

• For online submission please insert your corrections in the online correction form. Always indicate the line number to which the correction refers.

• You can also insert your corrections in the proof PDF and email the annotated PDF.

• For fax submission, please ensure that your corrections are clearly legible. Use a fine black pen and write the correction in the margin, not too close to the edge of the page.

• Remember to note the journal title, article number, and your name when sending your response via e-mail or fax.

• Check the metadata sheet to make sure that the header information, especially author names and the corresponding affiliations are correctly shown.

• Check the questions that may have arisen during copy editing and insert your answers/ corrections.

• Check that the text is complete and that all figures, tables and their legends are included. Also check the accuracy of special characters, equations, and electronic supplementary material if applicable. If necessary refer to the Edited manuscript.

• The publication of inaccurate data such as dosages and units can have serious consequences.

Please take particular care that all such details are correct.

• Please do not make changes that involve only matters of style. We have generally introduced forms that follow the journal's style. Substantial changes in content, e.g., new results, corrected values, title and authorship are not allowed without the approval of the responsible editor. In such a case, please contact the Editorial Office and return his/her consent together with the proof.

• If we do not receive your corrections within 48 hours, we will send you a reminder.

• Your article will be published Online First approximately one week after receipt of your corrected proofs. This is the official first publication citable with the DOI. Further changes are, therefore, not possible.

• The printed version will follow in a forthcoming issue.

Please note

After online publication, subscribers (personal/institutional) to this journal will have access to the complete article via the DOI using the URL: http://dx.doi.org/[DOI].

If you would like to know when your article has been published online, take advantage of our free alert service. For registration and further information go to: http://www.link.springer.com.

Due to the electronic nature of the procedure, the manuscript and the original figures will only be returned to you on special request. When you return your corrections, please inform us if you would like to have these documents returned. Fitting the spectrum using an adequate model will 58 give an estimate for the parameters characterizing the 59 plasma in the star's atmosphere. In this section, we con-60 sider the overall spectrum. The atmosphere is described 61 as a homogeneous slab, which absorbs the radiation 62 equation (e.g., [START_REF] Pomraning | The Equations of Radiation Hydrody-253 namics[END_REF]):

_####_ Page 2 of 4 Eur. Phys. J. D _#####################_
69

I ν = B ν exp(-κ ν L)+S ν (1 -exp(-κ ν L)), (1)
κ ν = n hν 4π N 2 B 2n φ n2 (ν), ( 2
)
82
where the sum is carried out over the upper levels presently no general formula accounting for all of these 96 effects and sufficiently simple so as to be applicable in a 97 fast fitting routine. In particular, the modeling of Stark 98 broadening is a major issue; we discuss it in the next 99 section. For the sake of simplicity in the fitting here, 100 we use a Lorentzian function with a width indepen-101 dent of n. The only effect retained here is the Zeeman 102 splitting, in such a way so as to reproduce the triplet 103 structure visible in the lines. Both linear and quadratic 104 Zeeman perturbations have been taken into account; for 105 the quadratic Zeeman effect, analytical formulas from 106 [START_REF] Garstang | [END_REF] have been used. Figure 2 shows the fitting result. 107 From the fit, the atmospheric temperature was esti-108 mated to be 19000 K, the magnetic field was found to 109 be 460 T, and a value of 3 × 10 20 m -2 was estimated 110 for N 2 L. The inferred values for the effective temper-111 ature and the magnetic field are in good agreement 112 with previous reports [7,8] (the magnetic field values 113 around 8 MG reported therein corresponds to the polar 114 field, which is approximately twice the observed mag-115 netic field). Note that there is no direct estimate of the 116 density N 2 through the model used here. An extension 117 providing the atomic density can be devised based on 118 a suitable line shape formula, accounting for collisional 119 broadening; an example of such a formula is provided 120 by the Lindholm impact theory and its variants [13]. 121 Databases providing line widths and shifts and appli-122 cable to isolated lines have precisely been developed in 123 this framework, e.g., see [14,15]. In white dwarf atmospheres, the electron density is 126 usually high enough so that the Stark effect due to 127 microfield provides a significant contribution to the line 128 broadening. Because of its strong dependence to the 129

i dU dt (t)=[H 0 + V (t)]U (t). ( 4
)
148 here H 0 is the Hamiltonian corresponding to the unper-149 turbed atom and the "unified theory" [START_REF] Capes | [END_REF]) or involving ad hoc stochastic 166 models (MMM [18,19], FFM [20,21]). Fully numerical 167 simulations can also be performed; in this case, they 168 may serve as a benchmark for models [22]. Essentially,

V (t)=-d • E(t) is the Stark effect

56 3 Fitting

 3 the overall spectrum 57

Fig. 1

 1 Fig. 1 Spectrum of the white dwarf considered in this work. The data were obtained from the SDSS database (www.sdss.org). The Balmer lines are split in triplets, which indicates a strong magnetic field

  70 where B ν is the Planck function, κ ν is the extinc-71 tion coefficient, L is the atmosphere thickness, and S ν 72 denotes the re-emission contribution. Assuming local 73 thermodynamic equilibrium (LTE), this term can be 74 approximated as S ν ≃ C ν exp(-hν/k B T ) where C ν 75 involves the photon absorption cross sections of bound-76 free and free-free transitions [11]. The frequency depen-77 dence of C ν is determined by Gaunt factors and is usu-78 ally weak; for simplicity, we will neglect this depen-79 dence in the analysis and leave C ν ≡ C as an adjustable 80 parameter. The extinction coefficient reads 81

83(

  n>2), N 2 is the density of absorbers, i.e., the density 84 of atom in the n = 2 level, B 2n is the Einstein coeffi-85 cient for absorption, and φ n2 (ν) is the normalized line 86 shape function relative to the 2 → n transition. The 87 stimulated emission has been retained through the for-88 mal substitution N 2 B 2n → N 2 B 2n [1 -exp(-hν/k B T )], 89 which holds at LTE. The line shape function φ n2 is 90 determined by the physical mechanisms that perturb 91 the energy levels, namely the Stark effect due to the 92 plasma microfield, the collisions with neutral atoms 93 (both resonant and van der Waals interactions), and 94 the Zeeman effect due to the magnetic field. There is

Fig. 2

 2 Fig. 2 Calculated best-fit spectrum for the white dwarf SDSS J124851.31-022924.73. A good agreement with the experimental white dwarf spectrum is obtained. See text for the fitting model

124 4

 4 Stark broadening modeling125

  150 term; E(t) is the electric field generated at the emitter's 151 location by the perturbers. Several approaches have 152 been developed for the evaluation of the spectrum (3), 153 usually based on the separation between the time scales 154 relative to the electrons and to the ions. For exam-155 ple, in the framework of the so-called standard model, 156 the ions are assumed motionless (static approximation) 157 and the electrons are described within a collisional pic-158 ture in terms of a collision operator (impact approxi-159 mation). Extensions to this model have been designed 160 in such a way so as to retain the time-dependence of the 161 ionic microfield (ion dynamics issue) and to account for 162 incomplete collisions with the electrons; such extensions 163 are usually based on kinetic theory approaches, involv-164 ing either particle distribution functions (e.g., through 165

169a

  simulation consists in numerically integrating the 170 time-dependent Schrödinger equation that governs the 171 dynamics of an atom perturbed by a fluctuating electric 172 field, itself being generated from a numerical integration 173 of the Newtonian equations of motion for the charged 174 particles moving in the vicinity of the emitter. We have 175 applied the computer simulation method reported in 176 [23] to the fitting of the shape of the Hβ absorption 177 line. One advantage of this line, with respect to the 178 other ones, is that it is far enough from its neighbors 179 (Hα and Hγ) so as not to overlap with them; it is also 180 deeper than Hα, which suggests this line is less subject 181 to re-emission. A focus has been put on Stark broad-182 ening; in particular, the broadening due to collisions

Fig. 3 A

 3 Fig. 3 A computer simulation method for Stark broadening modeling has been applied to the fitting of the Hβ line shape. The inferred density was Ne =1.5 × 10 17 cm -3

4 230

 4

  46 ure 1 shows the considered spectrum. The absorption 47 lines correspond to the Balmer series. They are well 48 resolved up to n = 6. The lines are split in triplets, 49 which indicates a strong magnetic field. This spectrum 50 has relatively weak noise near the lines, which makes it 51 a suitable candidate for analysis. Although there have 52 been many studies in white dwarf spectra, only a very 53 few discussed the presence of a strong magnetic field 54[7,8]. Here, we use a new fitting model which accounts 55 for this, including both Stark and Zeeman effects.

	1	1 Introduction
	2	White dwarfs are the post-end stage of stars of low or
	3	medium mass, typically M 10M ⊙ [1]. They are com-
	4	posed mostly of electron-degenerate matter, surrounded
	5	by a thin atmosphere of partially ionized plasma. Inves-
	6	tigations have shown that the majority of white dwarfs
	7	have an atmosphere of pure hydrogen as a result of
	8	gravitational setting, which removes helium and heavier
	9	elements from the atmosphere and moves them toward
	10	inner layers [2,3]. These atmospheres can be considered
	11	as hydrogen plasmas, which are similar to some that can
	12	be created in laboratory experiments [4]. As a rule, the
	13	physical properties of white dwarfs are inferred from
	14	spectroscopic analysis; the fitting of a spectrum can
	15	provide relevant parameters such as the surface tem-
	16	perature or the surface gravity, and modeling efforts
	17	are still ongoing in order to have accurate estimates
	18	[5,6]. The electron density in a white dwarf atmosphere is high enough (up to 10 17 cm -3 , and higher) so that
	20	the line shapes are affected by Stark broadening and,
		hence, can serve as a probe of the electron density N e .

21

Triplet lines, which are visible in some white dwarf spec-22 tra

[7,8]

, are the result of the Zeeman effect. The line 23 splitting is proportional to the magnetic field B and, 24 hence, provides a direct estimate of magnetic field. Fol-25 lowing previous works [9,10], we report on the analysis 26 of a white dwarf spectrum where the absorption lines 27 are subject to strong Zeeman splitting. A fitting model 28 that accounts for continuous background and Zeeman 29 effect is applied to the spectrum. We address the issue of 30 Stark broadening modeling by using a computer simu-31 lation technique. An estimate of the electron density in 32 the atmosphere is obtained based on Balmer line shape 33 a e-mail: joel.rosato@univ-amu.fr (corresponding author)

calculations. The article is organized as follows: Sect. 2 34 presents the spectrum investigated in this work, Sect. 4 35 presents the fitting model and the inferred values for the 36 temperature and the magnetic field, and Sect. 5 reports 37 on Stark broadening calculations that have been carried 38 out to infer the electron density. A discussion about the 39 fitting results is done. 40 2 Presentation of the white dwarf spectrum 41

Many white dwarf spectra have been recorded and 42 made available through dedicated databases, such as 43 the Sloan Digital Sky Survey (SDSS). In this work, we 44 investigate a spectrum observed in a white dwarf of DA 45 type, i.e., with an atmosphere rich in hydrogen. Fig-

Data Availability Statement This manuscript has no 231

density, this line broadening mechanism can be used as

Author Query Form

Please ensure you fill out your response to the queries raised below and return this form along with your corrections

Dear Author

During the process of typesetting your article, the following queries have arisen. Please check your typeset proof carefully against the queries listed below and mark the necessary changes either directly on the proof/online grid or in the 'Author's response' area provided below Query Details required Author's response 1.

Kindly check and confirm whether the corresponding author is correctly identified.

2.

Journal instruction requires a city for affiliations; however, these are missing in affiliation [START_REF] Shapiro | Black Holes, White 235 Dwarfs, and Neutron Stars: The Physics of Compact 236 Objects[END_REF]. Please verify if the provided city is correct and amend if necessary.

3.

Kindly check and confirm the updated publisher location is correct for the references [START_REF] Shapiro | Black Holes, White 235 Dwarfs, and Neutron Stars: The Physics of Compact 236 Objects[END_REF][START_REF] Pomraning | The Equations of Radiation Hydrody-253 namics[END_REF]. 4.

Please provide an explanation as to why there is no data or why the data will not be deposited. Your explanation will be displayed as 'Authors' comment'.


Mathias Géry
email: mathias.gery@univ-st-etienne.fr

Christine Largeron
email: christine.largeron@univ-st-etienne.fr

Franck Thollard
email: franck.thollard@univ-st-etienne.fr

UJM at INEX 2008: pre-impacting of tags weights

This paper 1 addresses the impact of structure on terms weighting function in the context of focused Information Retrieval (IR). Our model considers a certain kind of structural information: tags that represent logical structure (title, section, paragraph, etc.) and tags related to formatting (bold, italic, center, etc.). We take into account the tags influence by estimating the probability that a tag distinguishes relevant terms. This weight is integrated in the terms weighting function. Experiments on a large collection during INEX 2008 IR competition showed improvements for focused retrieval.

Introduction

The focused information retrieval (IR) aims at exploiting the documents structure (e.g. HTML or XML markup) in order to retrieve the relevant elements (parts of documents) for a user information need. The structure can be used to emphasize some particular words or some parts of the document: the importance of a term depends on its formatting (e.g. bold font, italic, etc.), and also on its position in the document (e.g., title terms versus text body).

Different approaches have been proposed to integrate the structure at the step of querying or at the step of indexing. Following [START_REF] Fuller | Coherent answers for a large structured document collection[END_REF], we propose to integrate the structure in the weighting function: the weights of terms are based not only on the terms frequencies in the documents and in the collection, but also on the terms position in the documents. This position can be defined by XML tags. This approach raises two questions: how to choose the structural weights? How to integrate them in the classical models? Some works propose to choose empirically the tags and their weights [START_REF] Rapela | Automatically combining ranking heuristics for html documents[END_REF] or to learn them automatically using genetic algorithms [START_REF] Trotman | Choosing document structure weights[END_REF]. These approaches use generally less than five tags. We propose to learn automatically the tags weights, without limit on the number of tags.

Concerning the integration of the structure weights, Robertson et al. suggests to preserve the non linearity of the BM25 weighting function by pre-impacting structure on the terms frequencies instead of impacting it directly on the global terms weights [START_REF] Robertson | Simple BM25 extension to multiple weighted fields[END_REF]. We propose to apply this approach in the context of focused XML IR.

The main contribution of this paper is a formal framework integrating structure, introduced in the next section. We present in section 3 our experiments and in section 4 our results in the INEX 2008 competition.

A structured document model

We consider in this paper the problem of extending the classical probabilistic model [START_REF] Robertson | Relevance weighting of search terms[END_REF] that aims at estimating the relevance of a document for a given query through two probabilities: the probability of finding a relevant information and the probability of finding a non relevant information.

Our model takes into account the structure at two levels. Firstly, the logical structure (e.g. tags section, paragraph, table, etc.) is used in order to select the XML elements that are handled at the indexing step. These elements are the only ones that can be indexed, ranked and returned to the user. Secondly, the formatting structure (e.g. bold font, italic, etc.) and the logical structure are integrated into the terms weighting function. For a given tag we can estimate if it emphasizes terms in relevant documents or term in non relevant part of documents. A learning step computes a weight for each tag, based on the probability, to distinguish relevant terms and non relevant ones. At querying step, the relevance of an element is estimated based on the weights of the terms it contains, combined with the weights of the tags labeling those terms.

Term based score of XML elements

The relevance of an element e j for a query Q is function of the weights of the query terms t i that appear in the element. We use the weighting function BM25 [START_REF] Robertson | Relevance weighting of search terms[END_REF]:

w ji = tf ji * (k 1 + 1) k 1 * ((1 -b) + (b * ndl)) + tf ji * log N -df i + 0.5 df i + 0.5 (1)
With tf ji : frequency of t i in e j ; N : number of elements in the collection; df i : number of elements containing the term t i ; ndl: ratio between the length of e j and the average element length; k 1 and b: classical BM25 parameters.

Tag based score of XML elements

The relevance of an element e j relatively to the tags is based on the weights, noted w ′ ik , of each term t i labeled by a tag b k . We used a learning set LS in which the relevant elements for a given query are known. Given the set R (resp. NR) that contains the relevant (resp. non relevant) elements, a contingency table can be built:

R NR LS = R ∪ N R t ik ∈ e j r ik nr ik = n ik -r ik n ik t ik / ∈ e j R -r ik N -n ik -R + r ik N -n ik Total R N R = N -R N
With R: number of relevant terms; N R: number of non relevant terms. r ik : number of times term t i labeled by b k is relevant; i r ik : number of relevant terms labeled by b k ; n ik : number of times term t i is labeled by b k ; nr ik = n ik -r ik : number of times term t i labeled by b k is not relevant.

Then, w ′ ik can be used to distinguish relevant terms from non relevant ones according to the tags that mark them. This is closely related to probabilistic IR model, but in our approach tags are considered instead of terms and terms instead of documents.

w ′ ik = P (t ik |R)(1 -P (t ik |N R)) P (t ik |N R)(1 -P (t ik |R)) = r ik × (N R -nr ik ) nr ik × (R -r ik ) (2)
Moreover, we hypothesize that the property for a tag to distinguish relevant terms does not depend on terms, i.e. the weight of a tag b k should be the same for all terms. We finally estimate for each tag b k a weight w ′ k :

w ′ k = ti∈T w ′ ik |T | (3)

Global score of XML elements

In order to compute a global score, we propose a linear combination f claw 2 between the weight w ji of a term t i and the average of the weights w ′ k of the tags b k that mark the term3 :

f claw (e j ) = t ik ∈ej /ti∈Q w ji × k/t ik =1 w ′ k |{k/t ik = 1}| (4)
In previous experiments [START_REF] Géry | Integrating structure in the probabilistic model for information retrieval[END_REF], f claw slightly improved recall but the results were not convincing. Even if the estimation of the tag weights must be carefully addressed, it appears that the way such weights are integrated into the final score is essential. Following [START_REF] Robertson | Simple BM25 extension to multiple weighted fields[END_REF], we take advantage of the non linearity of BM25 by pre-impacting the tags weights at the term frequency level. More precisely, tf is replaced by ttf4 in BM25:

ttf ji = tf ji × k/t ik =1 w ′ k |{k/t ik = 1}| (5)
3 Experiments

We have experimented these models during the INEX 2008 IR competition in a classic IR way (granularity: full articles) as well as in a focused IR way (granularity: XML elements). The English Wikipedia XML corpus [START_REF] Denoyer | The wikipedia XML corpus[END_REF] contains 659,388 strongly structured articles, which are composed of 52 millions of XML elements (i.e. 79 elements on average; with an average depth of 6.72). The whole articles (textual content + XML structure) represent 4.5 Gb while the textual content only 1.6 Gb. The original Wiki syntax has been converted into XML, using both general tags of the logical structure (article, section, paragraph, title, list and item), formating tags (like bold, emphatic) and frequently occurring link-tags.

Experimental protocol

The corpus enriched by the INEX 2006 assessments on 114 queries has been used as a training set in order to estimate the tags weights w ′ k . We have evaluated our approach using the 70 queries of INEX 2008.

Our evaluation is based on the main INEX measures (iP [x] the precision value at recall x, AiP the interpolated average precision, and M AiP the interpolated mean average precision [START_REF] Kamps | INEX 2007 evaluation measures[END_REF]). Note that the main ranking of INEX competition is based on iP [0.01] instead of the overall measure M AiP , in order to take into account the importance of precision at low recall levels.

Each run submitted to INEX is a ranked list containing at most 1 500 XML elements for each query. Some runs retrieve all the relevant elements among the first 1 500 XML returned elements, and some others retrieve only part of them. Note that a limit based on a number of documents (instead of e.g. a number of bytes) allows to return more information and therefore favors runs composed by full articles. We have calculated R[1500] (the recall at 1 500 elements) and S[1500] (the size of these 1 500 elements in Mbytes).

Tags weighting

We have manually selected 16 tags (article, cadre, indentation1, item, li, normallist, numberlist, p, row, section, table, td, template, th, title, tr ) in order to define the XML elements to consider. These logical structure tags will be considered during the indexing step and therefore those will define the elements the system will be able to return.

Regarding the other tags (namely the formatting tags), we first selected the 61 tags that appear more than 300 times in the 659,388 documents. We then manually removed 6 tags: article, body (they mark the whole information), br, hr, s and value (considered not relevant).

The weights of the 55 remaining tags were computed according to equation w ′ k in equation 3. Table 1 presents the top 6 tags and their weights, together with the weakest 6 ones and their weights. Their frequencies in the whole collection is also given. 4 Results: focused task

Our aim was firstly to obtain a strong baseline, secondly to experiment focused retrieval (i.e. elements granularity) against classic retrieval (i.e. full articles granularity), and thirdly to experiment the impact of tags weights in the BM25 weighting function. Table 2 presents the 3 runs that we have submitted to INEX 2008 Ad-Hoc in focused task. The structure is not taken into account in R1, where the documents are returned to the user (articles granularity) as well as in R2 where the elements are returned (elements granularity), while in R3 the tags weights are integrated in BM25 in a focused retrieval (elements granularity -TTF)

Parameters

The parameters of the chosen weighting functions (namely BM25) were tuned in order to improve classic retrieval (articles granularity) and focused retrieval (elements granularity). Among the parameters studied to improve the baseline, we can mention the use of a stoplist, the optimization of BM25 parameters (k 1 = 1.1 and b = 0.75), etc. Regarding the queries, we set up a better "andish" mode and consider or and and, etc . . . . Some specific parameters (e.g. the minimum size of the returned elements) were also tuned for focused retrieval.

Our baseline and all other runs have been obtained automatically, and using only the query terms (i.e the title field of INEX topics). We thus do not use fields description, narrative nor castitle.

INEX ranking: iP [0.01]

Our system gives very interesting results compared to the best INEX systems. Our runs are compared on the figure 1 against F OERStep, the best run submitted to INEX 2008 according to iP [0.01] ranking, on 61 runs yet evaluated in the focused task. This run outperforms our runs at very low recall levels. Our run R1 gives the best results at recall levels higher than 0.05. This is also shown by the M AiP presented in table 3. 3 presents the results of our 3 runs submitted to the track Ad-Hoc (focused task).

Articles versus elements

Our second aim was to compare classic retrieval of full articles versus focused retrieval of XML elements. We therefore indexed either the whole articles or the elements, and the parameters of the system were tuned also for focused retrieval.

It is interesting to notice that the BM25 model applied on full articles (R1) outperforms our focused retrieval results (R2) considering M AiP , despite the fact that BM25 parameter ndl is designed to take into account different documents lengths and thus documents granularities. Classic IR weighting functions, indexing and querying process, are undoubtedly not well adapted to focused retrieval. However, this is consistent with other results obtained during the INEX 2007 campaign where some top ranked systems only consider (and therefore return) full articles.

On the other hand, the focused run R2 returns a smallest quantity of information. Indeed, the total size of the 1 500 XML elements returned (for each query) is reduced to 51 Mb instead of 390 Mb for classic retrieval of full articles.

Pre-impacting of tags weights on terms weights

Finally, our third aim was to experiment the impact of tag weights in term weighting function in a focused retrieval scheme. In order to understand the pro and cons of our structured model, the weighting functions and the same parameters used for the baseline runs were also used with our structured model.

The figure 1 shows that our TTF strategy (R3) improves dramatically the focused retrieval at low recall levels (from 0.5697 to 0.6640 following iP [0.01] ranking). However, it does not improve focused retrieval enough to reach better results than classic retrieval.

These results confirm also that, according to Robertson and al. [START_REF] Robertson | Simple BM25 extension to multiple weighted fields[END_REF], it is important to keep the non linearity of the BM25 weighting function by "preimpacting" term position in the structure of document (in other terms, tags weights) on the terms frequencies (strategy TTF) instead of "post-impacting" it directly on the terms weights (strategy CLAW, cf. [START_REF] Géry | Integrating structure in the probabilistic model for information retrieval[END_REF]).

Conclusion

We proposed in [START_REF] Géry | Integrating structure in the probabilistic model for information retrieval[END_REF] a new way of integrating the XML structure in the classic probabilistic model. We consider both the logical structure and the formatting structure. The logical structure is used at indexing step to define elements that correspond to part of documents. These elements will be indexed and potentially returned to the user. The formatting structure is integrated in the document model itself. During a learning step using the INEX 2006 collection, a weight is computed for each formatting tag, based on the probability that this tag distinguishes relevant terms. During the querying step, the relevance of an element is evaluated using the weights of the terms it contains, but each term weight is modified by the weights of the tags that mark the term.

The baselines are rather strong as the score of the BM25 run on article (run R1) is ranked seven of the competition according to the iAP [0.01] ranking.

Our strategy TTF gives better results than focused retrieval (R2 ) and classic retrieval (R1 ) at low recall levels (iP [0.01]). That shows the interest of focused IR (R3 vs R1 ), and the interest of using structure (R3 vs R2 ). Pre-impacting the structure on terms frequencies (TTF, R3 ) gives also better results than "post-impacting" it on final terms weights (CLAW, [START_REF] Géry | Integrating structure in the probabilistic model for information retrieval[END_REF]). Actually, T T F changes significantly the performances of the methods when considering the iP [0.01] or the M AiP measure.

TTF (R3 ) gives also good recall results (M AiP = 0.2347; R[1500] = 0.6097). Focused IR eliminates more non relevant elements than relevant elements (R3 vs R1 ): R[1500] decreases by 16% while S decreases by 40%.

We have presented a document model integrating explicitly the structural information in the weighting function, and a learning process of tags weights. We reach the same conclusions than [START_REF] Robertson | Simple BM25 extension to multiple weighted fields[END_REF] about the interest of pre-impacting structure, with a very different collection, a more heterogeneous one that contains a much larger set of tags (> 1 thousand).

In previous experiments, a basic average function, that considers all the tags equally (CLAW), gives better results than other combining functions (multiplication, only the closest tag, etc.). But, we think that a finest combining function (e.g. taking into account the distance between terms and tags) should improve the results.

Fig. 1 .

 1 Fig. 1. Recall / Precision of 3 runs on 61 runs yet evaluated in the focused task

Table 1 .

 1 Weight w ′ k of the 6 strongest and 6 weakest tags

	Top strongest weights	Top weakest weights
	tag	weight	freq. tag	weight freq.
	h4	12,32		307 emph4 0,06	940
	ul	2,70	3'050 font	0,07 27'117
	sub	2,38	54'922 big	0,08 3'213
	indentation1 2,04	135'420 em	0,11	608
	section	2,01 1'610'183 b	0,13 11'297
	blockquote	1,98	4'830 tt	0,14 6'841

Table 2 .

 2 Our 3 runs submitted to INEX 2008 Ad-Hoc, focused task

	Run (name)	Granularity Tags weights
	R1 (JMU expe 136) articles	-
	R2 (JMU expe 141) elements	-
	R3 (JMU expe 142) elements	TTF

Table 3 .

 3 Our 3 runs compared to 61 "focused" runs

	Run (rank)	iP [0.01] M AiP R[1500] S[1500]
	FOERStep (winner)	0.6897 0.2076 0.4468 78
	R1: articles (14)	0.6412 0.2801 0.7871 390
	R2: elements (37)	0.5697 0.1208 0.2761 51
	R3: elements+TTF (9) 0.6640 0.2347 0.6097 234

This work has been partly funded by the Web Intelligence project (région Rhône-Alpes, cf. http://www.web-intelligence-rhone-alpes.org).

CLAW: Combining Linearly Average tag-Weights.

wji: the BM25 weight of term ti in element ej, cf. eq. 1.

TTF: Tagged Term Frequency. t ik = 1 means that ti is labelled by b k .


Mathias Géry
email: mathias.gery@univ-st-etienne.fr

Christine Largeron
email: christine.largeron@univ-st-etienne.fr

Christophe Moulin
email: christophe.moulin@univ-st-etienne.fr

UJM at INEX 2008 XML mining track

This paper 1 reports our experiments carried out for the INEX XML Mining track, consisting in developing categorization (or classification) and clustering methods for XML documents. We represent XML documents as vectors of indexed terms. For our first participation, the purpose of our experiments is twofold: Firstly, our overall aim is to set up a categorization text only approach that can be used as a baseline for further work which will take into account the structure of the XML documents. Secondly, our goal is to define two criteria (CC and CCE) based on terms distribution for reducing the size of the index. Results of our baseline are good and using our two criteria, we improve these results while we slightly reduce the index term. The results are slightly worse when we sharply reduce the size of the index of terms.

Introduction

The INEX XML Mining Track is organized in order to identify and design machine learning algorithms suited for XML documents mining [START_REF] Denoyer | Report on the xml mining track at inex 2007 categorization and clustering of xml documents[END_REF]. Two tasks are proposed: clustering and categorization. Clustering is an unsupervised process through which all the documents must be classified into clusters. The problem is to find meaningful clusters without any prior information. Categorization (or classification) is a supervised task for which, given a set of categories, a training set of preclassified documents is provided. Using this training set, the task consists in learning the classes descriptions in order to be able to classify a new document in one of the categories.

This second task is considered in this article. Moreover, even if the content information (the text of the documents), the structural information (the XML structure of the documents) and the links between the documents can be used for this task, we have only exploited the textual information. Indeed, this is our first participation to this track and our aim was to design a framework that could be used as a baseline for further works dealing with structured documents.

More precisely, we focus on the preprocessing step, particularly the features selection, which is an usual step of the knowledge discovery process [START_REF] Yang | A comparative study on feature selection in text categorization[END_REF][START_REF] Joachims | Text categorization with support vector machines: Learning with many relevant features[END_REF][START_REF] Forman | An extensive empirical study of feature selection metrics for text classification[END_REF]. On textual data, this step can be essential for improving the performance of the categorization algorithm. It exists a lot of words in the natural language, including stop words, synonymous, etc.. These words are not equally useful for categorization. Moreover, their distribution must also be considered. For example, words that appear in a single document are not useful for the categorization task.

So, we need to extract from the text a subset of terms that can be used to efficiently represent the documents in view of their categorization. In this paper, the documents are represented according to the vector space model (VSM [START_REF] Salton | Introduction to modern information retrieval[END_REF]). Our aim is to adapt some VSM principles, for example the measure of the discriminatory power of a term, to the categorization task. We propose two criteria based on terms distribution aiming at extracting the indexing terms from the training set corpora. After a brief presentation of the VSM given to introduce our notations in section 2, these criteria are defined in the following section. Our categorization approach is described in section 3 while the experiments and the obtained results are detailed in sections 4 and 5.

2 Document model for categorization

Vector space model (VSM)

Vector space model, introduced by Salton and al. [START_REF] Salton | Introduction to modern information retrieval[END_REF], has been widely used for representing text documents as vectors which contain terms weights. Given a collection D of documents, an index T = {t 1 , t 2 , ..., t |T | }, where |T | denotes the cardinal of T , gives the list of terms (or features) encountered in the documents of D. A document d i of D is represented by a vector d i = (w i,1 , w i,2 , ..., w i,|T | ) where w i,j represents the weight of the term t j in the document d i . In order to calculate this weight, TF.IDF formula can be used.

TF: term representativeness

TF (Term Frequency), the relative frequency of term t j in a document d i , is defined by: tf i,j = n i,j l n i,l where n i,j is the number of occurrences of t j in document d i normalized by the number of terms in document d i . The more frequent the term t j in document d i , the higher is the tf i,j .

IDF: discriminatory power of a term

IDF (Inverse Document Frequency) measures the discriminatory power of the term t j . It is defined by:

idf j = log |D| |{d i : t j ∈ d i }|
where |D| is the total number of documents in the corpus and |{d i : t j ∈ d i }| is the number of documents in which the term t j occurs at least one time. The less frequent the term t j in the collection of documents, the higher is the idf j .

The weight w i,j of a term t j in a document d i is then obtained by combining the two previous criteria:

w i,j = tf i,j × idf j
The more frequent the term t j is in document d i and the less frequent it is in the other documents, the higher is the weight w i,j .

Criteria for features selection

This VSM is widely used for text mining and information retrieval, as well for free format document like scientific articles as for semi structured document written in markup languages like HTML or XML.

But, in the context of categorization, even for limited collections, the dimensionality of the index can be exceedingly large. For example, in INEX collection, 652 876 non trivial words have been identified. This is a real problem for categorization since the terms belonging to this bag of words are not necessarily discriminatory features of the categories. So, we introduced two criteria (CC and CCE) in order to select a subset of T providing a description of the documents belonging to the same category. We consider that these terms must be very frequent in the documents of the category and, on the contrary, that they must be infrequent in the other categories.

Category Coverage criteria (CC)

Let df k j be the number of documents in the category C k where term t j appears, and f k j be the frequency of documents belonging to C k and including t j :

df k j = |{d i ∈ C k : t j ∈ d i }|, k ∈ {1, ...r} (1)
f k j = df k j |C k | (2)
The higher the number of documents of C k containing t j , the higher is f k j . On the other hand, the term t j could be considered as a discriminant term if most of the documents, where t j appears, belongs to the same category. Thus, a first criteria, noted CC (Category Coverage), is computed as follows:

CC k j = df k j |C k | * f k j k f k j CC k j = (f k j ) 2 k f k j
If the value of CC k j is high, then t j is a characteristic feature of the category C k .

Category Coverage Entropy criteria (CCE)

The frequency f k j considers the number of documents containing t j but it does not take into account the number of occurrences of t j in the category. It is the reason why we consider also p k j the frequency of t j in the category C k and a measure commonly used in information theory, called entropy, which evaluates the purity of the categories for the term t j . In the context of text categorization, it measures the discriminatory power of t j . Let n k j be the number of occurrences of t j in the documents of C k and p k j the corresponding frequency:

n k j = di∈C k n i,j p k j = n k j k=1,r n k j
The Shannon entropy E j of the term t j is given by [START_REF] Shannon | A mathematical theory of communication[END_REF]:

E j = - k=1,r (p k j ) * log 2 (p k j )
The entropy is minimal, equal to 0, if the term t j appears only in one category. We consider that this term might have a good discriminatory power for the categorization task. Conversely, the entropy is maximal if t j is not a good feature for representing the documents i.e. if t j appears in all the categories with the same frequency.

We propose a second criteria, denoted CCE (Category Coverage Entropy), combining f k j (from CC) and entropy. CCE is defined by:

CCE k j = (alpha * f k j ) + (1 -alpha) * (1 - E j M axE )
where alpha is a parameter and M axE is the maximal value of E. When the term t j is characteristic of the category C k , the value of the criteria is high.

For each category, a subset of the terms of T corresponding to the highest values of the criterion is built. Then, the index is defined as the union of these subsets.

Experiments

Collection INEX XML Mining

The collection is composed of 114 336 XML documents of the Wikipedia XML Corpus. This subset of Wikipedia represents 15 categories, each corresponding to one subject or topic. Each document of the collection belongs to one category. In the XML Mining Track, the training set is composed of 10% of the collection.

Preprocessing

The first step of the categorization approach that we propose, consists in a preprocessing of the collection. It begins by the construction of the list all the terms (or features) encountered in the documents of the collection. This index of 652 876 terms is build using the LEMUR software2 . The Porter Algorithm [START_REF] Porter | An algorithm for suffix stripping[END_REF] has also been applied in order to reduce different forms of a word to a common form. This operation reduces the index to 560 209 terms. However, it still remains a large number of irrelevant terms that could degrade the categorization, e.g.: numbers (7277, -1224, 0d254c, etc.), terms with less than three characters, terms that appear less than three times, or terms that appear in almost all the documents of the training set corpus. The index obtained at this stage is denoted I. In our experiments, its size is reduced to 161 609 terms on all the documents of the collection and to 77 697 on the training set.

Features selection

However, as explained in the previous section, the terms of I are not necessarily appropriated for the categorization task inasmuch they are not discriminatory for the categories. This is the reason why our criteria based on entropy and on frequency are used to select more suited features. The terms were sorted according to CC and CCE and only those corresponding to the highest values are retained. In our experiments, the top 100 terms by class and the top 10 000 terms by class were considered for each criteria to build four indexes, denoted respectively CC 100 and CC 10000 using CC and CCE 100 and CCE 10000 using CCE. Table 1 indicates the size of these indexes. 1. Indexes sizes Using one of these indexes, the content of a document is then represented by the tf.idf vector model described in the first section.

The second step is the categorization step itself. Two usual methods of classification are used: Support Vector Machines (SVM) and k-nearest neighbors. Only the most promising results obtained with the SVM were submitted. SVM was introduced by Vapnik for solving two class pattern recognition problems using Structural Risk Minimization principal [START_REF] Vapnik | The Nature of Statistical Learning Theory[END_REF]. In our experiments, the SVM algorithm available in the Liblinear library3 has been used. The results provided by this approach are presented in the next section.

Experimental results

This work has been done with a dual purpose: firstly develop a categorization text approach usable as a baseline for further work on XML categorization taking into account the structure, and secondly evaluate performances of this method using our selection features approach.

Global results

We have submitted 5 experiments using our 5 indexes presented in

Baseline results

Our baseline corresponds to the first experiment (expe 1), which was ranked 3th with a quite good recall: 0.7873.

Selection features improves results

When we select 10 000 terms for each class using CCE (expe 3) and CC (expe 5), we reduce the size of the index to respectively 77 580 and 75 181. This reduction is small compared to the size of the baseline index (77 697). However, it lets us to slightly improve our baseline to 0.7874 with CCE and 0.7876 with CC. These three runs obtained the three best results of the XML Mining challenge.

Selection features reduces indexes

The last two submitted runs correspond to the selection of the first 100 terms for each class using CCE (expe 2) and CC (expe 4). As presented in table 1, the size of the index is sharply reduced to 909 terms for CCE and 1 051 for CC. This reduction respectively correspond to 85% and 74% of the size of the baseline index. Even if the obtained results are lower than the results obtained with larger indexes, they are still relatively good. Indeed, the obtained recall is 0.6770 with CCE and 0.7230 with CC.

Conclusion

We proposed a categorization text approach for XML documents that let us obtain a good baseline for further work. For now we just used CC and CCE criteria as a threshold to select terms in order to build the index. For future work, we aim at exploiting the computed value of CC and CCE to improve the categorization. Moreover, we could use the structure information of XML documents represented by the links between document to improve even more the results.

table 1 .

 1 The global results of XML Mining 2008 are synthesized in table 2 (participant: LaHC).

	Rank Participant Run	Recall Documents
	1	LaHC submission.expe 5.tf idf T5 10000.txt 0.7876 102 929
	2	LaHC submission.expe 3.tf idf T4 10000.txt 0.7874 102 929
	3	LaHC submission.expe 1.tf idf TA.txt	0.7873 102 929
	4	Vries	Vries classification text and links.txt 0.7849 102 929
	5	boris	boris inex.tfidf.sim.037.it3.txt	0.7379 102 929
	6	boris	boris inex.tfidf1.sim.0.38.3.txt	0.7347 102 929
	7	boris	boris inex.tfidf.sim.034.it2.txt	0.7309 102 929
	8	LaHC submission.expe 4.tf idf T5 100.txt 0.7230 102 929
	9	kaptein kaptein 2008NBscoresv02.txt	0.6980 102 929
	10	kaptein kaptein 2008run.txt	0.6978 102929
	11	romero romero naive bayes links.txt	0.6813 102 929
	12	LaHC submission.expe 2.tf idf T4 100.txt 0.6770 102 929
	13	romero romero naive bayes.txt	0.6767 102 929
	14	Vries	Vries classification links only.txt	0.6232 102 929
	15	Vries	Vries classification text only.txt	0.2444 92 647

Table 2 .

 2 Summary of all XML Mining results.

This work has been partly funded by the Web Intelligence project (région Rhône-Alpes, cf. http://www.web-intelligence-rhone-alpes.org).

Lemur is available at the URL http://www.lemurproject.org

http://www.csie.ntu.edu.tw/ cjlin/liblinear/ -L2 loss support vector machine primal


A Kalitsov

M Chshiev

B Canals

C Lacroix

Voltage induced control and magnetoresistance of noncollinear frustrated magnets

Keywords: numbers: 72.15.Eb, 72.15.Gd, 75.10.Lp, 75.10.Jm, 75.47.-m, 75.47.Pq

Noncollinear frustrated magnets are proposed as a new class of spintronic materials with high magnetoresistance which can be controlled with relatively small applied voltages. It is demonstrated that their magnetic configuration strongly depends on position of the Fermi energy and applied voltage. The voltage induced control of noncollinear frustrated materials (VCFM) can be seen as a way to intrinsic control of colossal magnetoresistance (CMR) and is the bulk material counterpart of spin transfer torque concept used to control giant magnetoresistance in layered spin-valve structures.

The discovery of giant magnetoresistance (GMR) [1,2] in magnetic multilayered structures has generated a new field of spin-based electronics [3,4], or spintronics, which combines two traditional fields of physics: magnetism and electronics. A spin-valve concept [5] used in GMR structures allows controlling the magnetic configuration of its ferromagnetic layers by (i) application of relatively small magnetic fields or (ii) passing spin polarized currents using spin transfer torque (STT) [6,7]. These factors made them ideal systems for spintronic applications such as magnetic random access memories (MRAM) and magnetic field sensors used in read heads.

The advent of GMR has considerably increased an interest in related phenomenon in bulk materials, colossal magnetoresistance (CMR) [8,9], which is several orders higher than GMR and unlike the latter, can be viewed as an "intrinsic" property of material itself. To date, the CMR is typically observed in certain manganite compounds with the bulk magnetic configuration controlled by applying the magnetic field (similar to method (i) mentioned above for spin valves) but requires characteristic magnetic fields of several Tesla [9]. Such high fields make them inappropriate for use in spintronic applications where appropriate scale should be about Oersteds. However, one may expect the possibility of controlling the intrinsic magnetic configuration of the bulk materials (and thus of CMR) by passing spin polarized currents through them similar to STT mechanism (ii) mentioned above for spin-valves. Since the STT in the latter originates from non-collinearity of their adjacent magnetizations, the same requirement should hold for the bulk materials.

Here we promote for the first time magnetically frustrated bulk materials as a new paradigm for spintronic applications with high magnetoresistance which can be controlled with relatively small applied voltages and does not require spin polarized currents. This novel phenomenon may be viewed as the "bulk" counterpart of STT in layered spintronic structures (spin valves) and may represent a crucial interest both from applications and fundamental points of view. Below we demonstrate that the magnetic configuration of the bulk frustrated material is changed under applied voltage leading to the strong variation of its conductance.

The key mechanism at stake is a somewhat microscopic equivalent of spin torque, allowing for local spin flips within the microscopic spin texture. Those local moves requires non collinearity of the spin configuration, which would otherwise be insensitive to the current. Furthermore, in order to avoid a continuous response of the spin configuration, we require the ground states to be well defined minima, such that the switching from one to the other looks like a logical operation. These two aspects are well taken into account as soon as the underlying localized moments possess strong multiaxial anisotropies, in order to behave bitwise like and fulfill the non collinearity criterion. Several Rare-Earth-transition metal intermetallics systems are relevant candidates : in these systems frustration arises due to the competition of crystal field anisotropy, exchange and quadrupolar interactions (for a review see Ref. [10]). These intermetallics systems often show a non-collinear magnetic structure: this is the case for example of TbGa 2 [11], HoGe 3 [12] or Uranium compounds [13]. Naturally, pyrochlores in which frustration is due to the crystal structure are also relevant candidates, as it is well known that in such exotic 0.0 0.5 1.0 1.5

-7

E F =-2.5 eV E tot (eV) Voltage (V)
4in 2in2out

E F =-2 eV E tot (eV) 4in 2in2out (a) (b) 0.0 0.5 1.0 1.5 2.0 -7 -6 -8
-7

E F =-2.5 eV E tot (eV) Voltage (V)
4in 2in2out systems, non collinear low temperature magnetic phases are often stabilized [14], among those, few metallic compounds have been identified like Pr 2 Ir 2 O 7 [15,16] and the family of the pyrochlore Molybdates, R 2 Mo 2 O 7 [17].

E F =-2 eV E tot (eV) 4in 2in2out (a) (b)
For argument's sake, we have used a two-dimensional square lattice model in order to focus on the driving mechanism that is at stake and to describe its qualitative behavior. Realistic realizations would involve more complicated structures, like the one in the 3-D pyrochlores [14]. We consider a square lattice of classical localized moments S i with strong local on-site uniaxial anisotropy (D 0 ) along square diagonals (represented by unit vector n i ) and coupled through intersite exchange I ij . In addition, moments S i are coupled with conduction electrons through the local exchange (J 0 ) described here quantum mechanically using tight-binding model. Unlike collinear systems, non-collinearity provided by magnetic frustration is a key ingredient for switching phenomenon proposed here since its origin is due to STT mechanism acting locally by conduction electrons on localized mo-ments S i . The Hamiltonian of the system has the form:

Ĥ = - i,j I ij S i • S j -D 0 i (n i • S i ) 2 + t i,j,σ ĉ †σ i ĉσ j + H.c. -J 0 i ĉ †α i (σ αβ • S i ) ĉβ i , (1)
where I ij and J 0 are the intersite and the local exchange constants respectively, D 0 is the uniaxial anisotropy constant, t is the hopping integral between two neighboring sites, ĉ †σ i and ĉσ i are the creation and annihilation operators of the conduction electron with the spin σ on site i and σ αβ is the vector of Pauli matrices. For the chosen model, nearest neighbor interactions I 1 are irrelevant (because nearest neighbor spins are always orthogonal) so that only the second nearest neighbor one I 2 is taken into account. We emphasize that dealing with longer ranges interaction does not affect further results reported here since only energy differences between magnetic configurations matter.

In the limit of low temperature and D 0 → ∞ the localized moments S i are strictly collinear with n i (S i n i ) and the total energy of the system within a constant shift becomes

E tot = -I 2 i,j S i • S j + T r[ Ĥ ρ] = = -I 2 i,j S i • S j - i 2π j EG < jj (E)dE, (2)
where G < and ρ represent the "lesser" non-equilibrium Green function and density matrix respectively [START_REF] Lifshitz | Physical Kinetics, Course of Theoretical Physics[END_REF], and the Green function indices include spin index as well.

We adopt a conventional transport approach where the system is considered to consist of three regions: semiinfinite left (L) and right (R) regions connected to the middle scattering region (M ) which is assumed to be very long. Despite of all three regions being identical, such a subdivision is chosen in order to consider all nonequilibrium processes occuring in the scattering region M while the left and right regions are in thermodynamical equilibrium described with Fermi-Dirac distribution functions f L(R) = f (E -µ L(R) ), where µ L(R) the chemical potential in the left(right) region. The "lesser" Green function can be written as

G < ij = i(f L ψ L i ψ L * j + f R ψ R i ψ R * j ), (3)
where ψ L(R) j

is the single-electron wave function on site j ∈ M incident from the left (right). It is straightforward from the definition of retarded and advanced Green functions that

G a jj -G r jj = i(ψ L j ψ L * j + ψ R j ψ R * j ). (4)
Since the middle region is assumed to be very long, the finite potential drop eV = µ L -µ R results into infinitesmal change from one unit cell to another within the middle region M . In this case the reflections from the L|M and R|M boundaries are negligible and one can write for the left(right) wave functions ψ

L(R) j = c L(R) exp (+(-)k x aj),
where k x is the wave vector along x-axis and a is the lattice constant. One can show that the expression for the charge current which is proportional to (f

L -f R ) [19] is satisfied when |ψ L j | 2 = |ψ R j | 2 .
Indeed, the charge current can be expressed using the lesser Green function as [20]:

I ∼ [(G < j+1,j -G < j,j+1 )]dE.
Using Eq. ( 3) and aforementioned wave functions it yields the following form:

I ∼ sin k x a f L |c L | 2 -f R |c R | 2 dE which is proportional to (f L -f R ) only when |c L | 2 = |c R | 2 (and consequently |ψ L j | 2 = |ψ R j | 2
). Then, it straightforwardly follows from Eq. ( 3) and Eq. ( 4) that

G < jj (E) = i(f L + f R )ℑ(G r jj (E)) (5)
where ℑ notation for imaginary part is used. The total energy given by expression (2) of the system per a unit cell yields:

E tot = -I 2 i,j S i • S j - 1 2π (f L + f R ) j ℑ(G r jj (E))EdE, ( 6
)
where G r jj (E) is calculated as jth diagonal element of the matrix [E -Ĥ + iδ] -1 . Thus, the magnetic configuration corresponding to the minimum of the total energy is resulting from the interplay of intersite exchange interaction between localized moments on the next nearest neighbour sites (the first term in Eq. ( 6)) and the exchange interaction between the localized moment and the conduction electron described by the second term in Eq. ( 6).

For simplicity, among the possible magnetic configurations of the system (S i n i ) we select "4in" and "2in2out" ones which are the two configurations with the largest energy separation in this example (see Fig. 1). We choose for parameters J 0 = 2 eV , t = 1 eV and I 2 = 0.1t [21].

In Fig. 2 we show the dependence of the total energy for the two aforementioned selected configurations as a function of applied voltage. First, one can note that E tot strongly depends on the Fermi level position which is used to represent different materials. Indeed, when the system is in equilibrium, its "2in2out" ("4in") state is energetically more favorable for E F = -2 eV (-2.5 eV ) as shown in Fig. 2 strong variation of the total energy as a function of applied voltage causes the system to switch from an equilibrium state to another one at certain critical voltage (see Fig. 2). Moreover, further increase of the applied voltage may again reverse the state of the system (Fig. 2a). This demonstrates that the magnetic configuration of the system can be controlled by the applied voltage.

The mechanism of these dependences may be understood from the corresponding total density of states (DOS) [START_REF]The Density of States is equal to -1/π P j ℑG R jj (E), where summation includes all atoms in the unit cell and spins[END_REF] for "2in2out" and "4in" configurations represented in Fig. 3 where for convenience purposes only negative energy range is shown since DOS(E)=DOS(-E). The DOS for both magnetic states have sharp peaks and band gaps causing strong dependence of the preferable configuration on the Fermi level E F both in and out of equilibrium. Indeed, as on can see from the second term in Eq. ( 6), the total energy is defined by the sum of two products of the DOS with E integrated till µ L and µ R . In the absence of applied voltage, E tot is calculated with µ L = µ R = E F , i.e. defined by the position of the Fermi level indicated by black arrows in Fig. 3 corresponding to two values -2 eV and -2.5 eV used in Fig. 2(a) and (b), respectively.

One might expect that the conductances of two aforementioned magnetic configurations may strongly differ due to significant differences in their DOS, thereby leading to high magnetoresistance values defined as (σ 2in2out xx -σ 4in xx )/σ 4in xx ), where σ 4in xx and σ 2in2out xx are the conductances for "4in" and "2in2out" magnetic configurations, respectively. We calculate the linear response conductance using the Kubo formula [22]

σ xx = πe 2 L T r vx δ(E -Ĥ)v x δ(E -Ĥ) ,
where vx is the electron velocity operator along the applied voltage direction, L → ∞ is the length of the system and

δ(E -Ĥ) = - 1 π ℑ( Ĝr (E)) = i π ( Ĝr (E) -Ĝa (E))
The calculated conductance of the system for "4in" and "2in2out" magnetic configurations are represented in Fig. 4. As one can see, both σ 4in xx and σ 2in2out xx strongly depend on position of the Fermi energy E F and strongly correlate with DOS picture (see Fig. 3). For the two values of the Fermi energy used above, the MR ratio is found to be about 336 % for the case when E F = -2.5 eV and becomes infinite for E F = -2 eV .

One can note that similar behavior is expected to occur in systems with collinear antiferromagnetic structure: depending on the band structure and the band filling, the ferromagnetic state can be close in energy with the AF one, and an applied voltage could stabilize the ferromagnetic state. However in such a case, there is no spin torque since the structure is colinear, and no mechanism for spin reversal (except through spin waves). On the contrary, the mechanism proposed here neither require any defect, spin wave nor a fine tuning of any kind. The noncollinear nature of the magnetic configurations considered provides a robust mechanism for microscopic magnetic switching. Furthermore, this mechanism does not require spin polarized current injection.

In conclusion, we have shown that magnetoresistance of noncollinear frustrated bulk magnets may reach ex-tremely high values and their magnetic configuration may be controlled by applied voltage. The proposed phenomenon is the bulk material analog of spin transfer torque used in layered spin-valve structures. The model chosen here used to prove the concept, in combination with inverse approach of designing the crystal structure from predefined electronic structure properties [START_REF] Franceschetti | [END_REF], opens the way towards experimental realization of the proposed phenomenon.

This work was supported by Chair of excellence Program of the Nanosciences Foundation in Grenoble.

FIG. 1 :

 1 FIG. 1: (a) "4in" and (b) "2in2out" configurations for square lattice

FIG. 2 :

 2 FIG. 2: (a) and (b): total energy as a function of applied voltage for EF =-2 eV and EF =-2.5 eV, respectively, shown by arrows in Fig. 3.

FIG. 3 :

 3 FIG.3: Total density of states for "4in" (black solid) and "2in2out" (red dashed) configurations represented in Fig.1.

FIG. 4 :

 4 FIG.4: Conductance as a function of the Fermi level.