N O Popovych

O V Larina

S M Orlyk
email: orlyk@inphyschem-nas.kiev.ua.

P I Kyriienko

S O Soloviev

S Dzwigaj

DESIGN OF BIFUNCTIONAL CATALYSTS BASED ON BEA ZEOLITES FOR TANDEM PROCESSES WITH PARTICIPATION OF ETHANOL

The results from the creation of BEA zeolite catalysts, containing heteroelements of Ta(Nb) and modified with d-metals (Ag, Cu, Zn), for processes involving ethanol, are presented. The activity of Ag/BEA in the selective reduction of nitrogen oxides (NO x ) with ethanol is due to the presence of acidic and surface redox sites; the promoting effect of hydrogen is produced by the presence of Lewis acid sites and silver nanoclusters. The highest ethanol conversion rates and selective formation of the target product with a high (>60%) yield in the production of 1,3-butadiene from ethanol are achieved in the presence of Cu/Ta/SiBEA. Key words: BEA zeolite; SCR NO-C 2 H 5 OH (+H 2 ); Ag/BEA; Ta, Nb heteroelements; ethanol/acetaldehyde; 1,3-butadiene.

.

Zeolite catalysts are used in such processes because targeted control of their acid-base and redox characteristics is possible even at the synthesis step and particularly for zeolites with the BEA structure by dealumination followed by incorporation of the cations of certain heteroelements into the crystal lattice.

This work presents the results from the design of BEA zeolite catalysts based on determination of the role of the components of the catalytic compositions in revealing the bifunctional surface characteristics, the optimum combination of oxidation-reduction and acid-base sites for increased activity, selectivity, and productivity of the catalysts in processes involving ethanol -the selective reduction of nitrogen oxides (to molecular nitrogen) and the production of 1,3-butadiene.

Ag/BEA ZEOLITES IN THE SELECTIVE CATALYTIC REDUCTION OF NO BY ETHANOL

Compositions based on zeolites of various structural types (MFI, MOR, FAU, BEA) have been proposed for the selective catalytic reduction (SCR) of nitrogen oxides by organic compounds (hydrocarbons, oxygenates). Interest in silver-containing compositions increased significantly after discovery of the promoting effect of hydrogen on the SCR of NO x [5,6]. Considerable attention has been paid in the literature to investigation of the nature of the active sites of such catalysts, and Ag + ions [7][8][9], silver nanoclusters (Ag n d+ , n £ 8) [6], and surface structures of the Ag-O-Al type [10,11] have been examined as active sites for the SCR of NO x [7][8][9]. In order to determine the effect of aluminum in the composition of the Ag/BEA catalysts on the SCR of NO with ethanol zeolite catalysts with various degrees of dealumination, prepared by a two-step post-synthesis method, were investigated [12,13]. The preparation procedure included treatment with an aqueous solution of nitric acid at the first step in order to remove the organic template and aluminum atoms and impregnation of the zeolite with an aqueous solution of silver nitrate at the second step. Samples of zeolite catalysts with silicate modulus (Si/Al) of 12, 100, 200, and 1000 and silver content 1.5-4.0 wt.% were prepared.

The diffractograms of the initial AlBEA (Si/Al = 12), partly dealuminized AlSiBEA (Si/Al = 100 and 200), and dealuminated SiBEA (Si/Al = 1000) are typical of BEA zeolites and do not contain signs of amorphization or formation of new phases. The possibility of dealumination of BEA zeolite arises from the fact that in its structure already contains T atoms (i.e., vacancies of central aluminum atoms in the TO 4 tetrahedra of which the zeolite matrix is composed), and the zeolite lattice is therefore "flexible," which makes it possible to remove Al atoms without loss of crystallinity [13,14].

As a result of decrease in the volume of micropores (237®219 cm 3 /g) and increase in the volume of mesopores (52®55 cm 3 /g) treatment with nitric acid reduces the specific surface area of SiBEA and AlSiBEA somewhat in comparison with AlBEA (625®580 m 2 /g), as is demonstrated by low-temperature adsorption/desorption of nitrogen. With the introduction of silver the specific surface area and volume of the micropores decreases somewhat in line with its content (4.0% Ag/AlBEA -605 m 2 /g, 233 cm 3 /g; 3.0% Ag/SiBEA -567 m 2 /g, 210 cm 3 /g). The decrease of the volume of the mesopores (52®48 cm 3 /g and 55®49 cm 3 /g) may indicate that nanoclusters of silver are formed in the zeolite cavities, thereby blocking access to the adsorbate [12].

Figure 1a shows the IR spectra of the investigated samples in the region of the stretching vibrations of the OH groups. The absorption bands in the spectrum of AlBEA belong to the hydroxyl groups attached to the aluminum (3781 and 3667 cm -1 ), acid Al-O(H)-Si (3615 cm -1 ), and silanol Si-OH (3750 cm -1 with a shoulder at 3740 cm -1 ) groups [15]. After dealumination of the zeolite the absorption bands that belong to the Al-OH and Al-O(H)-Si disappear, indicating removal of the aluminum from the lattice, and new bands appear (3739, 3709, and 3529 cm -1 ), corresponding respectively to the isolated inner, terminal inner, and Si-OH groups with a hydrogen bond, located in the vacancies of the T atoms of the zeolite lattice [16]. For the partly dealuminated zeolites the intensity of the absorption bands of Al-OH and Al-O(H)-Si is reduced, indicating incomplete removal of the aluminum from the zeolite lattice. The dealumination of the zeolites is also confirmed by the X-ray fluorescence method and 27 Al MAS-NMR [14]. With the introduction of Ag into the dealuminated (partly and fully) zeolites the intensity of the bands corresponding to the OH groups in the vacancies of the T atoms decreases, indicating that the silver is located at these positions.

In order to determine the valence state and dispersity of the silver the samples of the zeolite catalysts (after treatment under the conditions of SCR with ethanol in the presence of hydrogen) were studied by diffuse reflectance spectroscopy in the ultraviolet and visible spectral regions (DR UV-Vis spectroscopy) (Fig. 1b). The absorption bands at about 225 nm are due to electronic transitions in the Ag + cations, those at 260-290 nm are due to charge transfer in the nanoclusters Ag n d+ (n £ 8), while the absorption bands at 415 nm indicate the presence of metallic silver nanoparticles in the samples [5][6][7]. The intensity of the band at 225 nm is higher for the Ag/AlSiBEA(100) sample. With increase in the degree of dealumination the intensity of the bands corresponding to the silver nanoclusters decreases, and the intensity of the absorption bands of the silver nanoparticles increases. The presence of highly dispersed oxidized silver in the zeolites was confirmed by XPS on the basis of the proximity of the bond energies of the respective valence electrons.

For analysis of the electronic and coordination state of the active sites of the Ag/BEA zeolites samples with various degrees of dealumination were studied by the Fourier IR spectroscopy with low-temperature adsorption of CO (Fig. 1c). In the case of AlBEA and AlSiBEA the absorption band at 2177 cm -1 is characteristic of CO combined with the bridging hydroxyl groups of Al-O(H)-Si zeolite. The carbonyl band at 2227 cm -1 is characteristic of extra-framework aluminum. The main band in the IR spectrum of SiBEA (2157 cm -1 ) corresponds to CO combined with silanol groups Si-OH in the vacancies of the T atom. As shown by the appearance of absorption bands at 2170-2195 cm -1 , the introduction of silver into the BEA zeolites leads to the formation of new Lewis acid sites. The frequency of the stretching vibrations of Ag + -CO depends on the chemical composition of the support. Thus, the carbonyls on Ag/SiO 2 are characterized by n(CO) at 2169 cm -1 [17], which is close to the value observed for CO adsorbed on the Ag/SiBEA sample. The absorption band at 2186 cm -1 with a shoulder at about 2195 cm -1 in the spectra of Ag/AlBEA and Ag/AlSiBEA corresponds to the symmetrical and asymmetrical vibrations of the Ag + -(CO) 2 dicarbonyls. This indicates at least two coordination vacancies for each Ag + cation, which exhibit the characteristics of Lewis acid sites (LAS).

In order to determine the amount and strength of the acid sites the zeolite samples were studied by IR spectroscopy with pyridine as probe molecule [12,13]. In the IR spectra of pyridine adsorbed on the Ag-containing zeolites there were new absorption bands at 1605 and 1450 cm -1 , which correspond to Lewis acid sites. The intensity of these bands in the spectra of the Ag/AlBEA samples increases with increase of the desorption temperature in comparison with the value in the spectra of the dealuminated (to various degrees) zeolites, which indicates an increase in the strength of the acid sites. The concentration of LAS was calculated from the integral intensities of the absorption bands at 1463, 1455, and 1450 cm -1 and the desorption temperature of the pyridine 150 °C (Fig. 1d). With increase of the Si/Al ratio of the zeolite matrix from 12 to 100 the concentration of Lewis acid sites increases and then decreases with further increase of the ratio. Thus, the silver in the Ag/BEA catalysts makes them bifunctional, forming both redox (Ag + and Ag n d+ ) and additional LAS.

Table 1 shows the values for the conversion of NO during selective reduction with ethanol and the temperature ranges in which they were reached. We note that the temperature dependences of NO conversion in the presence of Ag/BEA zeolite catalysts are described by typical dome-shaped curves [13,18]. A similar degree of transformation of NO was obtained for catalysts with different degrees of dealumination (100 and 200) and an identical Ag content at temperatures up to 350 °C. At higher temperatures the conversion of NO on Ag/AlSiBEA(200), unlike Ag/AlSiBEA(100), decreases. In the presence of samples with a smaller content of Brønsted and Lewis acid sites (Ag/SiBEA(1000) and Ag/AlSiBEA(200)) the rate of oxidation of the ethanol by oxygen evidently increases with increase of temperature, and rate of reduction of NO x decreases accordingly.

The role of the silver sites (Ag + cations and Ag n d+ clusters) involves activation of the alcohol adsorbed on the LAS containing silver by partial oxidation with the formation of more reactive intermediate products (enols) [18][19][20]. Over the strong Lewis acid sites on the surface of Ag/AlBEA the nitrate-nitrite adducts are stabilized over a wider range of temperatures and interact with reducing agent that is activated over the silver sites. For the dealuminated samples with a smaller content of LAS comparable rates of formation of the intermediates and their interaction are possible over a narrower range of temperatures.

The role of the LAS in the Ag/BEA zeolites was confirmed by investigating the catalysts by DRIFT spectroscopy in situ in a stream of NO + O 2 [18]. Bands corresponding to bridged, bidentate, and monodentate nitrates were identified in the spectra of the Ag/AlBEA (12) and Ag/SiBEA(1000) samples, and their intensities were twice as high as those for the first sample. When the temperature was increased from 150 to 400 °C a decrease in the intensity of the bands was observed for both samples. At 400-500 °C, however, the intensities of the bands for Ag/AlBEA (12) hardly change at all, in contrast to the significant decrease for Ag/SiBEA(1000). Thus, the nitrate adducts over the strong Lewis acid sites of Ag/AlBEA (12) are stabilized over a wider range of temperatures. The effect of hydrogen on the selective reduction of NO in the presence of silver-containing BEA zeolites depends on the Si/Al ratio of the samples. It should be noted that hydrogen is not a selective reducing agent for NO. (Reduction is not observed when there is no ethanol in the reaction mixture.) Figure 1d shows the increases of the conversion of NO after addition of hydrogen for samples with various Si/Al ratios. The addition of hydrogen to the reaction mixture (0.05% NO + 0.1% C 2 H 5 OH + 10.0% O 2 + 0.5% H 2 in He) increases the conversion of NO at 400 °C on the Ag/AlSiBEA(100) samples: from 24% to 50% with a 1.5% Ag content and from 35% to 56% in the case of 3.0% Ag. The effect of hydrogen is less significant in the presence of the catalysts based on AlBEA (12) and AlSiBEA(200) and is practically absent for Ag/SiBEA(1000). The promoting effect of H 2 is reversible; after removal of hydrogen from the gas mixture the conversion of NO is reduced to the value obtained in the absence of H 2 . The effect of H 2 on the SCR process in the presence of Ag/BEA varies in line with the change in the content of Lewis acid sites in the zeolites.

In [21] it was shown that activation of H 2 in the silver clusters deposited on q-Al 2 O 3 (110) occurs in the region of contact between the cluster and the Lewis acid-base pair (Al-O). After dissociation of the hydrogen during interaction with the oxygen hydroxyl and hydroperoxide particles, which are more active in the partial oxidation of the reducing agent, are formed [22]. The insignificant promoting effect of H 2 for the dealuminated Ag/SiBEA zeolites may be due to the fact that the small surface concentration of the nitrate adducts levels out the increase in the rate of activation of the reducing agent in the presence of hydrogen. This is confirmed by the composition of the products from SCR of NO by ethanol (+H 2 ): in the presence of Ag/SiBEA(1000) the content of CO and acetaldehyde (the products of partial oxidation of the reducing agent) is increased after the addition of hydrogen to the reaction mixture, while the degree of transformation of the NO remains unchanged. A larger content of Lewis acid sites on the surface of Ag/AlSiBEA(100) leads to the appearance of a promoting effect from the hydrogen. It can therefore be supposed that the presence of the LAS on the surface of the catalyst is essential for the appearance of the H 2 effect in the SCR of NO.

The role of the silver thus involves the formation of the oxidation-reduction and Lewis acid sites over which the reagents are activated. The wider temperature range for the SCR of NO by ethanol in the presence of Ag/AlBEA compared with the dealuminated (to various degrees) zeolites is due to the high concentration of strong Lewis acid sites, over which the nitrate adducts capable of interaction with the reducing agent activated over the silver-containing sites are stabilized. The promoting effect of the hydrogen on the selective reduction of the NO is ensured by the presence of Lewis acid sites and silver nanoclusters on the surface of the catalysts.

MODIFICATION OF BEA ZEOLITES FOR THE PRODUCTION OF 1,3-BUTADIENE FROM ETHANOL

Ethanol is a promising renewable feedstock for a series of valuable chemical products. One of them is 1,3-butadiene (BD), which is widely used in the production of polymeric materials [23]. It can be obtained from ethanol (the Lebedev method) and by conversion of a mixture of ethanol with acetaldehyde (AA) (the Ostromyslenskii method). Tantalum-silicate compositions Ta 2 O 5 -SiO 2 are catalysts for the production of BD by conversion of a mixture of ethanol and AA [24,25], in which the active sites are Lewis acid sites containing tantalum.

The use of BEA zeolite catalysts for the production of BD arises from the possibility of targeted control of their acid-base characteristics even at the synthesis step and in particular by dealumination followed by incorporation of certain heteroelements into the crystal lattice. The heteroelement can be present in the zeolite structure in two forms: built into the lattice (highly dispersed in the M n+ state) and outside the lattice (in the form of clusters and MO x nanoparticles) [26]. In the zeolite materials with the incorporated heteroelements it is possible to achieve high selectivity in the production of BD on account of the formation of a large number of active sites in one of the key steps of the process -aldol condensation of the AA.

The dealuminated BEA zeolites with incorporated tantalum atoms were prepared by a post-synthesis method [27,28] and were characterized by XRD, IR spectroscopy, and DR UV-Vis spectroscopy. As a result of dealumination the Brønsted acid sites [Al-O(H)-Si], which are active in the side reactions of dehydration of ethanol to diethyl ether (DEE) and ethylene, are removed from the zeolite. Here, during subsequent incorporation of the heteroelement (Ta) Lewis acid [Ta + (-O-Siº) 4 ] and basic [O -Ta(-O-Siº) 4 ] sites, which are active in the conversion of the ethanol/acetaldehyde mixture into 1,3-butadiene, are formed on the surface of the dehydrated Ta/SiBEA zeolite.

Samples of 1(3)% Ta/SiBEA were studied in the conversion of ethanol and its mixture with AA. The zeolites are not very active in the conversion of ethanol, and the main products are AA, ethylene, and DEE. The selectivity in the formation of the targeted product is no higher than 16%, and the yield of BD is low (~3%) (Fig. 2a). Increase of the tantalum content of the zeolite catalyst from 1 to 3 wt.% does not lead to any substantial increase in the yield of BD. At the same time, it is possible to achieve high selectivity in the conversion of the ethanol/AA mixture. The highest selectivity with respect to BD amounts to 87%-90% with 23%-28% conversion and AA fractional content of 0.31 in the mixture with ethanol, which is substantially higher than the values described in the literature [24].

The selectivity of formation of BD in the presence of Ta/SiBEA increases with increase of the AA content in the reaction mixture, and here the selectivity with respect to ethylene and DEE decreases. The amount of the products probably decreases on account of retardation of the side processes involving dehydration of the ethanol. The aldol condensation of AA and the dehydration of ethanol to ethylene and DEE occur over identical active sites (Lewis acid sites) and consequently compete. Since Meerwein-Ponndorf-Verley (MPV) reduction of crotonaldehyde (CA) also occurs over acid sites increase of the AA content of the initial reaction mixture can be expected to lead to a larger amount of CA in the reaction products. Crotyl alcohol is not found in the products of conversion both of ethanol and of the ethanol/AA mixture, and it can therefore be stated that the targeted process is not limited by dehydration of the crotyl alcohol.

The high selectivity and the yield of the product (butadiene) are not changed after 7-8 h continuous use of the Ta/SiBEA zeolite catalysts [29]. It should be noted that the crystal structure of the tantalum-containing catalysts is unchanged even after repeated cycles of use (total ~400 h) and subsequent regeneration (in a stream of 5% O 2 in Ar, 500 °C, 2 h).

As mentioned above, the selectivity of formation of BD in the conversion of ethanol in the presence of the Ta/SiBEA catalyst is less than 20%, which is due to the insufficient reactivity of the active sites of the zeolite compositions for the key step of the process -the dehydrogenation of ethanol to AA. It is known that doping of oxide compositions with the compounds of d metals and of silver, copper, and zinc in particular intensifies the dehydrogenation of ethanol [30][31][32][33].

To form surface active sites for dehydrogenation of ethanol to acetaldehyde the Ta/SiBEA zeolite catalyst was therefore doped with d metals (Ag, Cu, Zn), and this led to a substantial increase of the conversion of ethanol and the selectivity with respect to BD (Fig. 2c,d). The highest values for these properties (88% and 73% respectively) were observed in the presence of Cu/Ta/SiBEA zeolite. The selectivity of formation of BD at Ag/Ta/SiBEA was also fairly high (63% at 83% conversion of ethanol). In the presence of the Zn/Ta/SiBEA sample the selectivity for BD amounts to 43% at 52% conversion of ethanol. Here significant amounts of AA (selectivity of 23%), ethylene, and DEE (overall selectivity of 30%) are observed in the reaction products, indicating an insufficient amount of acid sites of medium strength for the aldol condensation and a significant amount of strong-acid sites that are active in the ethanol dehydrogenation side reactions.

Increase of the load on the Cu/Ta/SiBEA catalyst (ethanol weight hourly space velocity increased by three times -to 1.5 h -1 ) leads to a substantial decrease in the conversion of ethanol, the selectivity of formation of the targeted product, and the productivity of the catalyst with respect to BD -from 0.19 to 0.12 g/(g cat •h) [34]. Here the amount of dehydration products (ethylene and DEE) and acetaldehyde in the reaction products is increased significantly due to the kinetics of the multistage process.

In order to determine the role of the additions of d metals (Ag, Cu, Zn) in the formation of the active sites of the Ta/SiBEA zeolites the modified samples were investigated by IR spectroscopy with the use of probe molecules (pyridine and deuterochloroform, CDCl 3 ).

The IR spectra of pyridine adsorbed on the surface of the modified Ta/SiBEA samples contain bands that indicate the presence of Lewis acid sites (absorption bands at 1613, 1599, 1492, 1454, and 1446 cm -1 ) (Fig. 3a) [14,[35][36][37]. The absorption bands of pyridine coordinatively bound to the Lewis acid sites are significantly weaker for Ta/SiBEA. It can therefore be assumed that additional acid sites are formed in the matrix of the dealuminated BEA zeolite on account of entry of the silver, copper, or zinc ions into the vacancies of the T atoms. On the basis of the intensities of the bands in the IR spectra of the Fig. 3. IR spectra of adsorbed pyridine (a) and deuterochloroform (b) on samples 1% Ta/SiBEA (1), 1% Ag/1% Ta/SiBEA (2), 1% Cu/1% Ta/SiBEA (3), and 1% Zn/1% Ta/SiBEA (4) after vacuum treatment at 250 °C (pyridine) and 50 °C (CDCl 3 ). modified (Ag, Cu, Zn) Ta/SiBEA zeolites after desorption of the pyridine at 350 °C a larger amount of strong Lewis acid sites is observed for the zinc-containing sample. Weak LAS are mostly present in the Ag/Ta/SiBEA zeolite.

In the IR spectrum of deuterochloroform adsorbed on the surface of the Ta/SiBEA sample there are bands at 2558 and 2556 cm -1 , which indicate the presence of weak basic sites (Fig. 3b) [38]. In the IR spectra of the modified samples the intensity of these bands is significantly lower. Moreover, the absorption band of CDCl 3 adsorbed on the basic sites of medium strength (~2220 cm -1 ) is absent from the spectra, indicating that the introduction of silver, copper, or zinc into the composition of Ta/SiBEA leads to blockage of the basic sites of medium strength with a simultaneous decrease in the amount of weak sites.

Comparison of the results of IR spectroscopy of the modified zeolites with the characteristics of the targeted process of production of BD from ethanol in their presence indicates that the incorporation of silver, copper, and zinc into the Ta/SiBEA zeolite leads to the formation of additional LAS that are active in the aldol condensation of AA. Here a large concentration of strong Lewis acid sites, which accelerate the side reaction of dehydration of ethanol to ethylene, is observed for the copper-and zinc-containing samples. With the addition of silver predominantly weak Lewis acid sites, less active in the targeted process, are formed.

According to the results from DR UV-Vis spectroscopy, XPS, and TEM methods there are clusters and subnanoparticles of the metals (Ag, Cu, Zn) on the surface of the modified Ta/SiBEA compositions. Their stabilization in a highly dispersed state even after catalysis is promoted by the silicate matrix of the zeolite. Their role mainly includes the formation of active sites for dehydrogenation of ethanol to AA at a higher rate than over the basic sites of Ta/SiBEA. It can be supposed that the dehydrogenation of ethanol occurs over active sites formed with the participation of silver, copper, and zinc, while further transformation of the ethanol/AA mixture by aldol condensation, MPV reduction, and dehydration of the formed CA take place over the tantalum-containing acid sites of the zeolite. The production of BD from ethanol in the presence of M/Ta/SiBEA (M: Ag, Cu, Zn) can be represented schematically in the following way:

The effect of the nature of the niobium particles in the Nb/SiBEA zeolites on their catalytic characteristics in tandem processes and in the production of 1,3-butadiene from ethanol in particular was established in [26]. It was shown that the active sites formed with the participation of the various states of niobium in the structure of the Nb/SiBEA zeolite accelerate various steps of the production of BD from ethanol. Incorporation of the niobium into the crystal lattice of dealuminated BEA zeolite at the vacancies of the T atoms is confirmed by XRD, NMR, IR and DR UV-Vis spectroscopy; the ratio of niobium in the various states depends on its (Nb) amount added to the sample. According to the results from DR UV-Vis spectroscopy the Nb/SiBEA samples (0.7 and 2 wt.% Nb) contain highly dispersed niobium in the Nb 5+ state; the first sample contains a larger amount of incorporated Nb(V), whereas the second sample contains extra-framework Nb(V) and NbO x nanoparticles.

Table 2 gives data on the catalytic activity (selectivity) of the Nb/SiBEA zeolites in the conversion of ethanol and its mixture with AA. The main products from the conversion of ethanol are BD, AA, ethylene, and DEE, and small amounts (<1%) of crotonaldehyde and butene isomers were also detected. The 0.7% Nb/SiBEA catalyst is characterized by a high ethanol conversion rate and selectivity for the formation of BD, as a result of which the yield of BD is almost twice the yield for the 2.0% Nb/SiBEA sample. The samples also differ with respect to the dehydrogenation (AA) and dehydration (ethylene and DEE) products. The higher selectivity for the formation of AA on 2.0% Nb/SiBEA is probably due to the smaller number of aldol condensation sites on its surface for subsequent transformation of AA to BD. The higher selectivity of the 0.7% Nb/SiBEA sample for ethylene and DEE indicates a larger number of acid sites over which the ethanol dehydration side reaction occurs.

The transformation of the ethanol/AA mixture gives a higher degree of conversion and a higher yield of the targeted product (BD) in the presence of 0.7% Nb/SiBEA. On this basis it can be supposed that the aldol condensation takes place with a higher rate over the sites formed with participation of the Nb(V) incorporated into the crystal lattice of the BEA. At the same time, substantially higher selectivity for the formation of BD (70.6%) is obtained on the 2.0% Nb/SiBEA catalyst than on 0.7% Nb/SiBEA (54.8%).

The observed increase of the amount of ethylene (from 11.9% to 29.5%) in the products of conversion of the ethanol/AA mixture in the presence of the 2.0% Nb/SiBEA catalyst and the almost identical value of the selectivity for the formation of BD for the sample with 0.7% of Nb with increase of the temperature of the process from 325 to 350 °C give reason to suppose that the ethanol dehydration side reaction can occur over acid sites formed by extra-framework Nb(V) and NbO x .

The products of transformation of the mixture of ethanol with AA in the presence of 2.0% Nb/SiBEA contain a significantly larger amount of CA (10.5%-11.8%). The considerable difference in the selectivity of formation of CA may be due to the fact that its reduction by ethanol becomes a key factor in the production of BD.

To assess the effect of the location of the niobium in the BEA zeolite matrix on the activity of the Nb/SiBEA samples in the MPV reduction of CA by ethanol their activity in the transformation of the ethanol/CA mixture was tested [26]. A higher degree of conversion and a larger amount of BD produced (58.2%) were observed on the 0.7% Nb/SiBEA sample, and it can therefore be assumed that the MPV reduction takes place at a higher rate over the sites formed with participation of incorporated Nb(V).

The specific rate of formation of 1,3-butadiene (turnover frequency (TOF)) in the transformation processes of ethanol and ethanol/AA mixture over the 0.7% Nb/SiBEA catalyst, which only contains Nb(V) incorporated into the crystal lattice of BEA zeolite, is several times higher than the TOF for the 2.0% Nb/SiBEA sample with extra-framework Nb(V) and nanodispersed NbO x (Table 2). The rate of dehydrogenation of ethanol in the presence of 0.7% Nb/SiBEA is almost four times higher than for 2.0% Nb/SiBEA. The 0.7% Nb/SiBEA sample is also more active in the side reactions of dehydration of ethanol to ethylene and DEE, which take place over the acid sites on the surface of the catalyst.

Thus, the dehydrogenation of ethanol, the aldol condensation of AA, and its MPV reduction by ethanol, and likewise the tandem process for the production of BD from ethanol, take place more effectively over the Lewis acid sites formed by Nb(V) incorporated into the BEA zeolite matrix than over the sites formed with participation of extra-framework Nb(V) and NbO x nanoparticles.

To summarize the results of the investigations into the bifunctional nature of the active sites on the surface of zeolite catalysts with the BEA structure we draw attention to the following points. The appearance of SCR activity in Ag/BEA catalysts and the promoting effect of hydrogen on the selective reduction of nitrogen oxides (NO x ) by ethanol arises from the bifunctional surface characteristics, i.e., the presence of acid sites and redox sites. Effective catalysts for the production of 1,3-butadiene from ethanol are proposed as a result of investigation of the nature of the effect of the active components and modifying additives of Cu(Ag,Zn)/Ta/SiBEA zeolite compositions on the acid-base characteristics of the surface; the highest degrees of conversion and selectivity for the formation of the targeted product (butadiene) with high yields (>60%) are obtained in the presence of Cu/Ta/SiBEA. Incorporation of niobium into the structure of the BEA zeolite leads to the formation of a larger amount of Lewis acid sites per mole of Nb. This makes it possible to produce BD from ethanol (ethanol/acetaldehyde mixture) in the presence of catalysts with smallest content of active component.

Fig. 1 .

 1 Fig. 1. Characteristics of Ag/AlSiBEA samples: IR spectra (a), UV spectra (b), IR spectra of CO adsorbed on the samples at -173 °C and 10 -3 Pa (c); the concentration of Lewis acid sites (IRS, pyridine, 150 °C) and the increase of the conversion of NO (DX, %) during SCR with ethanol in the presence of H 2 (d).

Fig. 2 .

 2 Fig. 2. Characteristics of the conversion of ethanol into 1,3-butadiene in the presence of zeolite compositions 3% Ta/SiBEA (a) (T = 325 °C, W = 0.8 h -1 ) and 1% Ag(Cu,Zn)1% Ta/SiBEA (b, c) (T = 325 °C, W = 0.5 h -1 ).

TABLE 1 .

 1 Catalytic Activity of Ag-Containing BEA Zeolites in SCR of NO with Ethanol (0.05% NO + 0.1% C 2 H 5 OH + 10.0% O 2 (+ 0.5% H 2 ) in He; 24000 h -1 )

		Conversion of NO, %/attainment temperature, o C
	Ag/BEA(Si/Al) catalyst
		NO + C 2 H 5 OH + O 2	NO + C 2 H 5 OH + O 2 + H 2
	1.5% Ag/AlBEA(12)	26/280-380	32/280-370
	4.0% Ag/AlBEA(12)	34/300-440	50/350-385
	1.5% Ag/AlSiBEA(100)	24/400-450	48/390-420
	3.0% Ag/AlSiBEA(100)	35/400-500	57/400-435
	1.5% Ag/AlSiBEA(200)	23/340-380	27/355-400
	3.0% Ag/AlSiBEA(200)	26/350-390	39/365-410
	1.5% Ag/SiBEA(1000)	27/310-350	27/310-350
	3.0% Ag/SiBEA(1000)	31/320-360	34/310-335

TABLE 2 .

 2 Catalytic Characteristics of Nb/SiBEA Zeolite Catalysts in the Conversion of Ethanol and the Ethanol/Acetaldehyde Mixture (T = 325 °C, W = 0.8 h -1 )

					Selectivity, %			Specific rate, mol/(mol Nb•h)
	Catalyst	Conversion, %	BD	AA Ethylene DEE CA Others	Transformations of ethanol (ethanol/AA)	Formation of Ethylene + BD DEE
						Ethanol
	0.7% Nb/SiBEA	42.4	27.9 16.8	42.0	12.2 0.1	1.0	98	27	53
	2% Nb/SiBEA	27.1	26.0 29.6	28.0	15.6 0.1	0.7	22	6	10
					Ethanol/AA = 2.7
	0.7% Nb/SiBEA	38.5	54.8	-	36.8	2.4	2.5	3.5	89	49	35
	2% Nb/SiBEA	22.1	70.6	-	11.9	2.0 11.8	3.7	18	13


É Lichtfouse

B Elbisser

J Balesdent

A Mariotti

G Bardoux

Isotope and molecular evidence for direct input of maize leaf wax n-alkanes into crop soils

Keywords: n-alkane, soil, plant, waxes, Zea mays, 13 C

The contribution of plant carbon to crop soils can be followed by isotope labelling at natural abundance, such as growing a C 4 plant on a soil which was previously under C 3 vegetation. For this purpose, carbon isotope compositions and relative abundances of n-alkanes of maize leaf waxes and maize crop soils were compared. Isotope values of soil n-alkanes increased with time of maize cultivation as the result of maize carbon integration into soil organic matter. With increasing time of cultivation, the increase in isotopic difference between n-heptacosane (C 27 ) and n-nonacosane (C 29 ) is explained, at least partly, by a direct input of maize leaf n-alkanes. The amount of maize-derived carbon within each n-alkane has been calculated by isotopic means.

INTRODUCTION

Distributions of sedimentary n-alkanes have been widely used for source identification in sediments. However, the specificity of n-alkanes as source parameters might be limited because our knowledge of biological n-alkanes, or their possible precursors, is relatively poor, especially among microorganisms, as reviewed by Collister et al., 1994a. For example, the distribution of n-alkenes from the alga Botryococcus Braunii (A race) with a strong predominance of oddcarbon number homologues in the range C 27 to C 33 , is very close to n-alkane distributions in terrestrial plants (Lichtfouse et al., 1994a). Moreover, recent isotope studies have shown that extracted sedimentary n-alkanes can derive not only from higher plants n-alkanes but also from n-alkanoic acids [START_REF] Lichtfouse | Tracing biogenic links of natural organic substances at the molecular level with stable carbon isotopes: n-alkanes and n-alkanoic acids from sediments[END_REF] or from algal n-alkenes [START_REF] Freeman | Predictive isotopic biogeochemistry: hydrocarbons from anoxic marine basins[END_REF], Collister et al., 1994a[START_REF] Spooner | Stable carbon isotopic correlation of individual biolipids in aquatic organisms and a lake bottom sediment[END_REF], Lichtfouse et al., 1994a).

In soils, studies of possible sources of organic components seem more straightforward, as compared to sediments, because the input of phytoplankton is minor relative to the input of higher plants. However, numerous soil organisms are able to biosynthesize n-alkanes or their precursors (Kolattukudy, 1976[START_REF] Dinel | Soil lipids: origin, nature, content, decomposition, and effect on soil physical properties[END_REF]. For this reason, and also because food chain flows from plants to soil biomass are far from being understood, the direct input of leaf wax nalkanes into soil is still to be proved. The fate of plant carbon into a soil can be followed by isotopic labelling in natural conditions. Hence, growing a C 4 plant, maize (bulk δ 13 C -12‰) on a soil which was previously under C 3 vegetation (bulk δ 13 C -26‰) can cause the isotopic composition of bulk soil organic matter to increase from -26‰ toward -13‰ [START_REF] Balesdent | Natural 13 C abundance as a tracer for studies of soil organic matter dynamics[END_REF]. Because leaf wax n-alkanes from C 3 plants and C 4 plants are also isotopically different [START_REF] Rieley | Sources of sedimentary lipids deduced from stable carbon-isotope analyses of individual compounds[END_REF][START_REF] Rieley | Gas chromatography-Isotope Ratio Mass Spectrometry of leaf wax n-alkanes from plants of differing carbon dioxide metabolisms[END_REF], Collister et al., 1994b, Lichtfouse et al., 1994a), the isotopic composition of soil n-alkanes should be related to their source. In this preliminary report, we have evaluated the contribution of leaf wax n-alkanes to soil n-alkanes via molecular and isotope analysis. EXPERIMENTAL Detailed procedures will be described elsewhere [START_REF] Lichtfouse | [END_REF]. Zea mays was grown in an experimental field at La Minière, France, at one crop per year on a soil which has previously been cultivated only with C 3 plants. Soils and Z. mays leaves were extracted with CHCl 3 -MeOH 3/1 v/v. Alkane fractions were obtained from extracts by silica gel chromatography (Lichtfouse et al., 1994b). n-Alkanes were identified by gas chromatographymass spectrometry (GC-MS) and by co-elution of pure standards. Their purity (>99.9%) was carefully checked by GC-MS and gas chromatography-combustion-isotope ratio mass spectrometry (GC-C-IRMS) in order to give accurate isotopic values [START_REF] Lichtfouse | Enhanced resolution of organic compounds from sediments by isotopic gas chromatography-combustion-mass spectrometry[END_REF]. The overall deviation was ± 0.3‰. n-Alkane relative abundances (area/area) were determined by gas chromatography using a flame ionization detector. GC-C-IRMS conditions are given in Lichtfouse et al. (1994a). Carbon isotopic compositions are expressed in per mil. relative to the PDB standard: δ 13 C = [( 13 C/ 12 C sample -13 C/ 12 C std)/( 13 C/ 12 C std)] x 10 3 , where 13 C/ 12 C std = 0.0112372.

RESULTS AND DISCUSSION

Initial isotope ratios

Carbon isotope compositions and relative abundances of selected n-alkanes from maize leaf waxes and maize crop soils are reported in Table 1. Before maize cropping, the soil was cultivated with C 3 plants. Distributions of soil n-alkanes in the range C 25 to C 35 were very similar to those of higher plants, e.g. they display a strong odd-carbon number predominance [START_REF] Lichtfouse | [END_REF]. Before cultivation, the isotope compositions for n-heptacosane (C 27 ) and n-nonacosane (C 29 ) were -34.2‰ and -35.7‰, respectively (Table 1, Figure 1). These values are well within the range of those reported for C 3 plants [START_REF] Rieley | Sources of sedimentary lipids deduced from stable carbon-isotope analyses of individual compounds[END_REF][START_REF] Rieley | Gas chromatography-Isotope Ratio Mass Spectrometry of leaf wax n-alkanes from plants of differing carbon dioxide metabolisms[END_REF], Collister et al., 1994b, Lichtfouse et al., 1994a). Also, observation of distributions together with isotope values strongly suggests a direct contribution of C 3 plants waxes from previous C 3 vegetation, although a contribution from soil organisms which biosynthesize those n-alkanes on C 3 carbon is not excluded.

Isotopic variations with increasing time of cultivation

Isotope compositions of n-alkanes increased with cultivation time (Figure 1). Such variations are explained by the input of 13 C enriched carbon from maize into soil organic components. A direct input of maize leaf waxes into soil is possible since isotope compositions of maize n-alkanes are 13 C enriched: -19.1 ‰ for n-heptacosane (C 27 ) and -18.4‰ for nnonacosane (C 29 ). Also the progressive addition of maize n-alkanes during cultivation should increase the isotope values of the soil n-alkanes, as shown in Fig. 1.

From the isotopic trends, it is obvious that the amount of maize-derived carbon in n-heptacosane (C 27 ) increased faster than in n-nonacosane (C 29 ). This difference can be quantitated by two means. Firstly, it can be assumed that isotopic trends versus cultivation time follow a linear law on a short time period, though this is probably not true for longer experiments. The slope calculated by linear regression, which is representative of the isotopic increase versus cultivation time for each component, is indeed much higher for n-heptacosane (1.028 ‰.yr -1 ) than for n-nonacosane (0.300 ‰.yr -1 ). Secondly, each soil n-alkane can be regarded as a mixture of two distinct isotopic sources: "C 3 " n-alkane from the initial soil and "C 4 " n-alkane from maize waxes. At a given time of cultivation, the molar percentage of maizewax derived carbon M in each soil n-alkane can be calculated by M = 100.(δδ 0 )/(δ mδ 0 ) δ: isotope composition of soil n-alkane at a given time of cultivation. δ 0 : isotope composition of soil n-alkane before maize cultivation (time = 0). δ m : isotope composition of n-alkane for maize waxes. Calculated M values after 5 crops of maize (4.33 yrs.) are 31.1‰ for n-heptacosane and 8.1‰ for n-nonacosane. Those values show again that n-heptacosane incorporated more maizederived carbon than n-nonacosane, probably via direct input of leaf waxes. Indeed, this hypothesis is strengthened by consideration of n-alkane relative abundances, as discussed below.

n-Alkane ratios in soils and maize

In order to get more insight on the mode of integration of maize carbon into crop soils, n-alkane relative abundances from maize leaves and soils were compared using C 27 /C 29 ratios (Table1). Maize waxes had a much higher C 27 /C 29 ratio (0.916) than did the soil before maize cultivation (0.321). It is thus expected that the C 27 /C 29 ratio would increase with time of cultivation if a direct input were happening. C 27 /C 29 ratios indeed showed a regular increase from 0.321 to 0.411 during 5 years of maize cropping (Fig. 2), thus demonstrating that soil n-alkanes were derived, at least partly, by direct input of maize leaf waxes. Further evidence is given by comparison of C 27 /C 29 ratios with δ 27 -δ 29 in soils (Figure 3). The correlation observed confirms that variations of isotopic differences between n-heptacosane and n-nonacosane resulted, at least partially, from a higher input of n-heptacosane from maize leaf waxes.

CONCLUSION

The contribution of plant waxes to crop soils has been evaluated by in situ labelling under natural conditions. A simple calculation has been used to evaluate the contribution of isotopically distinct sources within each n-alkane. Futhermore, soil n-alkanes have been shown to derive, at least partly, from the direct input of plant leaf n-alkanes. Isotope studies of other nalkanes and organic fractions will be published elsewhere [START_REF] Lichtfouse | [END_REF]. Further work on soil alcohols and fatty acids is in progress.

Fig. 1 .

 1 Fig. 1. Carbon isotope composition of soil n-alkanes versus time of maize cultivation. Note the higher isotopic increase of n-heptacosane (C 27 ) versus n-nonacosane (C 29 ). Linear regression gives δ 27 = 1.028.time -33.66 (r = 0.96) and δ 29 = 0.300.time -35.74 (r = 0.94).

Fig. 2 .

 2 Fig. 2. n-Heptacosane (C 27 )/n-nonacosane (C 29 ) ratios in crop soils versus time of maize cultivation. Linear regression gives C 27 /C 29 = 0.021.time + 0.319 (r = 0.97).

Fig. 3 .

 3 Fig. 3. Abundance ratios versus isotopic difference for n-heptacosane (C 27 ) and n-nonacosane (C 29 ) in maize crop soils. Linear regression gives δ 27δ 29 = 0.026.C 27 /C 29 + 0.269 (r = 0.93).

Table 1 .

 1 Carbon isotope composition and relative abundance of selected n-alkanes from crop soils and maize leaf waxes.

	_________________________________
	Cultivation Carbon number C 27 /C 29
	time (years) 27	29
	0 *	-34.2 -35.7	0.321
	0.42	-32.5 -35.3	0.332
	1.00	-32.7 -35.5	0.335
	1.42	-32.1 -35.4	0.355
	1.75	-32.3 -35.5	0.338
	2.42	-31.1 -35.1	0.378
	3.42	-29.5 -34.7	0.389
	4.00	-29.7 -34.5	0.403
	4.33	-29.5 -34.3	0.411
	_________________________________
	Maize	-19.1 -18.4	0.916
	_________________________________

* Uncultivated soil previously under C 3 vegetation.

Acknowledgements-We thank Profs. John hayes from Indiana University and Geoffrey Eglinton from Bristol University for helpful discussions. Dr. May Balabane from the Institut National de la Recherche Agronomique at Versailles is gratefully acknowledged for sampling assistance. We are indebted to Prof. James Maxwell from Bristol University for a critical review of this manuscript.


Ashesh Shah

Jérôme Coste

Jean-Jacques Lemaire

Ethan Taub

W M Michael Schüpbach

Claudio Pollo

Raphael Guzman

Karin Wårdell

Erik Schkommodau

Simone Hemm-Ode

Intraoperative Quantitative Tremor Evaluation in Deep Brain Stimulation Surgery

Keywords: Deep brain stimulation, Intraoperative monitoring, Acceleration, Tremor, Parkinson's disease, Essential tremor

Deep brain stimulation (DBS) is a common neurosurgical treatment for the tremor of Parkinson's disease, essential tremor, and tremor of other causes. The outcome depends on optimal placement of the permanently implanted electrode. Many centers perform DBS surgery under local anesthesia in order to confirm the therapeutic effect with intraoperative stimulation testing. Visual inspection-the method generally used to rate changes in tremor during stimulation testingis subjective, and its accuracy depends on the evaluator's experience [START_REF] Griffiths | Automated assessment of bradykinesia and dyskinesia in Parkinson's disease[END_REF]. This study presents the results of quantitatively estimating improvement in tremor during intraoperative stimulation tests in 15 patients. In addition, its influence on identifying the final position of the permanently implanted electrode is described.

Method:

We designed a 3D acceleration sensor system that is attached to the patient's forearm during surgery [START_REF] Shah | A method to quantitatively evaluate changes in tremor during deep brain stimulation surgery[END_REF]. During intraoperative stimulation tests, at each different position, accelerometric data are synchronously recorded with the changing stimulation current amplitude. The method was applied in 15 DBS procedures in 2 centers (University Hospital Bern, Switzerland & University Hospital Clermont-Ferrand, France); the data were analyzed offline to assess improvements in tremor and to identify tremor-suppressing stimulation current-amplitudes. For correlation analysis, the quantitatively and visually determined improvements in tremor were categorized into: no improvement, low improvement, average improvement, high improvement and tremor arrest. The quantitatively identified tremor-suppressing current amplitudes were compared to those identified by visual inspection, in order to determine the influence these findings would have had on the position chosen for permanent electrode implantation if they had been used for intraoperative decision-making. As this was a purely observational study, the accelerometric measurements were not, in fact, allowed to alter the surgical procedures in any way.

Results:

A total of 359 evaluations were available for a comparison of the improvement in tremor identified by accelerometry vs. visual inspection. Of these evaluations, 156 (43.5%) were assigned the same categorie by both methods; 296 (82.5%) fell in the same or neighboring categories; and 63 (17.5%) were at least 2 categories apart. The quantitatively identified tremorsuppressing current-amplitudes were significantly lower than the visually identified ones (1.13 ± 0.8 mA vs. 1.7 ± 0.8 mA [mean ± SD]). Of the 26 finally chosen positions for permanent lead implantation, 15 would have been different had the accelerometric data been considered.

Discussion and Conclusion:

The improvement of tremor brought about by test stimulation was rated in the same category by visual inspection and by quantitative measurement (accelerometry) in only 43.5% of the evaluations that we made in this study. In some of the evaluations where there was only mild tremor at baseline, the stimulation-induced improvement in tremor was classified by visual inspection as tremor arrest, while accelerometry revealed a very mild residual tremor. This fact explains many of the instances in which the tremor ratings obtained by the two methods were only 1 category apart, but it cannot account for the 17.5% of evaluations that were 2 or more categories apart. The quantitative assessment of tremor as performed here yields different findings from assessment by visual inspection alone. The tremor-suppressing stimulation current-amplitudes are lower, and this, in turn, can often affect the chosen site for permanent electrode implantation. Thus, quantitative tremor assessment can affect, and perhaps improve, targeting in DBS without altering the routine surgical procedure. This tentative conclusion awaits confirmation by further studies. Moreover, aside from its potential direct clinical utility, quantitative tremor assessment DBS surgery might be a useful adjunct to the clinical testing of new types of DBS electrode


Elasticity: Thermodynamic Treatment

The elastic behavior of materials is usually described by direct stress-strain relationships which use the fact that there is a one-to-one connection between the strain and stress tensors ε and σ; some mechanical elastic potentials can then be defined so as to express the stress tensor as the derivative of the strain potential with respect to the strain tensor, or conversely the latter as the derivative of the complementary (stress) potential with respect to the former.

Nevertheless, these potentials are restricted to specific thermodynamic conditions. This purely mechanical treatment may be generalized and improved by integrating its temperature dependence and by connecting the elastic potentials with the classical thermodynamic functions within a consistent thermodynamic framework. This leads to a more general definition of thermoelasticity (Sect. 1), which can be used for a better understanding of the properties of the elastic moduli (Sect. 2), for a natural definition of rubber elasticity (Sect. 3), and for the prediction of the effective thermal expansion coefficients of heterogeneous materials (Sect. 4).

Thermodynamic Definition of Thermoelasticity

General Case

The thermodynamic framework relies on the two fundamental (conservation) classical principles, for which local expressions can be derived from the global ones (Germain et al. 1983). Infinitesimal strains as well as quasi-static evolutions only are considered. For any given body, V, subjected to purely thermodynamic processes, the first principle expresses the fact that, whatever the mechanical behavior, the variation of the total (i.e., internal plus kinetic) energy is equilibrated by the variation of the mechanical work of the applied forces and by the variation of the heat received by the body.

When thermal conduction only occurs through the boundary cV, with q the heat flux vector and n the outward unit normal, and with r denoting the specific heat production rate, the global heat rate, Q, is given by Q l V r dVk c V q :n dS l V (rkdiv (q )) dV. Taking into account the equilibrium equations, which are used to connect the stress field with the exterior loading, the local expression of the first principle finally reads ρec l σ : εc jrkdiv (q )

(1)

where ρ is the mass per unit volume and e the internal energy per unit mass. Here, variables with a dot above them denote time derivatives and, for second-order tensors a and b, the notation a : b represents a : b l

Σ ij a ij b ij l a ij b ij .
A similar treatment can be applied to the second principle that relates to the existence of an absolute temperature, T, and of a thermodynamic function, the entropy, obeying a constitutive inequality ensuring the positivity of the dissipated energy. This leads to the following local inequality:

ρsc jdiv E F q T G H k r T 0 ( 2
)
where s is the entropy per unit mass. Elimination of the thermal source, r, between Eqns. (1) and ( 2) and use of the free energy per unit mass, f l ekTs, lead to the so-called Clausius-Duhem inequality

Φ l σ : εc kρ ( f c jsT c )k q T :] T 0 (3)
which expresses the non-negativity of the dissipation, Φ, per unit volume. Note that the total dissipation, Φ, includes both an 'intrinsic' part, which is directly connected to the mechanical behavior, and a thermal part (k(q \T):] T ), associated with conduction.

Thermoelasticity

Thermoelastic behavior is defined by stating that the thermodynamic functions e and s as well as the stress tensor, σ, are one-to-one functions of the independent state variables ε and T. In this case, the intrinsic dissipation does not depend on ] T in Eqn. (3). This means that both the intrinsic and the thermal dissipation must be non-negative separately. From the definition

f c l cf cε ij εc ij j cf cT T c
it follows for any evolution (ε0 , TI ) from the equilibrium state (ε, T ):

E F σ ij kρ cf cε ij G H εc ij kρ E F sj cf cT G H T c 0 (4)
The independent variables ε0 and TI may have arbitrary values, including zero; this results in the following relationships:

σ l ρ cf cε (ε, T ), i.e., σ ij l ρ cf cε ij ; s lk cf cT (ε, T ) (5)
which are the constitutive equations of thermoelasticity (Salenc: on 1995). The free energy, f(ε, T ), from which the mechanical behavior can be derived completely for any thermodynamic conditions, can be considered as a thermo-dynamic potential, which is more general and powerful than the mechanical strain or stress potentials. The second expression of Eqn. ( 5) can be used, in association with the heat equation and boundary conditions, to derive the coupled thermomechanical evolution of a thermoelastic body. An alternative treatment, relying on the state description through the variables σ (instead of ε) and T, can be developed by using the stability of any deformed equilibrium state. It makes use of the potential f * defined by f * l (1\ρ) σ : εkf. The associated constitutive equations are then

ε ij l ρ cf * cσ ij (σ, T ), s l cf * cT (σ, T ) (6)

Elastic Moduli and Compliances

Symmetry

When the potentials f (ε, T ) or f *(σ, T ) are quadratic, the usual case of linear thermoelasticity is recovered: the constitutive equations (Eqns. ( 5) and ( 6)) yield a linear relationship between ε and σ through the classical elastic moduli C and compliances S, respectively. When this is not the case, the same equations can be used to connect linearly infinitesimal variations of stress and strain through ''tangent'' moduli L and compliances M. Their thermodynamic definition ensures the symmetry of these fourth-order tensors.

For example, for an isothermal evolution, L and M are simply given by the partial derivatives, at constant temperature, cσ\cε and cε\cσ. Specification of the indices and use of Eqns ( 5) and ( 6) lead to

L ijkl l cσ ij cε kl l ρ c#f cε ij cε kl , M ijkl l cε ij cσ kl l ρ c#f * cσ ij cσ kl (7)
Classical properties of second-order derivatives ensure that the isothermal elastic tangent moduli and compliances exhibit the diagonal symmetry

L ijkl l L klij , M ijkl l M klij (8)
This result is, of course, also valid for C and S in the case of linear elasticity.

Isothermal and Adiabatic Compliances

The elastic moduli and compliances depend on the thermodynamic regime. This can be illustrated by comparing isothermal (T constant) M iso and adiabatic (s constant) M ad compliances (Franc: ois et al. 1998).

They are defined by dε l M iso : dσ, dT l 0; dε l M ad : dσ, ds l 0 (9)

The isothermal compliances are given by Eqn. ( 7); the adiabatic ones result from Eqn. ( 6) and from the condition ds l 0. When the temperature dependence of ρ is neglected, the order of magnitude of the difference between M iso and M ad can be estimated by assimilating ε to the thermal strain, ε th l α∆T, so that cε\cT $ α, and cs\cT to C p \T, with α being the thermal expansion tensor and C p the specific heat at constant stress.

This results in

M iso ijkl kM ad ijkl $ T ρC p α ij α kl ( 10
)
(note that M ad also exhibits diagonal symmetry). For isotropic materials (α ij l αδ ij , where δ ij is the Kronecker symbol), and usual conditions, this proves that the relative difference between isotropic and adiabatic compliances hardly exceeds a few percent. However, this difference, which can only be understood from the thermodynamic treatment of elasticity, could have significant consequences in some extreme cases (dynamics, very high pressure or temperature variations, etc.).

Crystal vs. Rubber Elasticity

From Eqn. ( 5) and the relationship between f, e, and s, the constitutive equations of thermoelasticity can be put in the form

σ l ρ cf cε (ε, T ) l ρ ce cε kρT cs cε (11)
This expression clearly shows that thermoelasticity results from two sources: the strain dependence of the internal energy and that of the entropy. According to the materials under consideration, the first source or the second can be predominant: (i) some materials, especially crystals and polycrystals, cannot deform very much in the elastic range and there are very small changes of the entropy under deformation, whereas the strain sensitivity of the internal energy is quite high, owing to the nature of the interatomic interactions (so that σ $ ρce\cε);

(ii) other materials, such as elastomers (Mark and Lal 1982), can suffer very large reversible distortions of the macromolecular chains in the rubber elasticity regime, which strongly alters the entropy, almost without modification of the internal energy (so that σ $ kρT cs\cε).

For the first case, linear elasticity is generally an appropriate model; the elastic moduli are quite high and do not depend very much on temperature. On the contrary, rubber elasticity is characterized by much lower moduli (by a factor of 10$-10% on average) and a stronger temperature dependence. Unusual be- havior, which can be understood easily from a uniaxial version of Eqn. (11) (i.e., F $ kρT cs\cλ), is observed, such as shortening under load on heating or heat production by elongation. The entropic origin of rubber elasticity is easily proved by the temperature dependence of the force per unit area of the initial cross-section, F\S ! , needed to provoke a given elongation λ (new length over initial length).

Typical experimental data obey a two-slope scheme (Fig. 1). Below the glass transition temperature, T g , the force abruptly decreases with increasing temperature, as observed for crystalline materials; rubber elasticity is dominant for T T g , as proved by the increasing (almost linear) variation of the stress with temperature (the slope, given by kρcs\cλ, is positive since the entropy decreases with an increasing elongation) and by the very low value of the back-extrapolated stress at 0 K from the rubber-type response (kρce\cλ $ 0), which shows that the strain dependence of the internal is negligible.

Heterogeneous Materials (Linear Thermoelasticity)

This section is concerned with the prediction of the effective thermoelastic properties of heterogeneous materials, such as polycrystals or composite materials, from those of the constituents (grains or phases). Since the exact internal geometry of any representative volume element, Ω, of such materials is generally not known, these properties can only be estimated or bounded. For the sake of simplicity, restriction is made to linear behavior, isothermal conditions (with τ l TkT ! , a uniform change of temperature from the reference temperature, T ! , in the initial natural state), and homogeneous strain conditions, E. The derivation of bounds for the overall elastic moduli and thermal expansion coefficients especially benefits from a thermodynamic treatment of elasticity.

(i) Since potentials are additive quantities, the overall total free energy, F, directly derives from the local one by integration

ρ`F (E, τ) l f ρf (ε, τ) g l 1 Ω & Ω ρf (ε, τ)dΩ (12)
where ρis the overall mass per unit volume and f…g denotes a spatial average. Both at the local and the global levels, the free energy is quadratic and the constitutive equations are similar in form. From the local relationships

Σ l ρ cf cε l C : εkτφ, s l kρ cf cτ l φ : εjC ν τ, φ l C : α (13)
where C ν is the specific heat at constant strain, the following overall relationships can be derived (Stolz 1986):

l f σ g l ρ`c F cE l C eff : EkτΦ eff , ρ`s l f ρs g lkρ`c F cτ l Φ eff : EjC eff ν τ, Φ eff l C eff : α eff (14)
where the superscript ''eff'' refers to ''effective'' macroscopic quantities.

(ii) Owing to the linear character of the whole problem, the strain field, ε(x ), depends linearly on the loading parameters (E, τ), i.e.,

ε(x ) l A(x ) : Ejτη(x ) (15)
and a decomposition of the problem into two elementary ones, with E l 0 and τ l 0 respectively, can be carried out. This results in the following relationships:

C eff l f C : A g, Φ eff l f φ : A g, C eff ν l f C v gjf φ : η g (16)
These expressions show especially that the effective thermal expansion tensor, α eff , given by C eff : α eff l f (C : α) : A g, as well as the effective heat capacity at constant strain, C eff ν , do not reduce to the spatial averages of the corresponding quantities. This results basically from the fact that the thermal strain field is generally not compatible, which is responsible for the presence of a thermal residual (self-equilibrated) stress field and its associated stored elastic energy. The foregoing treatment can be developed in order to derive bounds for α eff and to perform a refined analysis of the thermal stresses, which cannot be obtained from a purely mechanical treatment of thermoelasticity.

Conclusions

The foregoing derivations have provided some illustrations of the benefits that can be gained from a thermodynamic treatment of elasticity. The main advantage lies in the fact that the temperature dependence of the elastic behavior is naturally taken into account so that the associated constitutive equations are valid for any thermodynamic conditions instead of being specified for particular ones (e.g., isothermal, adiabatic, etc.). Moreover, definite relationships can be derived between the elastic characteristics exhibited for different such conditions. This treatment yields a simple and meaningful interpretation of the entropic nature of rubber elasticity; in addition, it offers a comprehensive framework for the derivation of thermoelastic properties of heterogeneous materials, including the analysis of thermal stresses at the local level.

F

  Figure 1Temperature dependence of the force per unit area of the initial cross-section needed to maintain a given elongation for an elastomer sample.


Yves Michels

Étienne Baudrier

Loïc Mazo

RADIAL FUNCTION BASED AB-INITIO TOMOGRAPHIC RECONSTRUCTION FOR CRYO ELECTRON MICROSCOPY

Keywords: tomography, cryo EM, parametrization, unknown directions

A cryo Electron Microscopy dataset is composed of tomographic projections of an object (e.g. a macromolecule). The projection orientation information is unknown. The scope of this paper is the tomographic reconstruction of the observed object in the ab-initio case where the volume has to be estimated only from a raw projection dataset. A new approach based on a parametric model of the volume is presented. The description of the model and the search of the parameters are detailed. The accuracy and robustness of the proposed reconstruction method is shown on synthetic and real databases.

INTRODUCTION

In biology, the 3D structure of a protein plays an essential role in the determination of its biological functions. Cryo Electron Microscopy (Cryo EM) is one of the main imaging modalities used for macromolecule structure determination. It enables, after reconstruction to visualize the 3D electron density of large macro molecules (1,000 to 10,000 atoms). A Cryo EM dataset is composed of a large number (thousands to million) of images, each including a tomographic projection of the same molecule. Due to the acquisition principle, the projection direction (orientation) of each projection is unknown. Therefore, the construction of a 3D object from the Cryo EM dataset is a difficult problem with important biological consequences .

Because of the problem complexity, the computational cost of reconstruction methods is still an issue and can reach several hundreds of hours of CPU time for an ab-initio reconstruction. Recent softwares (SIMPLE [?], RELION [?] or Cryo-SPARC [?]) use massively parallelized implementations to reduce the execution time. Nevertheless, it is important to look forward the complexity reduction so as to treat bigger datasets (as heterogeneous datasets) or to reduce the computation time. We propose in this paper a method for abinitio reconstruction with a reduced number of operations that can be implemented in parallel.

The ab-initio reconstruction methods construct the density volume from a raw projection image dataset without prior volume information. Most of the ab-initio reconstructions involve two steps. In the first one, the orientations are estimated. Different approaches have been proposed in the literature: some are based on the common line property [?], some use non-linear dimension reduction to represent the projection images in a space where orientations can directly be estimated [?], and some use specific metrics as the moment based angular distance [?, ?]. In a second time, the volume is constructed. Several methods exist, using algebraic reconstruction [?], using iterative back propagation algorithms [?] or Fourier representation of the volume [?]. Other approaches propose to estimate jointly the orientations and the 3D volume. As in [?] where the construction is based on a cost optimization taking into account the common line property and the similarity between the projections of the estimated volume and the dataset. Facing the high computational cost of 3D computerized tomography, a reconstruction method based on a parametric model is presented in this study. The reduced number of parameters allows fast reconstruction resulting from an optimization.

The paper is structured as follows. In Section 2, the tomographic problem is defined and the reconstruction method is described. In Section 3, the implementation of our method including our optimization strategy is detailed and the evaluation of the proposed algorithm is given in Section 4. Finally, the conclusion of the study and some future works are given in the last section.

RECONSTRUCTION ALGORITHM

Let ρ : R 3 → R be a continuous function with compact support standing for the volume, and let R ∈ SO(3), the orientation, be a rotation matrix. Let P ρ (R) be the projection of ρ for the orientation R. The projection is given by the Radon transform with the formula:

P ρ (R)(y, z) = +∞ -∞ ρ R -1 (x, y, z) dx . (1)
We note π the projection given by the sampling, possibly noisy, of P ρ (R), and m the width of the projection image (and also the sampling resolution). Let Π be a finite set of sampled projections:

Π = {π k , k ∈ [[1, n p ]]}
where n p is the number of projections.

We propose to approximate the volume ρ by a sum of radial Gaussians ng l=1 ρ l where n g is the number of Gaussians and ρ l a radial Gaussian parametrized by its position, X l = (x l , y l , z l ), its width, δ l and its intensity ω l :

ρ l (x, y, z) = ω l ( √ 2πδ l ) 3 n g exp - R 2 l 2δ 2 l , (2)
where

R 2 l = (x -x l ) 2 + (y -y l ) 2 + (z -z l ) 2 .
We note S(n g ) the set of the n g × 5 parameters of the Gaussian sum and V (S(n g )) the Gaussian sum (abbreviated in V if no ambiguity). Note that the linearity of the Radon transform implies that the projection of the volume V is given by the sum of the projections of the Gaussians (ρ l ) l . For a given orientation, R, the projection of the volume V is given analytically by:

P V (R)(s, t) = ng l=1 ω l 2πδ 2 l n g exp - (s -s l ) 2 + (t -t l ) 2 2δ 2 l (3) where (r l , s l , t l ) t = R × (x l , y l , z l ) t .
Note that we choose the Gaussian functions nevertheless other integrable spherical functions can also be used.

The tomographic reconstruction is performed by jointly estimating the volume parameters, Ŝ, and the orientation parameters, { Rk , k ∈ [[1, n p ]]} that minimize cost function given by:

E = np k=1 m i=1 m j=1 (π k,i,j -πk,i,j ) 2 , (4)
where π k,i,j is the (i, j) pixel of the k th projection and πk the sampled projection of V ( Ŝ(n g )) with the estimated orientation Rk . Intuitively, for a Gaussian number n g small against the voxel number, this approach is faster than a voxel based volume as the optimization is done on a smaller set of parameters, e.g. for a 256 3 volume and n g = 50, there are 250 parameters for the proposed method in comparison to 256 [START_REF] Punjani | cryoSPARC: algorithms for rapid unsupervised cryo-EM structure determination[END_REF] for a voxel based one. In addition, each estimated projection is computed in O(m 2 ) instead of O(m 3 ) when line integrals are computed.

The optimization implementation is done with a stochastic gradient descent; it is easily parallelizable and is adapted to a multi-scale approach. The detail of the optimization is given in the next Section.

IMPLEMENTATION

The parameter estimation is done alternatively on the orientations and the volume, as described in Algorithm 1. Both optimizations are described in the next subsections. At each iteration, the gradient of E is computed using the expressions given in Equations (5, 6, 7).

∂E ∂ω l = 2 np k=1 m i=1 m j=1 (π k,i,j -π k,i,j ) 1 ω l πl,k,i,j , (5)
∂E ∂δ l = 2 np k=1 m i=1 m j=1 (π k,i,j -πk,i,j )× 2 δ l - (s l,k -s) 2 + (t l,k -t) 2 δ 3 l πl,k,i,j , (6)
∂E ∂X l = 2 m i=1 m j=1 π k,i,j -πk,i,j ∂ πk,i,j ∂X l , (7)
where

∂ πk,i,j ∂X l = Rk ×   0 s l,k -s t l,k -t   1 δ 2 l πk,i,j . (8)
The image πl,k is the discrete projection of V ( Ŝ(n g )) with the estimated orientation Rk . The position, width and intensity parameters are updated independently, due to their heterogeneity. To be independent from the normalization of the projection dataset, the intensities are updated with a constant step and the positions and the widths are updated using a global constant norm:

ng l=1 ||X (n) l -X (n+1) l || 2 2 = λ 1 , ng l=1 (δ (n) l -δ (n+1) l ) 2 = λ 2 ,
where • (n) correspond to the n th gradient iteration. The fixed λ 1 and λ 2 parameters depend on the desired precision of the parameter estimation. Note that our updating rules imply that passed a finite number of iterations, the resulting cost oscillate. Therefore, the volume optimization is terminated when the cost function stops decreasing.

Besides, to reduce computational cost and increase robustness to local minimums, each iteration is performed on a random subset of projections.

Projection orientation estimation

For a given volume parameter set, the search of the orientations is done independently for each projection. For a given projection, π k , the cost in function of the orientation is not a convex function. For the tested objects presented in 4.1, the cost function may have several dozens of local minimums as shown in Figure 1. Orientations are pre-estimated by an exhaustive search on a fixed grid containing 162 × 24 orientations regularly sampled on SO(3) then refined by convex optimization.

Due to the stochastic gradient descent, the whole optimization converge in a local minimum. Nevertheless only the resulting volume, V ( Ŝ(n g )), is important, not the parameter set, Ŝ(n g ). In addition, empirical evaluations show that the volume converges to an acceptable optimum that allows the orientations to converge to the global optimum.

EVALUATION

The reconstruction algorithm is performed on synthetic and real datasets to evaluate its ability to reconstruct finely the volume and to estimate the projection orientations.

Description of our datasets

Synthetic datasets are generated from Protein Data Bank files (PDB) where the centroids of the atoms which constitute the protein are given. The volume is generated with a sum of bounded Gaussians:

ρ a (x) = exp 1 -1/(1 - ||x|| 2 d 2 a ) ,
where d a is the average size of the amino acids. To reduce the complexity of the volume, each Gaussian is centered on amino acid centroids. An example of volume is shown in Figure 2. Our synthetic dataset is composed of 7 datasets of 5,000 projections with resolution 128 × 128. The datasets are generated from different proteins. Noisy dataset are noised by additive white Gaussian noise with the desired signal-to-noise ratio (SNR). The real dataset is composed of 1,000 projections from the TAF7 (TAFII55 complex), with resolution 160 × 160, averaged and corrected according to the contrast transfer function. Examples of projection images from different datasets are shown in Figure 3.

Results

The reconstructions are evaluated on the volume correlations:

ρ|V ||ρ|| ||V || ,
where •|• is the canonical scalar product. The estimated orientations are evaluated with the average angular difference between estimated and true orientations.

Two series of tests were performed, one to evaluate the ability to reconstruct finely the volume in the noiseless case, and one to evaluate the robustness of the reconstruction algorithm to the noise. The number of Gaussians used for the reconstruction are respectively 115 Gaussians and 35 Gaussians for the dataset with SNR +∞ and SNR 1/64. The results are given in Table The volume correlations in the noise free case are close to one. Therefore, our proposed method is able to reconstruct finely the volumes. In addition, the estimated orientation error is smaller than 0.02 rad (1.15 deg), except the first volume where 13 orientations error are higher than 1 rad due to the quasi symmetry of the volume. The results in Table 1 shows the robustness of our method to the noise, even with SNR 1/64 (-18 dB). Indeed, the constructed volumes have correlations higher than 0.9 except for one reconstruction. For totally asymmetric volumes, the orientation error is smaller than 0.05 rad. Our reconstructed volumes have been compared to the ones obtained by Cryo-SPARC for the 7 noisy datasets (SNR 1/64). The average correlation and angular difference are respectively 0.7881 and 0.7479 rad, in comparison to 0.9255 and 0.1355 rad for our reconstruction method. The estimated volume and orientations are more accurate for our ab-initio reconstruction that allow faster post processing refinement.

The reconstruction from real dataset is geometrically similar to a reference ab-initio volume obtained with a tilt acquisition as shown in Figure 4. Note that the reference volume is not finely reconstructed therefore the correlation coefficient is 0.858.

A test on 128 3 volumes where the number of Gaussians g is set to 50 and n p = 1.7 × 128, gives reconstruction times for our method more than ten times faster than the integral based reconstruction method of [?]. The existing demonstration algorithm will be developed in C++ and parallelized for GPU implementation. The optimal translation search will be added to the orientation estimation to improve the robustness of the reconstruction on real Cryo Electron Microscopy datasets. This method will be extended to heterogeneous datasets, arising from deformable volumes, where the deformations can be discrete or continuous.

Algorithm 1 : 3 . 1 .

 131 Reconstruction algorithmData: Π Result: Ŝ and { Rk , k ∈ [[1, n p ]]} initialization: random set Ŝ; while not converged do for k ∈ [[1, n p ]] do Rk = Orientations Estimation( Ŝ); end Ŝ = Volume Estimation( Ŝ, { Rk , k ∈ [[1, n p ]]});end Volume parameter estimation For a given set of orientations, { Rk , k ∈ [[1, n p ]]}, the volume parameters Ŝ are computed performing a gradient descent. The optimization step is initialized with the last estimated parameters.

Fig. 1 :

 1 Fig.1: Cost function for a projection from our real dataset in function of the angles θ and ϕ (on the sphere). For each orientation on the sphere, the in plane rotation is chosen to minimize the cost function.

Fig. 2 :

 2 Fig. 2: Synthetic dataset generated from a PDB.

Fig. 3 :

 3 Fig. 3: Example of projection images. (a) synthetic projection, (b) noisy realization of (a), (c) averaged projection from TAF7.

Fig. 4 :

 4 Fig. 4: Reconstructed volumes.

Table 1 :

 1 Reconstruction evaluation for synthetic datasets.

		Volume	angular	Volume	angular
	n	Corr.	difference (rad)	Corr.	difference (rad)
			SNR +∞		SNR 1/64
	1	0.9966	0.0368 ±0.1326	0.9554	0.2700 ±0.0458
	2	0.9662	0.0106 ±0.0065	0.9014	0.0275 ±0.0168
	3	0.9762	0.0137 ±0.0076	0.9587	0.0200 ±0.0115
	4	0.9663	0.0139 ±0.0087	0.8182	0.5249 ±0.9615
	5	0.9820	0.0147 ±0.0079	0.9329	0.0500 ±0.0345
	6	0.9903	0.0122 ±0.0064	0.9413	0.0224 ±0.0123
	7	0.9966	0.0190 ±0.0104	0.9703	0.0334 ±0.0241


H Maitournam

FINITE ELEMENTS APPLICATIONS NUMERICAL TOOLS AND SPECIFIC FATIGUE PROBLEMS

A systematic methodology for designing structures against high-cycle fatigue is developed. It relies upon: (i) computational methods for the calculation of the limit response of structures subjected to cyclic loading; (ii) numerical implementation of efficient high-cycle fatigue criteria. The first step in the prediction of high-cycle fatigue damage is the determination of stress cycles. In this paper we present three finite element procedures for the calculation of elastic-plastic structures subjected to cyclic loading (repeated moving contacts, small oscillatory contacts, etc.), namely the direct cyclic method, the stationary method and the simplified analysis of inelastic structures. These methods lead to easy determination of the possible stabilised response. Therefore, they avoid the lengthy repeated calculations performed with the classical finite element method and an incremental treatment of the loading history. The second step is the use of high-cycle fatigue criteria. We review some numerical tools for their implementations. Finally, some applications are presented. The first one is a problem of a cylinder subjected to an elliptical rotating pressure. The direct stationary method is used to determine the stabilised stress cycle and a macro-meso high-cycle fatigue criterion is used to detect crack initiation. The second application is the numerical simulation of a fretting fatigue map in relation with material fatigue properties.

INTRODUCTION

The systematic methodology for designing structures against high cycle fatigue we use, is summarized on figure 1. It consists in three uncoupled steps. The first one is the calculation of the mechanical stabilized state of the structures subjected to repeated loading ; the inputs needed are therefore the geometries of the structures, the materials they are made of, and the fatigue loading. The second step is the application of a crack initiation criterion with, as inputs, the intrinsic fatigue properties of the materials. And if necessary, a third step consisting in the use of propagation laws is performed. The mechanical stabilized state of the structures designed against high cycle fatigue is necessarily elastic. But in most of repeated contact problems, the first con-tacts generally generate plasticity setting up residual stresses, even in case of an elastic stabilized response (elastic shakedown). The quantification of those residual stresses is important in predicting the fatigue of such structures. It is well known that under a cyclic loading a rate independant elastic-plastic structure can have one of this four limit response : purely elastic, elastic shakedown or plastic shakedown, ratchetting. The possible mechanical stabilized state could be reached in few cycles or in an infinite number of cycles. It is therefore essential to have efficient numerical methods allowing a straightforward determination of this state. Classical finite element method with an incremental treatment of the loading history leeds to lenghty repeated calculations. The first part of this study is devoted to the presentation of three numerical proce dures for the direct calculation of the limit response of a structure subjected to cyclic loading : (i) the cyclic method (Dang Van et al.) for structures subjected to general cyclic loadings; (ii) the stationary method (Dang Van et al.) for structures subjected to moving loads; (iii) the simplified analysis of inelastic structures (Zarka et al.);

The second part of this paper is devoted to numerical analysis of fatigue. We describe how to numerically implement some fatigue criteria.

Finally, some applications are presented. The first one is a problem of rolling contact fatigue. It is represented by a cylinder subjected an elliptical rotating pressure. The direct stationary method is used to determinate the stabilized stress cycle. This limit state is necessarily elastic (or elastic shakedown) be avoid low cycle fatigue and incremental collapse. A macro-meso high cycle fatigue criterion is used to detect crack initiation. The second application is the numerical simulation of fretting fatigue tests on a particular experimental set-up considered as a structure. Fretting is the surface damage induced by small amplitude oscillatory displacements between components in contact. The direct cyclic method is used to evaluate the limit response of the structure. Then, we use macro-meso multiaxial fatigue criterion identified by classical fatigue tests to predict crack initiation.

DIRECT CYCLIC METHOD

As far as we are concerned by high cycle fatigue the limit response of the structure subjected to repeated loads is with elastic shakedown. The evolution towards this state is not important. The cyclic direct method allows the direct determination of the asymptotic response (i. e. stabilized mechanical state) of a structure subjected to a general cyclic loading without an incremental treatment of the whole loading history. The structure considered is made of von Mises elastoplastic material. It is subjected to a cyclic external loading composed of prescribed forces fd and displacements ud(t), which are time-periodic with a period T.

2.1

Principle of the method

We seek directly for the stabilized mechanical response of the structure by using the two following procedures : (i) large time incremental method developed Ladeveze et al. ([9, 10]), departing form step by step scheme and using a single "large" increment of time ; (ii) research of the solution in the space of periodic responses ( [START_REF] Akel | Determination of the limit response in cyclic plastic ity[END_REF]).

In the incremental method, the time interval [0, T] is discretized in n subintervals 0 < t1 < . . . < t n = T. For each loading increment ( !:J.fd, !:J.ud ) , the incremental quantities !:J.u, !:J.e, !:J.eP are evaluated such as to satisfy static, plastic and kinematic admissibility conditions. The loading path is followed step by step and cycle by cycle until the limit response is reached. The approach we propose here, is also iterative. But, it consists,

• (i) first, in calculating , for the whole cycle, the statically and kinematically ad missible responses ; practically, we calculate this solution at the discrete instants t; of the cycle; for each iteration j, at instant t;, we determine the stress tensor u(, which satisfies equilibrium for the loading (fd(t,)) and the strain tensor E{, which respects with the prescribed displacements ud(t,) corresponding to instant t,; these two conditions are equivalent to the resolution of equilibrium equations; for instance, when a constant stiffness is used, they are given by : in which F,, F,P'-1 and K are the vector of nodal external forces corresponding to f d (t,), the vector of plastic forces due to plastic strain Ef 1 -1 (assumed known) and the stiffness matrix of the linear elastic structure, respectively;

• in a second stage, in determining for the whole cycle, the plastically admissible response (the stress tensor u;: 1 , the plastic strain €�� 1 , the internal parameters at (•+I) ) ; at each time, the plastic strain €�� 1 is obtained by projecting the strain increment €; +1 -E{ on the elastic convex defined by the initial state ( u(, E/ 1 , at , ) p at time t, which is known; u; + l is then calculated from the constitutive law written in incremental form :

• The stabilized state is obtained when the plastic strain E1P' and the internal pa rameters (a� 1 ) are periodic i.e. their values at the beginning of the cycle are equal to those at the end of the cycle( €�1 and a� n ; otherwise, the same steps are performed in iteration j + 1 by initializing ( �:: f1+1, o:17 1 ) to ( E�1, o:1 "!; 1 ) and the pro cess is repeated until the obtention of elastic or plastic shakedown. Ratchetting is revealed by the non convergence of the procedure. This method is typically used for fretting problems in which we are concerned by small oscillatory contacts.

3

STATIONARY METHODS

3.1

Recall of the principle of the 'stationary methods'

The method is applicable in the case of structures subjected to repeated moving loads. The full details of the 'stationary methods' used in this paper and their finite elements implementation can be found in [START_REF] Van | Steady-state flow in classical elasto plasticity : application to repeated rolling and sliding contact[END_REF]. Here, we just recall the principle of the analysis. Let us consider a prismatic structure made of an elastic-plastic material (let's say a rail) subjected to loads (let's say normal and tangential pressures) moving with a velocity V= Vex (Fig. 2). We assume that there is a steady state in a reference frame moving with the loads ; for any material quantity A, it leads to the following relation :

A. = -VA, x (1)
) )

H u u ( = 1 if f = 0 (2) (3) (4)
(5)

( 6) (7)
(8)

(9) (10)
( . ) meaning positive part, and

( = 0 if f < 0 (11) (12)
A is the plastic multiplier, ( f = 0) defines the elastic domain, ke is the shear yield stress, a the coefficient of thermal expansion, (} the temperature elevation, u0 the ini tial stress, ak are internal parameters and Ak their associated forces, L tensor of elastic coefficients, Z and l are assumed known.

In these equations time derivatives have been replaced by space derivatives, ac cording to (Eq. ( 1)). Their integrations are performed in space, along the direction of the motion of the loads; so that, numerically (using the finite element method), neither the load, nor the structure are translated. The plastic behaviour can be very gen eral, for instance, in [START_REF] Van | Steady-state flow in classical elasto plasticity : application to repeated rolling and sliding contact[END_REF][START_REF] Van | Thermomechanical state near rolling contact area, Dissipative Processes in tribology[END_REF][START_REF] Van | Elastoplastic Calculations of the Mechanical State in Alternative Moving Contacts: Application to Fretting Fatigue, FRE T TING FA TIGUE[END_REF] the von Mises isotropic and/or linear kinematic hardening material, the non-linear kinematic hardening law proposed by Bower [START_REF] Bower | Cyclic hardening properties of hard-drawn copper and rail steel[END_REF] has been used.

The second idea is to use the stationarity of the limit state (in the absence of ratchetting the stress field and the internal variables are are periodic) in order to find its value directly. This condition is written as :

1 f.Pd t = 0 cycle and (13)
The two following numerical procedures for the determination of the stabilized state are derived [START_REF] Van | Steady-state flow in classical elasto plasticity : application to repeated rolling and sliding contact[END_REF].

(i) The pass-by-pass stationary method (PPSM) for the calculation of a single pass; so by computing the successive passes by this method, one can find all the features of the "stabilized" state (numbers of cycles before reaching it, residual stresses, plastic deformations, etc). This "stabilized" state can, of course be ratchetting; in such a case, the ratchetting rate is immediately deduced as the increment of plastic strains caused by a pass.

(ii) The direct stationary method (DSM) for the direct determination of the "stabilized" state if it is a shakedown (elastic or plastic); The ratchetting is indicated by a non-convergence of the algorithm.

Thus, these methods allow the description of strains and stresses due to cyclicly and alternatively moving contacts, the quick determination of the nature of the sta bilized state (elastic shakedown, plastic shakedown or ratchetting). Such analyses of repeated contacts using shakedown principles were carried out with semi-analytical methods, by Johnson [START_REF] Johnson | Contact mechanics and wear of metals[END_REF][START_REF] Johnson | The application of shakedown principles in rolling and sliding contact[END_REF][START_REF] Johnson | Contact mechanics[END_REF].

SIMPLIFIED ANALYSIS OF STRUCTURES

The Simplified Analysis of Inelastic Structures proposed by Zarka et al. [START_REF] Zarka | A new Approach to Inelastic Analysis of Structures[END_REF] leads to the direct calculation of the limit cycle (elastic or plastic shakedown) for an elastic-plastic material subjected to cyclic loading. It is based on the introduction of transformed variables linearly related to the internal parameters, and on simple rules for the deter mination of theses variables in the limit cycle. Elastic analysis are then only needed.

NUMERICAL IMPLEMENTATION OF FATIGUE CRITERIA

The first step is the numerical determination of the stabilized response of the structure (elastic state or elastic shakedown). The application of the criteria is essentially local at each point of the stress evaluations. For a point M, the stress cycle for a period T is discretized in time, t1, t2, ... tn = T, n is the number of instants. The stress tensor u is put in the 6-component vector form with, 0" 1 = 0" 11, 0" 2 = 0" 22, 0"3 = 0"33 , 0"4 = 0"1 2, 0"5 = cr23, cr5 = cr13. One defines : [START_REF] Bower | Cyclic hardening properties of hard-drawn copper and rail steel[END_REF] The stress deviator is S . It is defined in a 5-dimension space. (If for instance, we choose the stress deviator vector as 51 = (311 -522)/2, 52 = (311 + 822)vf3j2, 83 = 833,84 = 812, 0"5 = cr23 then J 2 = 11 S 11). .jJ:;: is the amplitude of the square root of the second invariant of the stress deviator.

It can be defined :

• as half of the longest chord, so it is simply determined by ;

( 1 5 ) and in a detailled form In this case (see the contribution of I. Papadopoulos) the determinatioi_t of the mean value ( p;;. ) can be ambiguous.

• as the elastic fluctuation : [START_REF] Johnson | Contact mechanics[END_REF] This is a classical min-max problem, which can be solved by using the solution proposed by Papadopoulos, or using other numerical procedures given for instance in [START_REF] Zhou | Non monotonic line search for minimax problems[END_REF].

One defines S* as the stress deviator for which min8 r max; 11 S ( t;) -sT 1 1 is obtained. So, p;;. = � m f x 1 1 S ( t ;)-S* 11 The first example is that of a rotation of a solid over a cylinder (Fig. 3). The contact surface is assumed to be elliptical with semi .. The application of fatigue criteria over the whole structure leads to a critical point located at the surface. On figure 5 we show in the ( T, p) diagram, loading paths for different depths below the contact.

material line

Fretting fa tigue proble m

This section is devoted to the numerical simulation of a fretting test. The experimental set-up is presented in [START_REF] Petiot | An analysis of fretting-fatigue failure combined with numerical calculations to pre dict crack nucleation[END_REF] and in [START_REF] Maouche | On a new method of eval uation of the inelastic state due to moving contacts[END_REF]. The prescribed oscillatory between the pads are linked to the prescribed oscillatory fatigue stress S(t) in the specimen. For a maximum stress Smax=500 MPa, the am plitude of displacement is 0.55 f.1m. The flexible beams are equiped with strain gauge in order to measure the clamping force P between pads and specimen and the friction force related to the displacements accommodation. The variations of the tangential force T( t) are recorded for each cycle fatigue and plotted as function of fretting fatigue stress S ( t) (fretting fatigue loops). By varying the operating parameters ( P, Smax ) , three regimes are established:

• Stick regime : Fretting fatigue loops keep a non evolutive closed shape. Loops are quite linear during the test. The macroscopic displacement between the contacting surfaces is mainly accommodated by elastic deformation in the near surface of the two components. No damage (wear or crack nucleation) appears during the 107 cycles of the test.

• Mixed stick-slip regime : Loops present an elliptical closed shape. There is partial slip and fatigue crack nucleation observed at the edges of the contact.

• Gross slip regime : Loops present a trapezoidal shape. Full slip occurs between the two contacting surfaces. In this regime, particle detachment is observed. The differents regimes are obtained for differents varying parameters ( P, Smax ) summarized in the map shown in Fig. 8.

7.2

Numerical prediction of wear and crack nucleation

The prediction of damage mechanisms carried in the map Fig. 8 requires first the calculation of the stresses history in the stabilized state. Secondly, Dang Van multiaxial fatigue criterion is applied if the stabilized state is elastic shakedown.

The direct cyclic method described in Sec. 2 is used here to simulate the set-up and to calculate the stresses history in the stabilized state. The specimen is modeled as a half space subjected to a constant normal force and a varying tangential force T(t) ami a fatigue stress S(t) varying linearly with T(t). The material is elastoplastic with a kinematical hardening (hardening modulus C=30 GPa) and with properties given in the previous subsection. Three cases are performed corresponding to the three regimes established in the map ( P and Smax are the prescribed parameters and T max is measured in the test).

• Stick regime P=140 N, Smax=350 MPa and Tmax=53 N (no damage) • Gross slip regime P=80 N, Smax=500 MPa and Tmax=64 N (wear)

In the case of loading in gross slip regime, the nonlinear stress-plastic strain cycles for an element located at the surface under the contact zone (x = -0.4a) are shown in Fig. 10 illustrating plastic shakedown. Wear is then associated with plastic sakedown and so to low cycle fatigue properties confirming the works of Ludema [START_REF] Kim | A correlation between Low cycle fatigue and scuffing properties of 4340 steel[END_REF] and John san [START_REF] Johnson | Contact mechanics and wear of metals[END_REF]. Owing to this numerical method, we obtain, better still, the plastic strain amplitude in the stabilized state. Since the response of material is respectively purely elastic and elastic shakedown in stick regime and mixed stick-slip regime, high cycle fatigue is concerned. The stress

Fig. 1 :

 1 Fig. 1: General methodology for designing structures against fatigue

Fig. 2 :

 2 Fig. 2: Structure under repeated moving load

3. 1 . 1

 11 Mechanical problemThe first idea consists in writing the equations governing the mechanical problem in the moving loads reference. Thus one obtains[START_REF] Van | Thermomechanical state near rolling contact area, Dissipative Processes in tribology[END_REF] : with -equations of motion -constitutive laws u = L: (e::: e-aBI) + Uo p -A 8 f -€ ,x -8 u f � 0, A 2: 0, Af = 0 f � 0, A 2: 0, Af = 0 1 8 f 8 f A = -( ( --8 :L:e::: ,x+aB,xtr(-

5 . 1 imple mentation of criteria based on stress invariants 5 . 1 . 1

 51511 Calculation of .jJ:;: and � :

  t;)-S 1 (ti))2 + (S 2 (ti)-S 2 (ti))2 + (S3(ti)-S3(ti))2+ ( 1 6 ) 2 (S4(t;)-S4(ti))2 + 2 (S5(t;)-S5(ti))2 + 2 (S6(t;)-S6(ti))2

5. 1 . 2

 12 Calculation of L.H ,a and L.H,m :'L,H -� ( �H -�H . ) ,a -2,max ,mzn1 �H = -( �H + �H . )

  -axes a et b. a is the semi-axe along the longitudinal direction b the semi-axe along the circonferential one. The ratio a/b is fixed to 5. The material constants are : E (Young modulus) = 210 GPa, v (Poisson coefficient) = 0.3, k (shear yield strenght) = 200 MP a. The pressure distribution is also assumed to be elliptical with maximum pressure P0 = 800 MPa. The limit state reached is an elastic shakedown. The residual stresses obtained are shown on figure 4 .

Fig. 3 :- 1 .Fig. 4 :

 314 Fig. 3: Cylinder under rotating contact 100 .------r-----.------,-----�------.------,

Fig. 5 :

 5 Fig. 5: Loading path for different depths at the contact center (po/ k = 4., f1 = 0 . )

  Two cylindrical fretting pads (diameter 10 mm) are clamped against the two surfaces of a flat uniaxial fatigue specimen tested under constant amplitude loading at a frequency of 20 Hz. The pads are made of 100C6 steel and the fatigue specimen is made of 3Cr-MoV steel. The mechanical properties are given below : Yield strength Tensile strength Young modulus Hardness (

  Fig. 6: loading path during fretting • Mixed stick-slip regime P=100 N, Smax=600 MPa and Tmax=80 N (crack nucle ation)

Fig. 7 :

 7 Fig. 7: loading path during fretting

Fig. 8 :Fig. 9 :

 89 Fig. 8: Fretting map (Petiot et al, [12])

5.2

Implementation of criteria based on critical planes One defines e• as the shear stress for which min er max; 11 e (t;) -eT 11 is obtained. So,

'

5.3

Macro-mesa approach i mplementation Dang Van et Papapoulos ( [START_REF] Ballard | High cycle fatigue and a finite element analysis[END_REF]) criteria implementation needs the calculations of the microscopic stresses u( t;) , which are the sum of the macroscopic stress �(M, t;) and the microscopic residual stresses p :

�(M, t;) = A(M)u(M, t;) + p(m, t;)

A is the localisation tensor.

The assumption of elastic shakedown for the microscopic stresses leads for time independant microscopic residual stresses ( S*). They are are evaluated by the following min-max problem :

� min max 11 5EL(M, i;)-ST 11= k* Papadopoulos gives a method for solving this problem. Another method, based on Faisible Sequential Quadratic Programming, can be found in [START_REF] Zhou | Non monotonic line search for minimax problems[END_REF].