Fabien Salaun

Pierre-Olivier Chapuis

Kumar Sourabh

Sebastian Saha

Volz
email: volz@em2c.ecp.fr

Polymer nanoparticles to decrease thermal conductivity of phase change materials

Keywords: polymer nanoparticles, phase change microparticles, microencapsulation

Microparticles including paraffin are currently used for textiles coating in order to deaden thermal shocks. We will show that polymer nanoparticles embedded in those microsized capsules allow for decreasing the thermal conductivity of the coating and enhance the protection in the stationary regime. A reasonable volume fraction of polymer nanoparticles reduces the conductivity more than predicted by Maxwell mixing rules. Besides, measurements prove that the polymer nanoparticles do not affect the latent heat and even improve the phase change behaviour as well as the mechanical properties.

Introduction

In the two past decades, microencapsulated Phase Change Materials (PCM) have drawn an increasing interest to provide enhanced thermal functionalities in a wide variety of applications [START_REF] Colvin | Microencapsulated Phase Change Materials for Storage of Heat[END_REF][START_REF] Kaska | Improvement of the Performance of Solar Energy of Waste Heat Utilization Systems Using Phase Change slurry[END_REF][START_REF] Mulligan | [END_REF][4]. When the encapsulated PCM is heated to the melting point, it absorbs heat as it goes from a solid state to a liquid state. This phase change slows down the temperature increase. It can be applied to clothes technology, building insulation, energy storage as well as to coolant liquids. On a more general basis, it can be used to design a broad variety of thermal transient regimes.

The PCMs used for ambient temperature related applications are carbohydrates with different chain lengths [5][START_REF] Pause | Development of Self-conditioning Sport Shoe with PCM[END_REF][START_REF] Leitch | [END_REF][START_REF] Bryant | Techtextil-Symposium, Fibers with enhanced, reversible thermal energy storage properties[END_REF] or paraffins [START_REF] Mccullough | Industrial fabric protective review[END_REF][START_REF] Nakahira | Heat-accumulating microcapsule dispersion[END_REF]. The PCMs are encapsulated in small spheres in order to be contained in a liquid state. The microcapsules possess approximate diameter of 1-10 µm and are resistant to abrasion, pressure, heat and chemicals according to the shell chemical compounds [START_REF] Jahns | Application of microcapsules as latent heat accumulators[END_REF][START_REF] Dietrich | [END_REF] Micro-sized capsules are required to provide a large contact area with the environment and ensure an optimal efficiency of the phase change.

Among a multitude of possible wall materials for microcapsules, amino resins, more especially melamine-formaldehyde, play a main role, more precisely in the patent literature [24]. Amino resins represent an interesting economical alternative, since these polymer raw materials have been produced on a large scale for many years, and since they have already been practised in many processes like phase separation and interfacial reaction (e.g. dicarboxylic acid dichlorides and di or triamines).

Preparation of microcapsules

The choice of the process of synthesis is generally dictated by the chemical nature of the polymer forming the membrane and by the desired mean diameter. The thermal properties of a microcapsule are influenced by three factors: its mean diameter, its expansion during the phase change process and its shell's chemical nature. In his study, Bryant noticed that the smaller the mean diameter was, the more the weight loss was significant at high temperature [16]. The microcapsule wall is greatly affected by the phase change process, since between the melting and the crystallization, the volume varies with an approximate magnitude of 10 %. This can induce the diffusion of the core material at high temperature. Most of synthesized microcapsules are intended for application fields at a temperature lower than 150 °C. Therefore, microcapsules containing phase change materials have yet been synthesized with urea-formaldehyde [START_REF] Yoshizawa | Membrane formation mechanism of crosslinked polyurea microcapsules by phase separation method[END_REF], urea-formaldehyde-resorcinol [START_REF] Soto-Portas | Preparation of oily polyamide microcapsules via interfacial polycondensation[END_REF], polyamide [START_REF] Uddin | Effects of cyclic operation on the characteristics of a microencapsulated PCM storage material[END_REF], gelatin-formaldehyde [START_REF] Bryant | Melt Spun Fibers Containing Microencapsulated http://patft[END_REF], polyurethane [START_REF] Zhang | Structure and Thermal stability of microencapsulated phase-change materials[END_REF], melamineformaldehyde [START_REF] Cho | Microencapsulation of octadecane as a phasechange material by interfacial polymerization in an emulsion system[END_REF] and urea-melamine-formaldehyde [START_REF] Fan | Super-cooling prevention of microencapsulated phase change material[END_REF] shell. In this last study, the factors improving thermal stability were the addition of cyclohexane to the noctadecane and the copolymerization of urea with melamine and formaldehyde.

The paraffins N-hexadecane, n-eicosane, tetraethyl orthosilicate and sodium phosphate dodecahydrate used as PCMs formulation were purchased from Acros organics. Melamine-formaldehyde resin (Arkofix® NM) used as shell-forming was kindly supplied by Clariant (France). Poly(vinyl alcohol) (95 % hydrolysed, average M.W. 95,000) and methylene diisocyanate (Suprasec 2030, Huntsman ICI) were employed to entrap inorganic salt. Nonionic surfactants, Tween® 20 (Acros Organics, France), Span® 85 and Polyethylene glycol 400 dioleate (Aldrich, France) were used as emulsifiers. For pH control, triéthanolamine and citric acid were used (Aldrich, France).

Our technique of microencapsulation for the simple core shell structure consists in a three step process. The first step is the crucial phase of emulsion. All physical properties of the microcapsules will be established at that time (e.g. size distribution, morphology, shape, …) [START_REF] Mc Quillan | Microencapsulation process factors which influence the sphericity of 1 mm o.d. poly(α-methylstyrene) shells for icf, general atomics report GA-A22850[END_REF][START_REF] Maa | [END_REF]. Indeed, the dispersion of the organic phase in the continuous phase is prepared by a fragmentation process where wide areas of oil are broken up by shear or elongation. Therefore the droplet sizes are highly affected by the physicochemical properties of the two immiscible liquids (e.g. viscosities, densities, pH, interfacial tension, temperature, …) and by the agitation system (especially stir speed and shape of stirrer) [15]. The preparation of our microparticles was carried out according to the method of in situ polymerisation of melamineformol. The inner phase was emulsified in an aqueous solution containing 4 g of Tween® 20 in 100 g of water and 9.2 g amino resin.

In the second phase, the pH was reduced to 4 by adding acid citric solution, at a stirring rate of 8000 rpm, at room temperature, during 3 minutes with a homogenizer.

The PCM core particles are generated in this phase.

Finally, the reaction mixture was heated at 60°C, stirring was continued using a blade stirrer at low speed (400 rpm, RW20, IKA, Germany) for 4 hours until the end of the polycondensation. The microcapsules were recovered by filtration, washed with methanol and water, and dried at room temperature during a night.

Figure 1 reports a Scanning Electron Microscope (SEM) image of the resulting coreshell microcapsules named 'ARD'. The size distribution appears as quite homogenous and the particle diameter ranges from 2 to 4 micrometers. The composition of each inner phase used in the preparation is listed in Table 1.

To produce polymer nanoparticles in the basic core shell structure, the first step consists in generating polymer spheres in paraffin. In this aim, we mixed 4 g of sodium phosphate dodecahydrate and 2 g of water. The solution is added to 0.5 g of mixture of nonionic sufactants in 7 g of n-alkane. After stirring during 15 min, the droplet size was reduced by homogenizing the emulsion during 15 min at 9500 rpm with an Ultra-Turrax homogenizer (Ika, Germany). In the same way, another emulsion was prepared by homogenizing 8 g of 95 % hydrolysed poly(vinyl alcohol) solution (5 wt.-%) in 8 g of n-alkane. The polymer spheres were prepared by shearing under high speed the two emulsions with 3 g of MDI to crosslink the shell at 50°C for 30 min. The resulting polymer nanoparticles were observed under SEM as illustrated in Figures 2. The mean diameter is about 50 nanometers and the scattering in size seems neglectible. The porous material formed by those particles is labeled '18' in the following.

In the second and last step, the previous mixture is introduced in the place of the basic PCM in the microencapsulation process described above. The surface morphology of the final product was observed with an optical microscope (Zeiss Axisokop) in the phase difference mode as shown in Figure 3 (top). The inner structure of the microcapsule is revealed in Figures 3 (bottom). A rather compact distribution of polymer nanoparticles is shown in both images.

Results and discussion

Latent heat

A conventional Different Scanning Calorimetry (DSC) was performed to prove that the polymer nanoparticles do not affect the latent heat. The thermal behavior of the microcapsules was recorded to analyse the influence of nanoparticles on latent heat storage, using a TA instrument type DSC 2920 piloted on PC with TA Advantage control software. The measurements were performed by using increase rates for the temperature ranging from 0.5 °C/min to 20 °C/min. The values of latent heat were Submitted to Polymer 6 found equal within the measurement accuracy for both pure paraffin microparticles and polymer nanoparticles based microcapsules. The value of the latent heat is equal to 176 J/g. The DSC data reported in Figure 4 however indicate a difference in the phase change behaviour. The temperature at which the solid-liquid transition starts is the same but the temperature interval is 50% larger for the NPs based PCM compared to the conventional mirocapsules. The phase change will hence be effective on a larger temperature range. We presume that interaction between nanoparticle surface and paraffin induces a specific molecular structure in paraffin.

An extra energy is necessary to break this molecular arrangement. This modification in the phase change is of course beneficial for textile application because temperature will be maintained under a larger heat flux amplitude.

Thermal conductivity

The thermal conductivity is estimated by using a scanning thermal microscope. It consists in a conventional atomic force microscope mounted with a hot wire probe as illustrated in Figure 5. The probe is a wollaston wire made of a platinum core 5 microns in diameter and a silver coating 70 microns in diameter. The silver coating is etched to uncover the platinum wire over a length of 2L=200 microns. This tip was studied in several of our previous works [START_REF] Bryant | Techtextil-Symposium, Fibers with enhanced, reversible thermal energy storage properties[END_REF]. A modulated (AC) electrical current is used to Joule heat the wire and the second harmonic of the temperature is measured [START_REF] Mccullough | Industrial fabric protective review[END_REF].

The wire temperature is related to the heat flux flowing from the tip to the sample.

And this heat flux between the tip and the sample is driven by a global thermal conductance G eq including the contributions of the contact conductance and the one of the sample G s .

We perform a two-step approach to identify the ratio between two thermal conductivities λ S1 and λ S2 of two different samples. Measurements are performed in ambient and in vacuum when varying contact force F between sample and tip. Those measurements provide the force derivatives of the sample conductance and the one of G eq . We assume that the contact conductance and the contact hardness do not vary significantly when probing the two dissimilar samples. In those conditions, it can be shown that the ratio between both thermal conductivities arises as follows:

λ S1 λ S 2 = dG S dF 1 dG S dF 2 ⎛ ⎝ ⎜ ⎞ ⎠ ⎟ 1/ 2 dG eq dF 2 dG eq dF 1 ⎛ ⎝ ⎜ ⎜ ⎞ ⎠ ⎟ ⎟ 1/ 4
(1)

Note that Eq. ( 1) does not rely on the determination of a specific mechanical or thermodynamic property. including polymer. Eq. ( 1) provides the ratio of λ 18 / λ ARD =0.31. From the reference data λ ARD =0.26 W.m -1 .K -1 and the very low thermal conductivity of λ 18 =0.08 W.m - 1 .K -1 is found to be only three times the one of air.

Figure 7 presents the evolution of the thermal conductance G eq as a function of the electrical power in the probe. The hotter the probe is, the more resistive samples are.

The linear behavior suggests that the samples might melt and become more resistive.

The relative values of the samples thermal conductivities are reported in the insert of to 70% more insulating than paraffin. This confirms the impact of the very insulating polymer nanoparticles on the PCM thermal conductivity.

Samples E2 and H/E include 68% of paraffin and 32% of polymer. Consequently, Maxwell-Garnett mixing rules predict a thermal conductivity decrease of 23%, which is about threefold less the experimental data of 60%. The thermal resistances at interfaces between polymer nanoparticles and paraffin might explain such a discrepancy. A volume fraction of 10% also corresponds to an average inter-particle distance of 0.35 time the particle diameter, i.e. about 15nm. We therefore suggest that percolating networks of polymer nanoparticles might be responsible for the drastic reduction in heat conduction.

This property implies a considerable gain, especially in the field of textile because the temperature drop between the inner and the outer face is increased threefold.

Bulk modulus

The Atomic Force Microscope allows for probing the mechanical properties of the surface. Figure 8, provides the force applied on the tip versus the cantilever altitude.

Submitted to Polymer

The stiffer the material is, the larger the slope is. Sample 18 is the stiffest material because it does not include any paraffin. The E2 sample has a lower melting point than H/E. The corresponding slope is the same as the one of sample 18 because we presume that the polymer contribution predominates because of its structural network included in a solid matrix. In sample H/E, paraffin is in the liquid phase at ambient which might explain the slight decrease of the slope. Isolated polymer nanoparticles might indeed slide among the aggregates. The ARD sample has the same slope as the H/E one in the range of the small forces because they include the same paraffin. The data point corresponding to the largest force seems to indicate that the capsule shell was broken. This mechanism appears only in the case of sample ARD because the polymer structure does not enforce the capsule.

Conclusion

We have fabricated core shell PCMs including polymer nanoparticles and have characterized their structural configuration via Scanning Electron Microscopy. The diameter of polymer nanoparticles is in the vicinity of 50nms and the shell size is of a few microns. We have shown that the addition of polymer nanoparticles has a very beneficial impact on thermal properties. Firstly, we have proven by using usual DSC that the temperature interval of the phase change is augmented by 50%. This significantly improves the role of thermal barrier of the PCM. We suppose that this effect is related to the ordered structure of paraffin molecules around the polymer nanoparticles. This configuration might require a larger amount of energy to start for melting.

Secondly, the thermal conductivity of polymer based PCM is decreased by 60% in comparison with the pure paraffin sample. This reduction is threefold larger than the one predicted by conventional Maxwell-Garnett mixing rules. We have highlighted that the average separation distance between particles is about 15nm only and aggregation occurs as proven by the SEM snapshots. The network of polymer nanoparticles might break the percolation of paraffin and hence explains this very special behavior.

  Furthermore melamine-formaldehyde microcapsules prepared by in situ polymerisation have impermeable shell.Up to now, most studies have targeted microcapsules, but encapsulation of micronanospheres remains an innovative field[25]. In this work, we show how polymer nanoparticles (NP) embedded in PCM microcapsules significantly decreases the thermal conductivity to improve the thermal barrier effect. We also prove that the modified microcapsules have better mechanical properties and even an enhanced phase change behavior. Submitted to Polymer 3 The fabrication of the simple core-shell PCM and of the polymer NPs based microcapsules are presented in the second section. The thermal measurements are described in the third paragraph. The phase change behaviour was characterized by conventional Differential Scanning Calorimetry to measure the latent heat and the temperature range of the solid-liquid transition. Both materials exhibited the same latent heat but the phase change occurs on a wider temperature interval in the case of the NP based PCM. The thermal conductivity was measured with a Scanning Thermal Microscope because the conventional methods are difficult to apply. Contact technique such as hot guarded plates introduce thermal contact resistances which are difficult to estimate and the low optical absorption of porous materials makes the application of optical techniques intricate. The value of the thermal conductivity of the polymer NPs was found to be extremely low, i.e. only twice the one of air. The addition of a reasonable fraction of NPs in the PCM capsule results in a decrease by a factor of three of the effective thermal conductivity. The outcome is a significant improvement of the thermal barrier in the steady state. The mechanical properties were estimated on a qualitative basis by using an atomic force microscope. The damage threshold was reached in the basic microcapsule before the one of the modified materials. The results are described and commented in the last part.

Figure 6

 6 Figure 6 reports on the thermal conductance G eq against force for samples 18 and ARD. The slopes are not clearly different showing that the contact conductance is predominant over the sample conductances. But the difference in levels indicate that the paraffin in sample ARD is significantly more conductive than sample 18

Figure 7 .

 7 Figure 7. Microparticle based samples, that are referred as E2 H/E and 18, are 60%

CaptionsFigure 1 :

 1 Figure 1: Scanning Electron Microscope image of the ARD sample. The spherical shell structure has a size of 1 to 3 microns.

Figure 2 :

 2 Figure 2: Scanning Electron Microscope image of sample 18. The polymer nanoparticles have diameters of about 50nms. The high particle density allows for aggregation.

Figure 3 :

 3 Figure 3: Optical (top) and Scanning Electron Microscope (bottom) images of sample E2 proving the high density and the configuration of polymer nanoparticles inside the shell. Although shells include 68% of paraffin, the nanosized polymer particles generates a very dense and thermally insulating structure.

Figure 4 :

 4 Figure 4: DSC characterization of the E2 sample showing the extension in the temperature interval of the phase change. Several increase rates for the temperature are displayed. Each of them leads to the same latent heat.

Figure 5 :

 5 Figure 5: Sketch of the thermal probe. A wollaston wire is formed as a tip. A mirror is stuck on the top face to allow for the detection of the tip deflection. The electrical resistance of the wire is measured to deduce its temperature and its dissipated heat flux.

Figure

  Figure 6a,b: (top) Thermal conductance G eq related to the heat flux exchanged between the thermal probe and the sample versus contact force between tip and sample. (bottom) Thermal conductance G s of the sample versus contact force between tip and sample. The force is proportional to the tip deflection, which is measured by the photodiode presented in Figure 5. The force units are hence provided by the photodiode signal in NanoAmperes.

Figure 7 :Polymer 12 Figure 8 :

 7128 Figure 7: Thermal conductance G eq related to the heat flux exchanged between the thermal probe and the sample versus the input electrical tip voltage squared. The insert represents the normalized thermal conductivities of the four samples.

Figure 2 Figure 3

 23 Figure 2

Figure 6 400

 6 Figure 6

Figure 7 405

 7 Figure 7

Table 1 Typical composition of inner phase for preparation of microparticles

 1

	Nomenclature	Compound wt.-%
	ARD	n-eicosane / n-hexadecane / tetraethyl orthosilicate 48 / 48 / 4
	18	PVA-MDI microsphères* 100
	E2	n-eicosane / PVA-MDI microsphères* 68 / 32
	H/E	n-eicosane / n-hexadecane / PVA-MDI microspheres* 34 / 34 /32

Figure 1

Submitted to Polymer


Guillaume Galliero
email: guillaume.galliero@univ-pau.fr

Sebastian Volz

Thermodiffusion in model nanofluids by molecular dynamics simulations

published or not. The documents may come

Thermodiffusion in model nanofluids by molecular dynamics simulations

Introduction

Among the transport phenomena, thermodiffusion (Soret effect), which couples mass and heat fluxes [START_REF] De Groot | Non-Equilibrium Thermodynamics[END_REF] , is probably the one for which the modeling of the corresponding transport property is the less satisfying in condensed phase despite recent improvements in simple fluids [START_REF] Kempers | [END_REF][3][4] , polymers 5,6 or colloids [7][8][9] . It should be noted that this improvement is largely due to the fact that experimental measurements of this property, which are hard to handle because of thermal convection, have been strongly improved during the last decade and provide now consistent values using different experimental techniques 10 .

One possibility to improve the modeling and the knowledge of the microscopic mechanisms of this cross effect is the use of molecular dynamics (MD) simulations on model systems.

Such an approach has already proves its efficiency, at least concerning the microscopic mechanism responsible of thermodiffusion on rather academic systems, Lennard-Jones binary mixtures [11][12][13][14] and ternary ones 15 , molecular fluids [16][17][18][19] , associative mixtures 20 , reactive mixtures 21 , ionic systems 22 , polymer 23 , and fluids in porous media [24][25][26] . Concerning the modeling of thermodiffusion using MD results, we can mention in particular the work of Galliero et al. 27 on the mass effect in Lennard-Jones (LJ) mixtures and the work of Artola and Rousseau 28 which has shown the pertinence of some thermodynamic models to predict thermodiffusion in peculiar LJ mixtures.

Nevertheless, in dilute systems, thermodiffusion (or equivalently thermophoresis) was hardly accessible by molecular dynamics despite an obvious interest. This is due mainly because of the too poor statistics in such systems when using usual equilibrium or non-equilibrium MD algorithms 11,17,29 . In addition, when the solute particle is large compared to the solvent, as it is often the case in systems studied experimentally, the mass diffusion process becomes particularly low. So the simulation duration needed to perform the computation of thermodiffusion using MD becomes too large to be reasonably accessible despite the continuous increase in computational power.

Therefore, the main aim of this work is to provide a new MD procedure to allow the computation of infinite dilution (single particle) thermodiffusion in model systems. In particular, this new scheme will be applied to model nanofluids, i.e. dense fluids containing nanoparticles described as in Vladkov and Barrat 30 , which have raised a large interest in thermal science due to their abnormally large thermal conductivity 31,32 . It should be nevertheless mentioned that the results provided here deal with small non-metallic nanoparticles (diameter up to 4nm) in a LJ fluid and not "true" nanofluids for which thermal conductivity has exhibited this special behavior.

The new MD scheme relies on the estimation of the thermophoretic force on a solute particle induced by a fluid subject to a thermal gradient. This enables, assuming a Stokesian behavior of the system, the direct knowledge of the single particle thermal diffusion factor (characteristic of the thermodiffusion amplitude). The approach proposed here is noted Single Particle Thermodiffusion Algorithm (SPTA) in the following. In addition, through the MD computation of the mass diffusion coefficient, the single particle thermal diffusion coefficient is deduced as well.

In a first part is provided the thermodiffusion formalism as well as the way the fluid and nanoparticles are modeled. In addition, the new molecular dynamics SPTA proposed is described in details together with the technical details of the simulations. Then, in a second part, the SPTA is validated through a comparison with what can be extrapolated from usual boundary driven Non-Equilibrium MD (NEMD) results for finite concentration, first on simple LJ mixtures and second on nanofluids. Finally, the SPTA is applied on various nanofluids for which the nature of the solvent (quality and viscosity) and of the solute nanoparticle (size, mass and internal stiffness) is changed to study the various influences of these physical parameters on the thermodiffusion in such model systems.

Theory and Modelling

Thermodiffusion

The diffusive mass flux of component 1, J 1 , in a binary mixture experiencing a temperature gradient and in mechanical equilibrium is well described, close to equilibrium, by 1 :

( ) [ ] T D w w w D T ∇ - + ∇ - = 1 1 1 12 1 1 ρ J ( 1
)
where ρ is the density, w 1 the mass fraction of component 1, T the temperature, D 12 the mutual mass diffusion and D T the thermal diffusion coefficient (also termed thermophoretic mobility).

When the steady state is reached, the two contributions to the mass flux J 1 , due to mutual diffusion and to thermodiffusion, balance each other and the Soret coefficient can be defined as:

( ) ( )

1 1 12 1 1 1 1 1 1 1 1 T T D w x S D w w T x x T ∇ ∇ = = - = - - ∇ - ∇ (2)
where x 1 is the molar fraction of component 1. To quantify the amplitude of the thermodiffusion effect relatively to that of mass diffusion it is more convenient to use the dimensionless thermal diffusion factor defined as :

12 D D T TS T T T = = α (3)
In this work, the solute is defined as component 1 and the solvent as component 2. Thus, a positive α T means that the solute particles tend to migrate, relatively to the solvent, towards the cold areas.

Interactions potential

Solvent

To model the solvent-solvent interaction, the usual truncated Lennard-Jones 12-6 (LJ) potential is used on spheres: 12 6 4

LJ ij ij ij ij U r r σ σ ε         = -                 ( 4
)
where σ is the distance at which the potential is equal to zero (the "atomic diameter"), ε ij the potential depth and r ij the distance between atoms i and j. The cutoff radius, r c , has been set to 2.5σ. LJ molecular parameters have been taken to be those of an Argon like fluid: σ=3.405 Å, ε=996 J.mol -1 . In addition, except when stated, a molar mass m s =0.04 kg.mol -1 has been used.

Nanoparticles

The nanoparticles are quasi-spherical and composed of atoms, having the same molecular parameters than those of the solvent, which are distributed on a FCC crystal as in the work of Vladkov and Barrat 30 . Interactions between nanoparticle atoms are described by a LJ potential, Eq. ( 4), and those within the nanoparticules are linked to their nearest neighbours through a Finite Extensible Nonlinear Elastic (FENE) bonding potential 33 :

2 2 0 0 0.5 ln 1 ij FENE ij r U R R κ       = - -        (5)
where R 0 is a finite extensibility and κ a spring constant. In this work, except when stated, the FENE parameters were set to R 0 =1.5σ and κ=30ε/σ 2 , which allows simulations with relatively large time steps. It should be noted that, with such parameters, the distance between bonded atoms is peaked on 0.97σ.

In addition, to avoid clustering effects, for atoms belonging to different nanoparticles, a cutoff radius equal to 2 1/6 σ (WCA potential [START_REF] Hansen | Theory of Simple Liquids[END_REF] ) has been taken which means that nanoparticles are purely repulsive to each other.

To modify the quality of the solvent (or equivalently the wetting properties between the nanoparticles and the solvent), the energy parameters of the cross interaction is expressed by:

ij ij i j k ε ε ε = (6)
where k ij is a parameter taken between 0.5 (bad solvent) and 1.5 (good solvent) which is of great importance on the amplitude 14 , and even the sign 28 , of thermodiffusion in simple LJ fluids.

To estimate the diameter of a nanoparticle, as it presents some FCC facets, the radius of gyration, R g , was estimated during the simulation using:

( )

p N i m c i g N r r R p ∑ = - = 1 2 .
. 2 (7) where N p is the number of atoms in the nanoparticle and r c.m. the position of the nanoparticle centre of mass. Assuming a homogeneous mass distribution, the effective outer diameter of the nanoparticle, d NP , can be estimated using:

2 / 1 2 3 20         = g NP R d (8)
In this work, various nanoparticle sizes have been simulated, see table I. In real units, as σ=3.405 Å, the range of nanoparticle diameter tested goes from 0.8 to 4 nm which remains small compared to usual nanofluids. In fact, it should be noted that a more usual diameter of 10 nm (which remains nevertheless a small value), would imply around 2.5 10 4 atoms in the nanoparticle (

3 / 1 p NP N d ∝
) and at least 2.5 10 5 in the fluid. This would lead to an inaccessible CPU time cost of at least 10 5 hours on a recent supercomputer to obtain one value of the thermal diffusion factor by MD simulation.

Measurement of mass diffusion

In order to deduce the thermal diffusion coefficient from the thermal diffusion factor, one needs the mass diffusion coefficient, see Eq. (3). This coefficient has been obtained during separate equilibrium MD simulations from the mean-squared displacement of the centre of the mass of all fluid molecules, . .m c i r ∆ , using:

( ) ( ) t r N m x m x m x x x D m c i 6 1 1 2 . . 2 2 2 2 2 2 1 1 2 1 12 ∆         + = ( 9
)
where N is the total number of particles, t the time and m i the mass of particle i.

Computation of the thermodiffusion: usual Non Equilibrium MD

To compute thermodiffusion in concentrated systems the Non Equilibrium Molecular Dynamics approach of Müller-Plathe and Reith 29 has been used. In this scheme, a biperiodical heat flux is imposed to the simulation, and after a transient state, the thermal diffusion factor can be deduced using eqs. (2-3) and the measured temperature and molar fraction profiles. An example, on a half simulation box, of the profiles obtained is shown on figure 1. It is interesting to note that, as expected from Eq. ( 2), the molar fraction gradient can be non linear in nanofluids when the nanoparticles molar fraction is low and α T is large as in figure 1. In addition, as a by-product, the effective thermal conductivity can be deduced using the Fourier's Law and the temperature profile.

To impose the heat flux, the simulation box is divided into N s slabs, of identical thickness and volume. Slabs 1 and N s are defined as the "hot" ones and slabs N s /2 and N s /2+1 as the "cool" ones. The heat flux is generated by an exchange, at a given swapping frequency for both half simulation boxes, of the solvent particle that has the lowest kinetic energy (from the "hot" slab) with the one that has the highest kinetic energy (from the "cold" slab). The four slabs where the exchanges are performed, as well as their first neighbors, have been discarded to measure the profiles.

New scheme to compute the single particle thermal diffusion factor.

When the concentration of the solute becomes too small, the usual MD approaches to estimate thermodiffusion, such as the one presented in the previous section, are inefficient mainly because of the too poor statistics. In addition, when the solute particle is large compare to the solvent molecule, the mass diffusion process (associated to the thermodiffusion characteristic time [START_REF] De Groot | Non-Equilibrium Thermodynamics[END_REF] ) is very slow and the MD simulations need to be performed for a very long duration.

Therefore, in order to study the thermodiffusion in dilute nanofluids using molecular dynamics a new approach is proposed to deduce the single particle (SP) thermal diffusion factor, SP T α , which corresponds to the infinite dilution limit.

Let us consider a single particle (solute) in a fluid (solvent) subjected to an established thermal gradient. This particle experiences a thermophoretic force, F T , which in turns will induce a velocity drift of the particle, v T , proportional to the temperature gradient in dilute systems:

T D SP T T ∇ - = v (10)
where SP T D is the single particle thermal diffusion coefficient.

If we consider that this velocity drift is small enough (the particle Reynolds number is small in the cases treated here) and that the time of the analysis is long enough (typically on the order of a few nanoseconds), inertial forces can be neglected. In that case, the hydrodynamic force related to friction (Stokes drag) acting on the particle, F h , is linearly related to the velocity of the sphere, v, through:

v F ξ - = h (11)
where ξ is the friction coefficient (usually expressed as proportional to the solvent viscosity and the hydrodynamic radius of the sphere [START_REF] Bird | Transport Phenomena[END_REF] ).

The force balance implies

0 = + T h F F
which induces that the velocity drift due to the thermophoretic force writes as:

ξ T T F v = (12)
Moreover, if we assume that a Stokes-Einstein law [START_REF] Cussler | Diffusion: Mass Transfer in Fluid Systems[END_REF] applies between the single particle mass diffusion, SP D 12 , and friction, we have:

ξ T k D B SP = 12 (13)
Finally combining Eq. (10) with Eqs. (12-13), we obtain that the single particle thermal diffusion factor, SP T α , in such systems can be expressed as:

T k D TD B T SP SP T SP T ∇ - = = F 12 α (14)
Such an expression implies that the measurement of the thermophoretic force acting on the particle for a given thermal gradient provides a straightforward estimation of the thermal diffusion amplitude in such dilute systems. In addition, if the single particle mutual diffusion coefficient, SP D 12 , is known, the single particle thermal diffusion coefficient can be deduced from:

T T k D D B T SP SP T ∇ - = F 12 ( 15
)
In order to measure the thermophoretic force, F T , acting on the particle using molecular dynamics simulations, a simple scheme is proposed.

First, an initial system composed of two particles and the fluid is constructed. The two particles that are centered around L x /4 and 3L x /4, L x being the size along x of the simulation box, are attached to the reference frame of the simulation box through a harmonic potential with a spring constant: 2 30 σ ε κ

C h = ( 16
)
where C is a numerical constant below unity and set to ≈0.1 in order to weakly perturb the system. Then, the NEMD approach described in the previous section is applied to generate a bi-periodical linear thermal gradient to the simulation box in the direction x. After a transient state, the location of the centres of mass of the two particles will be displaced relatively to their point of fixation because the thermophoretic force applies on them. The measure of this displacement for both particles, ∆x, provides an estimation of the amplitude of the thermophoretic force, F T , that applies on each particle through:

x F h T ∆ = κ (17)
To be consistent with the sign convention of the thermal diffusion factor, F T is counted positively when the displacement goes against the thermal gradient. It should be noted that this scheme is similar to the experimental approach used in Albanese et al. [START_REF] Albanese | [END_REF] to study thermal radiation forces produced by heat flow through solid slabs in nonisothermal liquids.

Technical details of the simulations

The simulation box is cubic. A reduced timestep, δt * =δtε 1/2 m -1/2 σ -1 , varying from 0.001 to 0.004 has been used. To integrate the equation of motion the Velocity Verlet algorithm has been applied. In a first step, to obtain the desired temperature and pressure, Berendsen thermostat and barostat were used 38 during at least 10 5 timesteps. Then, nonequilibrium molecular dynamics runs have been performed to collect the results. To avoid nonlinearity in the temperature profile and phase transitions during the NEMD simulations, a weak reduced energy flux, e.g. in figure 1,

2 / 1 2 / 3 3 * m J J U U - = ε σ
, ranging from 0.02 to 0.006, has been applied (which corresponds to an exchange period of 50-150 timesteps).

Simulations on nanofluids have been performed on systems composed of 4000 to 21794 atoms (solute +solvent) depending on the size of the nanoparticle to avoid as much as possible finite size effects and interactions between nanoparticles in the case of the new scheme proposed. All nanofluids simulations have been performed at T*=k b T/ε=1 and P*=Pσ 3 /ε=1 which corresponds to a dense liquid state (ρ*=Nσ 3 /V≈0.8 in a pure LJ fluid for such conditions). After the transient state, runs of 2 to 4 10 7 timesteps were used to obtain the thermal diffusion factor for classical NEMD simulations. Concerning the SPTA, four independent runs of 5-20 10 6 each have been carried out.

To generate the initial nanofluid systems a FCC bulk arrangement of atoms is constructed.

Then, random atoms are chosen as centre of nanoparticles and, avoiding the overlap, the atoms located within the radius of the nanoparticles are linked to their first neighbours by the FENE bond.

Concerning simulations on simple fluids mixtures, as previously done 14 , they have been performed on systems composed of 1500 LJ particles for duration varying from 2-8 10 7 nonequilibrium timesteps.

For the computation of mutual diffusion, four independent runs of 10 6 timesteps have been performed at equilibrium for each system. As expressed in Eq. ( 9), the slope of the meansquared displacement versus time was calculated during the linear regime.

Results

Preliminary study

Simple fluids

The first point was to analyse the validity of the proposed method to evaluate the single particle thermal diffusion factor, SP T α , in simple dense LJ binary mixtures. In fact, this approach is based on "macroscopic" liquid state assumptions and is therefore, a priori, limited to large particle immersed in a liquid. Nevertheless, it has been shown that the definition of a Stokes-Einstein like law at this scale makes sense [39][40] with an appropriate friction factor, ξ. In addition, the final relation of Eq. ( 14) derived for SP T α only implies the thermal gradient and the thermophoretic force and so does not need an explicit formulation of the friction factor ξ.

The first mixture tested is composed of species having the same molecular parameters (σ 1 =σ 2 and ε 1 =ε 2 ) except that m 1 /m 2 =10 (isotope like mixtures) i.e. ideal ones in the thermodynamic sense. Two thermodynamic conditions have been studied, one corresponding to a dense liquid at T * =1 and ρ * =0.8 and the second one to a dense supercritical gas at T * =1.686 and ρ * =0.477 that was studied previously 14 . In addition, results have been compared with the correlation recently proposed 27 to predict thermodiffusion in such systems.

Results shown on figure 1 clearly reveal that, the proposed approach to compute the infinite dilution values of α T is consistent with what could be extrapolated from the values obtained by classical NEMD simulations. It is interesting to note that, despite the underlying assumptions see section 2.5., this new approach is able to provide reasonable result for both molar fraction limits (x 1 →0 and x 1 →1) and for both states. Besides, the results provided by the correlation provided in Ref. 27 are in agreement with the MD ones despite non-negligible deviations (up to 15 %).

The second mixture is a n-decane/methane (modelled by LJ spheres) one at T * =1.686 and ρ * =0.477 previously studied 15 . By using the SPTA, it has been found that Besides, it is worth noting that the possibility of obtaining the SP T α coefficient for such fluid mixtures, for both infinite dilution limits, could be valuable as well because recently proposed models to estimate thermodiffusion [41][42] uses these data to predict the thermodiffusion on the whole molar fraction range.

Nanofluids

To test the validity of the SPTA when applied on nanofluids, which are the main concern of this work, nanofluids have been simulated for various nanoparticles volume fraction and for the two smallest sizes studied in this work, d NP /σ=2.4 and 4.03 at T*=1 and P*=1 (approximately T=120 K and P=42 MPa for the molecular parameters chosen). For such conditions, we have computed that the solvent viscosity η s =1.8±0. 1 10 -4 Pa.s and the solvent self diffusion D s =3.6±0.1 10 -9 m 2 .s.

Results provided on Fig. 3 confirm that the proposed SPTA is able to provide an estimation of SP T α which is consistent with what could be extrapolated from classical NEMD results. This is true even it seems that the SPTA values slightly underestimate extrapolated values as in the previous section. This may come from the fact that for very low concentration α T becomes independent of concentration as particles do not see each other as noticed experimentally by Ning et al. 43 or because of a weakness in the underlying theory as claimed by Brinquier and Bourdon 44 .

Besides, whatever the size and the volume fraction, the nanoparticles tend to migrate towards the cold area; α T being always positive and rather important as seen in Fig. 3. In addition, α Τ decreases strongly with the nanoparticles volume fraction. More precisely, it seems that the values tend towards zero for large nanoparticles concentrations, i.e. >30% in volume. It was not possible to see whether or not a sign inversion may occur for larger volume fraction because such concentrations were not accessible by the scheme adopted to construct initial systems. The noted α T evolution with respect to the nanoparticles concentration is consistent with the experiments and model proposed by Ning et al. 43 to describe the hard-sphere suspensions and the work of Rauch et al. 5 on polymer suspensions. Nevertheless the decrease noted here is larger than what can be expected from both models.

As a by-product of the NEMD simulations, the thermal conductivity can be obtained for the nanofluids studied. Figure 4 clearly shows that thermal conductivity tends to decrease with the nanoparticles volume fraction, this decrease being slightly more pronounced for the smallest nanoparticle. These trends can be explained by two antagonist effects, the thermal shortcuts induced by the presence of nanoparticles that are more conductive than the fluid 45 (volume effect), and the Kapitza resistance occurring at the interface nanoparticle-fluid 30 (surface effect). For the studied systems, the surface effects predominate over the volume ones.

It should be noted that the decrease of the thermal conductivity λ with the nanoparticles volume fraction is in contradiction with what is generally found in experimental reports on nanofluids [31][32] whereas most of MD simulations confirm our trend 30,46 . Nevertheless, the here studied systems are far from being similar to those that are experimentally analysed. In particular the nanoparticles are not metallic ones (and are hence less conductive), the nanoparticles are smaller than in experiments and the fluid is by far simpler than water. For instance, it does not involve association.

Besides, it is worth noticing that the thermal conductivity measured here is an effective one.

Thus, it takes into account the transport of energy through the Dufour effect (the symmetric effect to the Soret one). Therefore, if D T is positive, this effect tends to slightly decrease the effective thermal conductivity 1 when the concentration gradient is established, but no more than a few percents.

Single particle thermodiffusion

In this section all simulations using the SPTA have been performed so that T*=1 and P*=1 in the fluid.

Influence of the size

Among the parameters that may affect the single particle thermodiffusion is the size (i.e the radius) of the nanoparticle, as already noticed in section 3.1.2. Simulations have been performed for d NP /σ ranging from 2.4 to 11.59 as reported in Table I.

From the values provided in the appendix, it first appears that the sign of SP T α is always positive whatever the size (i.e. nanoparticle tends to migrate towards the cold areas) and the obtained values are from one to two orders larger than the one for non interacting particles 8 (α T =1). In addition, SP T α amplitude increases with the diameter in a monotonic way and roughly linearly with d NP .

If SP T D is deduced from the computation of SP D 12 and the results obtained on SP T α , further information can be gathered as provided in Fig. 5. The most striking result is that SP T D seems to be independent of the size of the nanoparticle, with a value that is always close to 0.7 10 -10 m 2 .s -1 .K -1 . It should be noted that such behaviour is consistent with findings on dilute polymer [47][48] and nanodroplets 49 . Besides, as shown by figure 5, SP D 12 is roughly proportional to the inverse of the radius of the nanoparticle which is consistent with a Stokes-Einstein law, Eq. ( 13), as ξ is proportional to the hydrodynamic radius.

Influence of the nature of the solvent

In this section and the following, we have focused our work on a nanoparticle composed of 55 atoms (d NP /σ=4.03).

To quantify the effect of the nature of the solvent on the single particle thermodiffusion, we have first evaluated the influence of the quality of the solvent by varying k ij from 0.5, which corresponds to a bad solvent, to 1.5, i.e. a good solvent, appearing in Eq. (7). In a second step, we have varied the viscosity of the solvent, without modifying its affinity with the solute, through the modification of m s from m s /m NP =1 to 10 -2 , as solvent viscosity scales [START_REF] Boon | Molecular Hydrodynamics[END_REF] as m s 1/2 .

Figure 6 shows that the quality of the solvent affects both thermodiffusion and mass diffusion.

SP T

D tends to increase with the quality of the solvent (as SP T α , see appendix). This is in agreement with previous MD works on different systems 23,28 and with a lattice model described in Ref. 51. But, in the range of the k ij values tested, no sign inversion has been noticed.

SP D 12 has an opposite behaviour as proved by figure 6. It decreases, nearly linearly, with the quality of the solvent. Such behaviour can be understood by the fact that increasing k ij will imply that more solvent particle will be structured or adsorbed around the nanoparticle. As a consequence the hydrodynamic radius will be larger and the mass diffusion coefficient will decrease in consequence.

As clearly shown on Fig. 7, both SP T D and SP D 12 strongly decrease with the viscosity of the solvent. More precisely they are both nearly proportional to η s -1 (or equivalently to D s ). It should be noted that, concerning SP T D in dilute polymer solutions, such dependence to the viscosity with a η s -1 behaviour, has been emphasized in a recent experimental/theoretical work 48 .

These behaviours are consistent with the fact that increasing the solvent viscosity leads to an increase of the friction coefficient ξ that is usually assumed [START_REF] Bird | Transport Phenomena[END_REF] to be proportional to η s . Hence, this increase leads to a decrease of SP D 12 and SP T D , Eqs. (13,15), if we assume that the thermophoretic force, F T , is weakly dependent on m s . Concerning this last point, the fact that SP T α slightly decreases with m s , as shown in the appendix, implies that F T is not perfectly independent of m s for the system studied. This may be explained by the fact changing m s , apart from modifying η s , does affect the thermal conductivity of the solvent (λ s is proportional 50 to m s -1/2 ). Hence, we can suspect that F T is affected by a modification of the difference (or the ratio) between λ s and λ NP , as explicitly done in some theoretical models to predict thermodiffusion 8,52 .

Influence of the nature of the particle

To analyse a possible influence of the internal dynamic of the nanoparticle, two parameters were varied while keeping those of the solvent constant: the mass of the nanoparticle atoms, m NP (from m NP /m s =1 to 50), and the strength of the intramolecular interaction through the modification of the "spring" constant κ (ranging from 10ε/σ 2 to 90ε/σ 2 ) appearing in the FENE potential, Eq. ( 5).

Contrary to what occurs in simple LJ fluids 27 , results shown on Fig. 8 indicate that SP T D seems to be independent of m NP , whereas SP D 12 decreases weakly with m NP . This last result indicates that the particles are too small to behave like Brownian ones (for which SP D 12 is independent of m NP ) as noted as well in previous works for similar systems 40,53 .

Another interesting point is that by changing the mass m NP , λ NP is modified because, as the FENE potential is weakly anharmonic, λ NP should roughly be proportional to (κ/m NP ) 1/2 .

Therefore the ratio between λ NP and λ s evolves in the same way than in the previous section, i.e. proportionally to (m NP /m s ) -1/2 . This may explain why SP T α evolves versus m NP /m s in the same way than in the previous section, as shown in the appendix.

As shown in Fig. 9 and contrary to the mass m NP , the nanoparticle internal stiffness, κ, which

Summary and Conclusions

In this work, we provide a new simple scheme to compute, using Non Equilibrium Molecular Dynamics simulations, the Single Particle (infinite dilution) thermodiffusion in dense fluids which is not accessible in a reasonable amount of time using standard MD scheme.

In a first part, by a comparison with extrapolations from finite concentration results using standard NEMD simulations, this scheme is shown to be efficient in ideal and non ideal Lennard-Jones binary mixtures as well as in model nanofluids. These last being described by a LJ fluid + quasi-spherical nanoparticles composed of atoms distributed on a FCC crystal and interacting through LJ potential plus FENE bonding with the nearest neighbours. For these systems it appears that the thermodiffusion and thermal conductivity amplitudes decrease with nanoparticles concentration for both nanoparticles sizes tested.

Then, using the proposed scheme, the influence of the nature of the solvent and of the particle on the single particle thermodiffusion is estimated in these model nanofluids. It appears that, in all cases studied here, the nanoparticle tends to migrate towards the cold area and that both the nature of the solvent and the nature of the particle affect the thermodiffusion amplitude.

The most striking result is that the thermal diffusion coefficient is shown to be independent of the size of the nanoparticle on the range tested (d NP = 0.8-4 nm), whereas the mass diffusion coefficient is roughly proportional to the inverse of the nanoparticle diameter (consistent with the Stokes-Einstein law). Concerning the influence of the nature of the solvent, it appears that

15 .

 15 Thus, on the LJ mixtures studied (ideal or not, liquid or supercritical), the SPTA is able to provide a consistent estimation of SP T α .

affects 45 λ 12 ,

 4512 NP , has a non negligible influence on SP T D . The larger the spring constant κ is, the larger SP T D (and λ NP ). In addition SP T α increases with κ, see Appendix. This result clearlyshows that the internal degrees of freedom of the nanoparticle affects thermodiffusion whatever is the coefficient to quantify it, its behaviour with κ is similar to the one of SP T D in a less pronounced manner as illustrated by Fig.9. In fact, increasing κ generates two effects: the nanoparticle is more conductive and the radius of the nanoparticle slightly decreases, e.g. d NP /σ=4.26 when κ=10ε/σ 2 and d NP /σ=3.79 when κ=90ε/σ 2 . Thus, the hydrodynamic radius decreases with κ and so SP D 12 increases with κ when assuming a Stokes-Einstein like behaviour.

Figure captions:Figure 1 :Figure 2 :

 12 Figure captions: Figure 1 : Stationary temperature and nanoparticles molar fraction profiles in a half simulation box (d NP /σ=4.03 and nanoparticles volume fraction equal to 7.31 %)

Figure 3 :

 3 Figure 3: Thermal diffusion factors versus nanoparticles volume fraction for two nanoparticle sizes, d NP /σ=2.4, circles, d NP /σ=4.03, squares. Open symbols correspond to usual NEMD results, full symbols to SPTA ones.

Figure 4 :

 4 Figure 4: Thermal conductivity in nanofluids for two nanoparticles sizes, d NP /σ=2.4, circles,

Figure 5 :

 5 Figure 5: Single particle thermal diffusion coefficient and mutual diffusion versus the inverse

Figure 6 :

 6 Figure 6: Influence of the quality of the solvent on SP T D , circles, and SP D 12 , squares.

Figure 7 :

 7 Figure 7: Influence of the viscosity of the solvent on SP T D , circles, and SP D 12 , squares.

Figure 8 :

 8 Figure 8: Influence of the mass of the atoms constituting the nanoparticle on SP T D , circles,

Figure 9 :

 9 Figure 9: Influence of the spring constant of the FENE potential, κ, on SP T D , circles, and

  Figure 1

  Figure 6

  Figure 7

  Figure 9

  Single particle thermal diffusion factor results for the various nanofluids studied in this work.

	d NP /σ	m NP (kg.mol -1 ) m s (kg.mol -1 )	κσ 2 /ε	k ij	α	SP T
	2.4	0.04	0.04	30	1	8.4±1.5
	4.03	0.04	0.04	30	1	20.1±4.1
	5.98	0.04	0.04	30	1	35.4±7.4
	7.73	0.04	0.04	30	1	52.8±10.5
	9.75	0.04	0.04	30	1	62.5±15
	11.59	0.04	0.04	30	1	91±50
	4.03	0.04	0.02	30	1	19.8±3.9
	4.03	0.04	0.008	30	1	24.0±6.4
	4.03	0.04	0.004	30	1	23.5±7.6
	4.03	0.04	0.002	30	1	26.1±7.8
	4.03	0.04	0.0008	30	1	27.6±9.1
	4.03	0.08	0.04	30	1	22.2±4.3
	4.03	0.2	0.04	30	1	23.5±4.5
	4.03	0.4	0.04	30	1	23.4±6.8
	4.03	0.8	0.04	30	1	25.5±6.5
	4.03	2.0	0.04	30	1	28.1±6.5
	SP T D (and SP D 12 ) is nearly proportional to the inverse of the solvent viscosity and increases 4.03 0.04 0.04 10 1 13.1±2.5
	with the quality of the solvent. Besides, it is shown that SP 4.03 0.04 0.04 20 D is unaffected by the mass of the 1 18.5±4.5 T 4.03 0.04 0.04 60 1 23.75±4.2
	nanoparticle atoms whereas SP D 12 slightly decrease with m NP (i.e. nanoparticles studied here 4.03 0.04 0.04 90 1 25.6±6.3
	are not Brownian particles). Finally, results indicate that both SP T D and SP 4.03 0.04 0.04 30 0.5 D 12 increase with the 8.5±2.9 4.03 0.04 0.04 30 0.75 16.7±4.5
	value of the internal stiffness of the nanoparticle (through the modification of the spring 4.03 0.04 0.04 30 1.25 23.4±4
	constant κ). 4.03	0.04	0.04	30	1.5	25.4±3.5

Appendix :

Table : Table I :

 :I Relation between numbers of atoms inside a nanoparticle and its effective diameter for the various sizes tested.

	N p	13	55	177	381	767	1289
	d NP /σ	2.4	4.03	5.98	7.73	9.75	11.59

Acknowledgments

We thank the TREFLE laboratory in Bordeaux and the CINES in Montpellier for having provided a large part of the computer time needed for this work. F. Montel is acknowledged for fruitful discussions.