Françoise Fenouillot
email: f.fenouillot@insa-lyon.fr

Françoise Méchin
email: francoise.mechin@insa-lyon.fr

Fernande Boisson

Pierre Alcouffe

Tomasz Pokropski

Tasnim Kallel

Mahdi Mnif

Coarsening of nanodomains by reorganization of polysiloxane segments at high temperature in polyurethane / , -aminopropyl polydimethylsiloxane blends

Keywords: polyurethanes, polysiloxanes, melt reaction, copolymer

The reaction in bulk at high temperature of ,-aminopropyl oligodimethylsiloxane and thermoplastic polyurethane (TPU) allowed observing interesting behaviour. Mixing at 200°C first involved dissociation of urethanes and splitting of polyurethane chains followed by reaction of the released isocyanates with amino end-groups of the oligosiloxane. At this stage, a copolymer was formed which morphology consisted of a very fine dispersion of the polysiloxane domains at the nanoscale (20 nm) with a narrow size dispersity. The polymer blend was perfectly transparent. Increasing the reaction time resulted in a significant coarsening of the morphology and a consequent loss of transparency. The reason for such a morphology evolution has been elucidated. The progressive formation of alkyl-alkyl urea linkages at the expense of aryl-alkyl bonds obtained earlier in the process caused an increase in the average number of successive polydimethylsiloxane (PDMS) blocks that organized in larger domains. The average number of consecutive polysiloxane segments was found to evolve from ~1.5 in the first 10-15 minutes to a value of 3-5 at the end of the reactive process.

Introduction

A large number of works is dedicated to the association of polyurethanes with polysiloxanes, either in simple blends prepared by a classical melt mixing technique, or using polydimethylsiloxane (PDMS) oligomers as reactive precursors in the synthesis of segmented polyurethane or polyurethane-urea copolymers [START_REF] Knaub | New reactive polymer blends based on poly(urethaneureas) (PUR) and polydisperse polydimethylsiloxane (PDMS): control of morphology using a PUR-b-PDMS block copolymer[END_REF][START_REF] Hill | Development of wear-resistant thermoplastic polyurethanes by blending with poly(dimethyl siloxane). I. Physical properties[END_REF][START_REF] Bremner | Development of wear-resistant thermoplastic polyurethanes by blending with poly(dimethyl siloxane). II. A packing model[END_REF]. The main objective is usually to take advantage of the specific surface properties linked to the introduction of silicone known to be highly hydrophobic. The applications are materials with minimal biological adhesion used in medical devices where the antithrombogenicity of the surface is needed (artificial heart or vascular grafts, implants…), or in coatings where the settlement of marine organism [START_REF] Wynne | Poly(dimethysiloxane)-ureaurethane copolymers. Synthesis and surface properties[END_REF] must be inhibited. Also, better hydrolytic stability can be achieved in the presence of silicone, as well as flame retardancy with a halogen-free component. Finally, gas permeation can be tuned depending on the structure of the material [START_REF] Park | Gas separation properties of polysiloxane/polyether mixed soft segment urethane urea membrane[END_REF][START_REF] Dolmaire | Modification of a hydrophilic linear polyurethane by crosslinking with a polydimethylsiloxane. Influence of the crosslink density and of the hydrophobic/hydrophilic balance on the water transport properties[END_REF].

Copolymers are usually obtained from the reaction of hydroxyalkyl-or aminoalkylterminated oligosiloxanes with diisocyanates and diols. The reaction is most often carried out in solution in mixtures of solvents that must be carefully selected to maintain the homogeneity of the reactive medium. Indeed, the large solubility difference between the non-polar silicone and the highly polar urethane segments tends to induce a macrophase separation. Therefore the mechanical properties obtained when using only PDMS as soft segment are often rather poor, also because of possible end-group instabilities of the hydroxy-terminated PDMS oligomers. In order to improve the miscibility of the reactive systems and the material final properties, the use of mixtures of PDMS oligomers with a second macrodiol such as poly(propylene oxide), poly(tetramethylene oxide) or preferably poly(hexamethylene oxide) has been recommended by several authors [START_REF] Martin | Polydimethylsiloxane/polyether-mixed macrodiol-based polyurethane elastomers: biostability[END_REF][START_REF] Sheth | Structure-property behavior of poly(dimethylsiloxane) based segmented polyurea copolymers modified with poly(propylene oxide)[END_REF][START_REF] Hernandez | Microstructural organization of three-phase polydimethylsiloxane-based segmented polyurethanes[END_REF][START_REF] Hernandez | A comparison of phase organization of model segmented polyurethanes with different intersegment compatibilities[END_REF][START_REF] Choi | Influence of soft segment composition on phaseseparated microstructure of polydimethylsiloxane-based segmented polyurethane copolymers[END_REF].

The particular surface properties of poly(urethane-urea) copolymers containing silicone are obtained through the segregation of incompatible and less polar soft silicone domains toward the surface to minimize the interfacial tension with air [START_REF] Chen | Surface composition of a series of dimethylsiloxaneurea-urethane segmented copolymers studied by electron spectroscopy for chemical analysis[END_REF][START_REF] Zhuang | Determination of the distribution of poly (dimethylsiloxane) segment lengths at the surface of poly[(dimethylsiloxane)-urethane]segmented copolymers by Time-of-Flight Secondary Ion Mass Spectrometry[END_REF]. The bulk compositional features of the copolymers are less understood. A complex microphase separation of PDMS segments from a mixture between a second soft segment (polyoxypropylene) and the hard segments was observed but not clearly explained [START_REF] Kinning | Bulk, surface, and interfacial characterization of silicone-polyurea segmented copolymers[END_REF]. Actually, the immiscibility of the three types of segments should involve a ternary structure but this point was not always elucidated in the literature.

Recently, Hernandez et al [START_REF] Hernandez | Microstructural organization of three-phase polydimethylsiloxane-based segmented polyurethanes[END_REF] investigated the microstructural organization of PDMS-based polyurethanes synthesized from mixed macrodiols, especially by SAXS and DMA; they found that the polyurethanes under study indeed organized into three phases: hard domains, silicone domains, and a mixed phase composed of the second macrodiol, PDMS end-group segments, and dissolved hard segment sequences. In polyurethanes, the size of the hard segment domains is often estimated to be 3-5 nm [START_REF] Cho | Phase structure of poly(dimethylsiloxane-ureaurethane) segmented copolymers as observed by solid state nuclear magnetic resonance spectra[END_REF][START_REF] Choi | Microstructural organization of polydimethylsiloxane soft segment polyurethanes derived from a single macrodiol[END_REF][START_REF] Mclean | Tapping-mode AFM studies using phase detection for resolution of nanophases in segmented polyurethanes and other block copolymers[END_REF]. On the other hand, the size of the silicone-rich domains varies from a few tens of nanometers to more than 100 nm depending on the PDMS content and molar mass [START_REF] Wang | Synthesis and characterization of organosiloxane modified segmented polyether polyurethanes[END_REF][START_REF] Kozakiewicz | Investigations of morphology of moisture-cured silicone-urethanes[END_REF].

By contrast to the typical solution method usually employed in the literature, here the copolymers will be synthesized in the melt state starting from a thermoplastic polyurethane and amino-functional oligodimethylsiloxanes. More precisely, this work focuses on the reactions between the polyurethane and ,-aminopropyl polydimethylsiloxane at high temperature (200°C). A special emphasis is put on the evolution of the morphology of the copolymer during reactive processing. For comparison, blends obtained by mixing the same TPU with a commercial polydimethylsiloxane-urea copolymer have been prepared and characterized.

Experimental

Materials

The thermoplastic polyurethane (TPU) elastomer was a commercial polyester-based TPU, Elastollan 685A (named E685), kindly supplied by BASF, with poly(butylene adipate) soft segments and hard segments based on 4,4'-diphenylmethane diisocyanate (MDI), and 1,4butanediol (BDO) (Figure 1). The ,-aminopropyl oligosiloxanes were supplied by Wacker Chemie. Three different molar masses were used referenced respectively as PDMS15D, PDMS40D and PDMS130D. They were all liquid at room temperature. A thermoplastic silicone elastomer, Geniomer 80 (G80), given by Wacker was also used. It was a polydimethylsiloxaneurea copolymer obtained by the stoichiometric reaction of an ,-aminopropyl polydimethylsiloxane with n M =2.75.10 3 g/mol ( 1 H NMR determination) with isophorone diisocyanate (IPDI) (Table 1, Figure 1). The molar masses of G80 and E685 are also given in Table 1.

Blending method

The thermoplastic polyurethane was dried at 90ºC for at least 4 hours. Geniomer 80 was dried under vacuum (200 mbar) at 30ºC for at least 4 hours. PDMS oligomers were dried during 9 hours at 90°C. The extrusion was performed in a DSM 15 cm 3 twin-screw mini-extruder. The extruder was charged with 13.5 g of materials corresponding to the maximum volume (15 cm 3 ) of the chamber.

Blends of polyurethane E685 with ,-aminopropyl polydimethylsiloxane oligomers.

Half of the calculated amount of TPU pellets was charged into the mini-extruder heated up to 200 ºC at a screw speed of 100 rotations per minute. The mini-extruder drive was then stopped, the silicone building blocks were injected into the chamber using a syringe, and the second half of the pellets was charged. Samples were extruded after different mixing times.

2.2.2.Blends of polyurethane E685 with polydimethylsiloxane-urea copolymer G80

The thermoplastic silicone and TPU pellets were dry blended and then fed at the same time in the extruder. In this case, the residence time in the extruder was set to 10 minutes.

The extrudates were always immediately quenched in cold water at the exit of the die.

Nuclear Magnetic Resonance Spectroscopy (NMR)

High-resolution liquid NMR spectroscopy was carried out with a Bruker DRX 400 spectrometer operating at 400 MHz for 1 H. Spectra were obtained with a 5 mm-QNP probe at 298K or 323K. Deuterated tetrahydrofuran was used as solvent. Chemical shifts values are given in ppm with reference to internal tetramethylsilane.

Transmission Electron Microscopy (TEM)

Transmission electron microscopy analyses were carried out at the ''Centre Technologique des Microstructures de l'Université Claude Bernard Lyon 1" on a Philips CM120 microscope operating at 80 kV. The samples were ultra-microtomed at low temperature (-100°C) by a LEICA ultra-microtom, to have thin sections with a thickness around 80 nm.

Size Exclusion Chromatography (SEC)

The blends were analyzed by SEC using a Viscotek device equipped with a Viscotek VE 1121 pump and using 3 Waters columns HR 5E. The eluent used was tetrahydrofuran (THF) at a flow rate of 1 mL/min. The detectors included differential refractometer Viscotek VE3580, and a Viscotek Dual T60 detector (combining a right angle light scattering detector and a viscosimeter). Low molar mass components were analyzed using special columns (Viscotek G3000H and G2000H) with refractive index detector. Average molar masses were calculated using 3 detectors calculation method. As amino-terminated oligodimethylsiloxane could not be directly detected and analyzed, both because they partially adsorb on the columns and because their refraction index is quite close to that of THF, they were initially modified using ptolylisocyanate in THF solution. The same method was used for the reactive blends whenever these oligomers could be present.

Results and discussion

Polyurethane/, -aminopropyl oligosiloxane blends

3.1.1.Reactive process

The mixing process of the polyurethane with ,-aminopropyl oligosiloxane at high temperature does not consist only in simple mixing but is likely to involve several chemical reactions. First, the thermal reversibility of urethane bond above a temperature depending on its chemical structure (usually above 150-170°C) is indeed well-known [START_REF] Yang | Thermal degradation of urethanes based on 4,4'diphenylmethane diisocyanate and 1,4-butanediol (MDI/BDO)[END_REF][START_REF] Joel | Thermal dissociation of urethanes studied by FTIR spectroscopy[END_REF] and was shown for example to induce a strong broadening of the molar mass distribution when heating an initially monodisperse polyurethane [START_REF] Camberlin | Synthesis and DSC study of model hard segments from diphenylmethane diisocyanate and butanediol[END_REF]. The stability scale for urethane bonds is given as follows:

aryl-NHCOO-aryl < alkyl-NHCOO-aryl < aryl-NHCOO-alkyl < alkyl-NHCOO-alkyl 120°C 180°C 200°C 250°C
This reversibility phenomenon occurs under classical TPU processing conditions, although it was shown from theoretical kinetics calculations that only a tiny fraction of the urethane bonds were actually split, for example about 1% at 235°C [START_REF] Endres | The kinetics of the thermal decomposition of thermoplastic polyurethane elastomers under thermoplastic processing conditions[END_REF]; however a classical commercial TPU typically bears 3-400 urethanes per chain and even such a small proportion of bond scission can have non negligible consequences on the average molar mass.

Moreover, if the TPU is processed in the presence of amino-terminated oligomers, depending on the experimental conditions various NCO/OH or NCO/NH 2 recombinations can occur, and the formation of more or less randomized copolymers could therefore be expected. In the present work at temperatures above 170 ºC, urethanes can split into isocyanate-terminated and hydroxyl-terminated chains. The reaction of isocyanates with amino groups of the oligosiloxane is favored rather than recombination with hydroxyl-terminated chains. This assumption is supported by a quite complete study involving model urethanes and various reactive monofunctional compounds by Lu et al [START_REF] Lu | Reactivity of common functional groups with urethanes: models for reactive compatibilization of thermoplastic polyurethane blends[END_REF] where the authors demonstrated that among all used functions, the amines (primary or secondary) were the most suitable to allow reactive compatibilization of TPUs with other functional but non miscible polymers. Therefore, under the used mixing conditions the formation of urethane-polysiloxane block or segmented copolymers is expected with a final structure depending on the composition of the blend (Scheme 1).

Because of such differences in reaction rates and stoichiometric imbalance, some excess of free hydroxyl groups should remain in the final blends, whereas amino-building blocks are expected to be ultimately entirely incorporated inside the chains, leaving no remaining NH 2 groups. Ideally, the final amount of released OH groups should be equal to the initial amount of NH 2 groups. The reaction of amino groups with isocyanate being very fast, the polysiloxane blocks should be rapidly incorporated along the TPU chains through the formation of aryl-alkyl urea groups (Scheme 1). Finally, the introduction of low molar mass silicone blocks and the resulting chain rearrangements, together with possible chain scissions during extrusion [START_REF] Endres | The kinetics of the thermal decomposition of thermoplastic polyurethane elastomers under thermoplastic processing conditions[END_REF] should result in an overall decrease in the average molar mass of the polymer blends compared to the neat TPUs, and this decrease should be even sharper when higher amounts of low molar mass species are initially added. This hypothesis was indeed confirmed experimentally by size exclusion chromatography (see below, § 3.2)

The axial force of the extruder screw is related to the viscosity of the molten mixture, thus it gives qualitative indications about mixing and reaction progress (Figure 2a). The injection of liquid oligosiloxanes is done at time t=0. After injection, a lubrication phenomenon induced by the addition of low viscosity immiscible silicone to the high molar mass viscous polyurethane leads to a drastic drop of the screw force. Afterwards, the force rises strongly because of two competing phenomena, physical dispersion of the liquid silicone and reaction of the oligomer amino end-groups with reversible urethane bonds. A maximum force is reached once the aminoisocyanate reaction is completed. Samples of the blend were taken off after 3 minutes mixing and then 2, 12 and 17 minutes after reaching the maximum of the force. The extrudates appearance is very different depending on the time of sample withdrawal. As shown in Figure 2, with PDMS40D, the first sample was totally opaque while the second was perfectly transparent.

Then this transparency is progressively lost by increasing the mixing time (the third sample is translucent), ending with a fully opaque blend. Explanations for such a behavior will be given in the next section devoted to the morphology analysis of the blends.

As it is classically done for reactive processing in a batch mixer, the apparent rate of reaction between the silicone amino end-groups and the TPU can be roughly estimated from the evolution of the force recorded on the extruder screws. The strong increase in the force enables to set approximately the completion of the reaction and the end of blend homogenization. Using this criterion it was found that the reaction time for total conversion was directly related to the proportion of oligosiloxane and to a lesser extent to its viscosity (see supplementary data). As an example, when the weight fraction of oligosiloxane increases, the screw force maximum is shifted to longer times. This feature is obviously not related to differences in the actual kinetics of the amine-isocyanate reaction, but rather to the influence of the mixing process on the reaction. Indeed, the amine-isocyanate reaction has to occur at the TPU-silicone interface. At the early stage of the mixing process, the silicone is dispersed in the shape of large liquid drops in the TPU matrix. As soon as the reaction proceeds, a copolymer is formed and located at the interface that has to be continuously renewed for the reaction to proceed. This mechanism takes time, which explains that higher reaction times are observed when for the same PDMS oligomer the weight fraction of silicone increases (although the initial NH 2 concentration also increases).

In contrast no clear variation was observed when using equal weight fractions of oligomers with increasing molar mass (i.e. decreasing NH 2 concentration). The influence of the viscosity of the silicone is less clear and quite similar apparent reaction rates are observed depending on the molar mass of the used oligomers(see supplementary data).

Morphology and reactions

As a preliminary comment, it is interesting to notice that the liquid ,-aminopropyl oligosiloxane is highly immiscible and incompatible with the E685 polyurethane as can be seen in Figure 3 where the micrograph of a blend prepared in solution is presented. In this blend prepared at ambient temperature, no chemical reactions occurred between the polyurethane and the silicone, and the morphology is coarse with silicone domains of 2-10 m diameter, typical of an immiscible mixture where interfacial interactions are very weak.

When the blend is prepared at high temperature, the extrudate collected after only 3 minutes mixing is highly inhomogeneous. It contains unreacted liquid silicone badly dispersed in the melt. These domains appear as white holes in the polymer slice cut for microscopy (Figure 2A). They are however quite small (0.5-2 m) and significantly smaller than those of the blend prepared at ambient temperature. This demonstrates that even at this short reaction time some copolymer is already formed at the interface and is compatibilizing the dispersion. The material also contains very small dark domains visible on the TEM images (zoom in Figure 2A) that could be composed of reacted polysiloxane finally dispersed.

The second sample, taken after 12 minutes mixing, just after the maximum of screw force is shown in Figure 2B. Very small spherical silicone domains of about 20 nm diameter are then present in the TEM image, they appear throughout the sample and larger diameter structures are no longer visible. The distribution of the diameters is very narrow. This morphology is related to the total reaction of the ,-aminopropyl oligosiloxanes leading to their complete insertion into the polyurethane macromolecules as will be demonstrated later. This very fine morphology explains that the extrudate is perfectly transparent since such small heterogeneities do not diffuse visible light. However, this morphology is not stable as seen in Figure 2C and D. The nanodomains of silicone coarsen significantly upon longer mixing time (22 and 27 minutes).

This coarsening is not caused by physical coalescence of droplets since they are no longer composed of individualized silicon chains, but rather of PDMS segments from the PU-co-PDMS copolymer. The coarsening is thought to be linked to a change in the structure of the macromolecules and more precisely to an increase in the average length of the PDMS blocks. At this point, it is necessary to get a deeper insight into all the different chemical reactions that are likely to occur in the TPU/silicone reactive medium.

 Formation of aryl-alkyl urea

As discussed before, the amine and the urethane can react directly, or the urethane can split into alcohol and isocyanate. The latter is able to react very rapidly with the amine and thus form an aryl-alkyl urea during the early stages of the process [START_REF] Lu | Reactivity of common functional groups with urethanes: models for reactive compatibilization of thermoplastic polyurethane blends[END_REF]. To assert this comment, Figure 4 depicts the 1 H NMR spectrum of the blend of E685 with 10wt% PDMS40D after 12 minutes mixing. A new signal labeled "b" from CH 2 NH(C=O) protons clearly shows the formation of the aryl-alkyl urea while the signal labeled "a" from CH 2 OH protons reflects the formation of hydroxyl-terminated TPU. In parallel, the peak of the amine end-groups of oligosiloxane decreases to reach a non detectable value. The conversion of these amine groups is already complete after 12 minutes of reaction with 10 wt% PDMS40D and almost complete after 10 minutes with 10wt% PDMS15D (Table 2). These mixing times are in good agreement with those estimated from the maximum force of the screws observed at 10 and 7 minutes for PDMS40D and PDMS15D, respectively.

 Formation of alkyl-alkyl urea

The formation of TPU/PDMS aryl-alkyl urea is not the only chemical reaction that can be encountered. By reacting their diurethane model with an aliphatic primary amine, Lu et al. have

shown that after longer mixing times, alkyl-alkyl urea were formed and aryl-aryl urea were formed even later [START_REF] Lu | Reactivity of common functional groups with urethanes: models for reactive compatibilization of thermoplastic polyurethane blends[END_REF]. Those functions come from the cleavage of unsymmetrical urea groups into alkyl isocyanate and primary aromatic amine. The reaction of the former with primary alkyl amine borne by the oligosiloxane leads to alkyl-alkyl urea while reaction of aromatic amine with urethane functions gives aryl-aryl urea. The reactions are outlined in Scheme 2.

As evidenced by 1 H NMR, under our conditions, a low proportion of alkyl-alkyl urea is formed during the early stages of the process (signal labeled "c" from CH 2 NH(C=O) protons) (Figure 4) and the concentration increases as the reaction proceeds (Table 2).

In the material, PDMS segments can be ended with three different functions (amine, arylalkyl urea and alkyl-alkyl urea, see Scheme 2) which ratio depends on the mixing time. These ratios are given in The overall contribution of polysiloxane terminal groups does not always complete 100% with respect to the signal at 0.59 ppm, which means that other functions coming from side reactions of amines should coexist in the blend. One small signal is actually observed at 10.9

ppm and could suggest the formation of some allophanate functions through the addition of aromatic isocyanate on TPU urethane groups. The NMR analysis of such allophanate functions have already been reported: their NH NMR signals were observed at 10.85 ppm in acetone-d6 solution [START_REF] Heintz | Effects of reaction temperature on the formation of polyurethane prepolymer structures[END_REF] and 10.65 ppm in DMSO-d6 [START_REF] Lapprand | Reactivity of isocyanates with urethanes: conditions for allophanate formation[END_REF].

 Formation of aryl-aryl urea

Concerning aryl-aryl urea functions, their formation is difficult to confirm as their NH NMR peaks are expected to be very close to the NH urethane signals. Nevertheless three different couples of signals were observed in the region of aromatic protons for all the sampled blends. As the main signals at 7.34 ppm and 7.02 ppm are attributed to protons of aromatic ring bearing urethane functions, those slightly shifted to higher fields, at 7.28 ppm and 6.96ppm respectively, were assumed to correspond to aromatic protons close to aryl-alkyl urea groups. The smaller signals at 6.81 ppm and 6.46 ppm could then belong to aromatic rings of aryl-aryl urea functions.

These signals were present in all the spectra of sampled blends. They are shown in the 1 H NMR spectrum of the 10wt% PDMS15D blend after 25 minutes mixing (Figure 5).

The evolution of the proportion of amino groups from the oligosiloxane, aryl-alkyl groups and alkyl-alkyl groups of the TPU-co-PDMS copolymer are collected in Table 2 and the example of E685/PDMS15D 90/10 wt% blend processed at 200°C is depicted in Figure 6. To summarize briefly, the amino groups of the PDMS oligomers react rapidly with the isocyanates issued from the splitting of the urethanes and form aryl-alkyl ureas. Then alkyl-alkyl ureas are formed from the splitting of the aryl-alkyl groups.

These chemical observations now allow to go back over macromolecular structure, and its consequences on the morphology evolution displayed in Figure 2. The insertion of polysiloxane moieties into TPU chains is produced very early in the process through the reaction involving urethane groups and amino-terminated siloxane. Thus, regarding the large excess of urethane links compared to available amine, most of the siloxane units are initially individualized along the chain as depicted schematically in Scheme 3. As a consequence these units segregate in the shape of nanosized spherical domains (about 20 nm diameter, see Figure 2B). Later on, increasing amounts of alkyl-alkyl urea groups are formed. Each of them implies the presence of two consecutive polysiloxane segments in the copolymer chain, and thus the formation of alkyl-alkyl ureas leads to an increase in the average length of the PDMS blocks in this copolymer, as shown in Scheme 3. These blocks can then segregate in the shape of larger domains (about 100 nm or more depending on the composition of blend, see Figure 2C).

The NMR data collected and summarized in Table 2 have been used to evaluate an average number of consecutive polysiloxane segments, N. Indeed the proportion of aryl-alkyl and alkyl-alkyl ureas brings this information if we neglect the side reactions and consider the situation where there is no residual amine. N is expressed as:

2 1 alkyl alkyl aryl alyl N      (1)
Where  aryl-alkyl and  alkyl-alkyl are the molar ratio of respectively aryl-alkyl and alkyl-alkyl groups relative to the total (aryl-alkyl + alkyl-alkyl).

The average number of concomitent polysiloxane segments reaches values around 3-5 at the end of the process (Table 2).

Polyurethane/ thermoplastic polydimethylsiloxane-urea copolymer blends

Finally, a blend of the same polyurethane, E685, with a high molar mass thermoplastic silicone has been extruded to be compared with the previous E685/amino-terminated PDMS blends. The thermoplastic silicone is a polydimethylsiloxane-urea copolymer (Geniomer G80) obtained by the stoichiometric reaction of an ,-aminopropyl polydimethylsiloxane with isophorone diisocyanate (Figure 1). It was chosen because its structure is very similar to that hypothesized for the PDMS segments.

The morphology of the E685/Geniomer80 90/10 w% blend depicted in Figure 7a is showing spherical silicone domains with 50 to 250 nm diameter, which is actually very small for an immiscible polymer blend. By comparing Figure 7a with Figure 7b one can observe that the morphology of the E685/Geniomer 80 is quite similar to that obtained by blending during 24 minutes the large molar mass ,-aminopropyl oligosiloxane, PDMS130D, with E685.

Such a fine structure is linked to the reactions that possibly occur during the mixing process (Scheme 4). Under the used experimental conditions, the urethane and urea linkages from respectively polyurethane and PDMS-urea chains once again are reversible and prone to dissociation releasing isocyanate, alcohol and amine end-groups that can further react to produce block and segmented copolymers. Such copolymers are well known to compatibilize immiscible blends and stabilize their morphology. The compatibilization seems to be highly efficient since very small droplets of silicone are dispersed in the TPU matrix. Those blends of TPU and silicone thermoplastic elastomer G80 were investigated by NMR.

Neither 1 H nor 13 C NMR spectra of the TPU / G80 blend (20 wt% of G80) exhibit significant signals suggesting the formation of new urea or urethane groups; e.g. no new carbonyl signal from urethane groups could be detected in the 13 C spectra. From the 1 H spectrum it is difficult to show the formation of new alkyl-alkyl urea functions since they are already present in the spectrum of G80. Moreover the amount of reactive groups from Geniomer 80 that could react with the TPU is low. When 20 wt% of Geniomer 80 is used, the ratio of the total amount of isocyanate from TPU (MDI) to the amount of isocyanate from G80 (IPDI) is around 10:1. If only some part of the urethane and especially of the urea groups split during extrusion and give free reactive groups [START_REF] Endres | The kinetics of the thermal decomposition of thermoplastic polyurethane elastomers under thermoplastic processing conditions[END_REF], such a low content does not allow to confirm or to evaluate the extent of the reactions by means of NMR. An indirect proof of the formation of TPU-PDMS copolymer can nevertheless be found with the presence in the 1 H spectrum of the TPU/G80 blend of a multiplet at 3.5 ppm from hydroxyl end groups (-CH 2 -CH 2 -OH), this signal is analogous to the one observed for the TPU/,-aminopropyl oligosiloxane blends. These end groups can only result from the splitting of urethane bonds and the free isocyanate was therefore able to recombine with some amine function originating from G80.

Attempts were finally done to monitor the extent of the reaction by means of size exclusion chromatography, SEC, both with oligosiloxanes and with the thermoplastic silicone elastomer. Initially one aim of the experiment was, especially in the first case, to detect if any non-reacted amino blocks remained in the final blends. But for mixtures with small amounts of oligomers, and even though the blends were modified with p-tolylisocyanate, it was in fact very difficult to detect the peak of these amino-blocks, and therefore SEC finally could not provide reliable data for the precise determination of remaining species. However it gave other valuable information. Using the amino-terminated oligomers, the main peak of the TPU globally shifted towards lower molar masses, which was another proof that a reaction occurred (scissions of TPU chains by low molar mass silicone blocks and rearrangements leading to a decrease in the overall average molar mass). As shown in Figure 8, this was not observed with the thermoplastic silicone elastomer, for which only a slight increase in the average molar mass was noted at the very beginning (also observed with the oligomers), followed by a rather constant value of this molar mass. As the reaction proceeds, and given that the TPU and the elastomer G80 have initially comparable molar masses, the rearrangement reactions induce no noticeable change in these M w values. Such result proves that the decrease in molar mass observed with the oligomers cannot be attributed to some thermal degradation, but definitely results from rearrangement reactions. If thermal degradation occurred, it would also be observed for the blends with G80 that are produced exactly under the same conditions, and it would result in a decrease in the overall molar mass in this case as well.

Conclusion

Generally, PDMS/TPU copolymers are obtained via a reaction involving TPU precursors and functional PDMS oligomers and most often the reaction is done in solution at moderate temperature. Here ,-aminopropyl polydimethylsiloxane (PDMS) and thermoplastic polyurethane (TPU) have been reacted in bulk and the elevated temperature allows observing specific behaviour. At 200°C the splitting of the urethane group of the TPU is releasing hydroxyl and isocyanate terminated TPU chains. The splitting is followed by the fast reaction of the amino-siloxane with the isocyanate. This reaction forms aryl-alkyl ureas and promotes the insertion of individual PDMS segments inside the TPU chains, forming TPU-co-PDMS copolymer with a very fine morphology (~20nm). However the morphology is not stable and coarsening is observed for longer reaction times. This morphology evolution had to be explained.

Actually, after the formation of the aryl-alkyl ureas several reactions proceed. One is the cleavage of unsymmetrical aryl-alkyl urea groups into alkyl isocyanate and primary aromatic amine. The reaction of the former with primary alkyl amine borne by the oligosiloxane leads to alkyl-alkyl ureas that have been identified with the help of NMR. In terms of macromolecular structure, the formation of one alkyl-alkyl urea group implies that two polysiloxane segments are now consecutive in the copolymer chain. Thus, the formation of alkyl-alkyl ureas leads to an increase in the (average) number of successive polysiloxane segments that can reach value of 3-5. The longer PDMS blocks segregate in the shape of larger domains (~100 nm or more) depending on the composition of the material.

Scheme captions

H N O TPU O H N O TPU N C O H N O TPU HO TPU O Si O Si O Si O n NH 2 O Si O Si O Si O n H N H N O TPU H N O O H N O H N

Scheme 1 .Scheme 2 .Scheme 3 . 4 .

 1234 Scheme 1. Expected chemical reactions during the blending of E685 polyurethane with ,-

Figure 1 .Figure 2 .Figure 3 .

 123 Figure 1. Structure of the used polymers. a) Thermoplastic polyurethane elastomer E685, b)

Figure 4 . 1 HFigure 5 . 1 H

 4151 Figure 4. 1 H NMR spectrum of E685/PDMS40D blend 90/10 wt% after 12 minutes mixing time

Figure 6 .Figure 7 .Figure 8 .

 678 Figure 6. E685/PDMS15D 90/10 wt% blend processed at 200°C. Evolution of primary amine

Scheme 2 .Scheme 4 .Figure 1 .Figure 2 .Figure 3 .

 24123 Scheme 2. Reactions of urethane function and amino group of the silicone depicting aryl-alkyl,alkyl-alkyl and aryl-aryl ureas. The scheme should be read from left to right with the plus/minus symbols referring to the direction of the neighboring arrow.

Figure 4 . 1 HFigure 5 . 1 H

 4151 Figure 4. 1 H NMR spectrum of E685/PDMS40D blend 90/10 wt% after 12 minutes mixing time (THF-d8, 25°C). (a) CH 2 OH from hydroxyl-terminated TPU, (b) CH 2 NH(C=O) from aryl-alkyl urea, (c) CH 2 NH(C=O) from alkyl-alkyl urea

Figure 6 .Figure 7 .Figure 8 .

 678 Figure 6. E685/PDMS15D 90/10 wt% blend processed at 200°C. Evolution of primary amine

  Table 2 for several blends as a function of mixing time. They were calculated by dividing the integral values of signals at 2.60 ppm (CH 2 NH 2 of residual amine), 3.15 ppm (CH 2 NH(C=O) of aryl-alkyl urea ) and 3.06 ppm CH 2 NH(C=O) of alkyl-alkyl urea) by the integral value of signal at 0.59 ppm (CH 2 groups linked to PDMS).

Table 1 .

 1 Number average molar mass of the polymers determined by SEC in tetrahydrofuran.

		Elastollan 685A	PDMS15D	PDMS40D	PDMS130D	Geniomer 80
	n M (g.mol -1 )	136000	1240	3040	8520	176000

Table 2 .

 2 Chemical groups ratio in E685/PDMS blends after different mixing times ( 1 H NMR determination).

			CH 2 NH 2 of	CH 2 NH(C=O) of	CH 2 NH(C=O) of	Total	Estimated
		Mixing time	residual amine	aryl-alkyl urea	alkyl-alkyl urea	(mol%)	number of
		(min)	(mol%)	(mol%)	(mol%)		adjacent
							polysiloxane
							segments, N
			2.60 ppm	3.15 ppm	3.06 ppm
	PDMS40D 10 wt%	12	0	64	28	92	1.9
		22	0	37	43	80	3.3
		27	0	33	60	93	4.6
	20 wt%	27	0	47	35	82	2.5
		37		37	56	93	4.0
		42	0	30	51	81	4.4
	30 wt%	21	33	36	17	86	-
		36	0	43	59	102	3.7
	PDMS15D 10 wt%	10	16	70	20	106	-
		15	0	68	31	99	1.9
		18	0	60	37	97	2.2
		25	0	54	53	107	3.0

Acknowledgements

Wacker Chemie is gratefully acknowledged for having supported this work and special thanks are addressed to Dr. Jürgen Pfeiffer.


M Kamali-Nejad

T Azib

Member, IEEE G Remy

Member, IEEE O Bethoux

Member, IEEE C Marchand

Optimization of a Power Electronic Structure for Hybrid Fuel Cell/Ultracapacitors Vehicle

This paper deals with the optimization of a parallel hybrid Fuel Cell (FC)/ Ultracapacitors (UCs) power source for automotive applications. The aim of this hybridization and its control are to fulfil the load requirements as well as to comply with the component constraints (high efficiency, reduced weight and cost, etc.). First, a classical hybrid architecture using FC/UCs and a two-converter structure is presented. After that, power requirements for automotive applications are analyzed, and the load power demands for the FC/UCs sources are deduced. Secondly, the model-based design approach is used for the optimization, and a selection of the main components to be optimized is presented. Thus, the simulation model and the control strategy are detailed. This model has been validated experimentally. Finally, an optimization algorithm is designed using Parallel Computing and the Genetic Algorithm toolbox of Matlab/Simulink. The retained criterion is based on the reduction of the total volume of the system. I.

INTRODUCTION

Nowadays, with more and more severe emission standards and the reduction of global fuel supplies, electric/ hybrid vehicles are a very viable technology. It offers a significant potential to meet both fuel economy and emission cleaning demands, particularly the Fuel Cell Hybrid Electric vehicle FCHEV [START_REF] Phuc | Control Strategy for Hybrid Electric Vehicles Based on Driver Vehicle Following Model[END_REF]. This system combines a fuel cell power system with an energy storage system (ESS) as a power source. The ESS (e.g., batteries or ultracapacitors) helps the fuel cell to supply peak power to the electric drive system during acceleration and recovers energy during deceleration, which consequently improves the efficiency and dynamic characteristics of the power system [START_REF] An | performance analysis according to the combination of energy storage system for fuel cell hybrid vehicle[END_REF]. On the other hand, the ultracapacitor FCHEV has three major advantages over a battery FCHEV, satisfying almost the same packaging and cost requirements: higher energy efficiency, better accelerating performance and extended fuel cell life [START_REF] An | performance analysis according to the combination of energy storage system for fuel cell hybrid vehicle[END_REF], [START_REF] Barrade | The use of supercapacitors for energy storage in traction systems[END_REF]. However, the energy storage system increases the complexity of the overall power plant, the control design having an important effect on system performance.

As well as the control design, it is very important for component sizing to be taken into account in order to reduce installation investment costs and to achieve good overall performance, which is the main objective of this work.

Firstly, in this paper, a classical hybrid architecture using FC/UCs and a two-converter structure is presented. Then, power requirements for automotive applications are analyzed using a frequency decomposition of the load power demands for the FC/UCs sources. Then, a selection of the main components to be optimized is presented using an energy repartition calculation. As the optimization used a modelbased design approach, the simulation model and the control strategy are detailed. This model has been validated experimentally using a Nexa Ballard PEMFC as main source and BOOSTCAP Maxwell ultracapacitors as auxiliary source. Finally, an optimization algorithm is introduced and configured for Parallel Computing using the Genetic Algorithm toolbox of Matlab/Simulink. The retained criterion is based on the reduction of the total system volume.

II. HYBRID VEHICLE STRUCTURE

For the powertrain of hybrid vehicles using fuel cells as main power source, different hybrid fuel cell system configurations have been suggested and investigated (series, parallel, and cascade) [START_REF] Hajizadeh | Intelligent power management strategy of hybrid distributed generation system[END_REF], [START_REF] Cacciato | A critical evaluation and design of bi-directional DC/DC converters for supercapacitors interfacing in fuel cell applications[END_REF], [START_REF] Jiang | Design and,experimental tests of control strategies for active hybrid fuel cell/battery power sources[END_REF]. It has also been proven that the parallel structure is most advantageous [START_REF] Cacciato | A critical evaluation and design of bi-directional DC/DC converters for supercapacitors interfacing in fuel cell applications[END_REF], [START_REF] Jiang | Design and,experimental tests of control strategies for active hybrid fuel cell/battery power sources[END_REF]: fewer component constraints, easy energy management and good reliability.

A. Two-converter parallel structure Among the major topologies of parallel structures for FC/SCs hybrid power system, the two-converter parallel structure (Fig. 1) is the most used [START_REF] Jiang | Design and,experimental tests of control strategies for active hybrid fuel cell/battery power sources[END_REF]- [START_REF] Azib | Supercapacitors for Power Assistance in Hybrid Power Source with Fuel Cell[END_REF]. Some of the advantages are: no pollution, no noise, constant DC bus, simple control, improved FC lifetime, reversible energy conversion.

B. Power load reference profile

The needed car power is mainly due to speed variations, tyre friction dissipation, aerodynamics dissipation and mass elevation. This power, Pmotor, can be expressed as:

2 1 cos( ) sin( ) 2 motor r x dV P V C Mg Mg M SC V dt α α ρ   = + + +     (1)
Where:

V and M are the vehicle speed and mass (in m.s -1 and kg); α is the road angle with an horizontal line (in rad); Cr and Cx are the friction and aerodynamic coefficients; ρ is the air density (in kg.m -3 ) and g is the g-force (in m.s -2 ) S is the front surface area (in m 2 ). M=1000; g=9.81; α=0; Cr=0.01; Cx=0.30; ρ=1.225; S=2.5

These values correspond to a mid-size car adapted to urban transportation. Thus, the vehicle power demand can be determined by the driver's requirements. Indeed, European light duty vehicles have to face the New European Driving Cycle (NEDC) which represents the typical usage of a car in Europe. The NEDC consists of repeated urban cycles (called ECE-15 driving cycle) and an Extra-Urban driving cycle, or EUDC. Fig. 2 shows the ECE-15 cycle with the speed and the power demand of a car following a flat road.

Sudden power changes can be noticed each time the driver requires a speed change. In this example, the car's average power is only of 0.72 kW, whereas the peak power reaches roughly 10 kW, which means a 13.7 (Pmax/Paverage) ratio. Finally, the power source optimization of the hybrid sources has to fulfill both requirements. These requirements are included in a model for optimization purpose.

III. MODEL DESIGN FOR OPTIMIZATION

Nowadays, the Model-Based Design approach is widely used to design systems from simulation to real-time application [START_REF] Nicolescu | Model-Based Design for Embedded Systems[END_REF]. The optimization strategy uses this Model-Based Design approach to accelerate the design of hybrid vehicle regarding various requirements [START_REF] Marco | The Use of Model Based Design Techniques in the Design of Hybrid Electric Vehicles[END_REF].

In this paper, we focus on the design of a FC/UC hybrid vehicle. The main principle of such a design approach is to optimize the parameter values of a model that has been previously validated using experimental results: to be more precise, a Matlab/Simulink model of the hybrid vehicle is designed including the control strategy. The model parameters are identified using a test bench of about 1kW, and simulation results are compared to the experimental results. Then, the 1kW Model is extended to optimize the system parameters of a 10kW car's power regarding the NEDC requirements. Of course, the parameter validity domain has to be carefully considered. Nevertheless, such method allows for greatly reducing the size of the test bench.

In the following section, the model and the control used are described. Then, some experimental validations of the presented system model are depicted.

A. Model description

The model-based design approach is used to facilitate the global optimization of the system. Indeed, the simulation includes the component models and their associated controls, and it allows to evaluate the system performances for a requested load on each component and to estimate the system's stability.

Ultracapacitors model

Various ultracapacitor models can be found in literature, especially for hybrid systems [START_REF] Ehsani | Modern Electric, Hybrid Electric, and Fuel Cell Vehicles Fundamentals, Theory, and Design[END_REF], [START_REF] Bauman | A comparative study of fuel cell-battery, fuel cell-ultracapacitor, and fuel cell-battery-ultracapacitor vehicles[END_REF]. Classically, an ultracapacitor theoretical model uses a transmission line with the voltage (VSC) dependent on distributed capacitance [START_REF] Rafik | Contribution to the sizing of supercapacitors and their applications[END_REF].

However, to take into account the supercapacitor behavior during charge and discharge and to preserve a sufficient accuracy, a RC model is utilized to describe the nonlinear electrode behavior:

-The resistor R models the supercapacitor ohmic loss, usually called equivalent series resistor (ESR). -The capacitor C represents the supercapacitor capacitance during charging and discharging effects. The RC parameter values are extracted directly from the manufacturer's datasheets.

Fuel Cell model

In low dynamic conditions, the fuel cell characteristics can be considered as a voltage source with ohmic, kinetic and mass transfer resistances [START_REF] Barbir | PEM Fuel Cells-Theory and Practice[END_REF]. Consequently, the fuel cell can be connected either in series with a current source, which is the case in the two-converter structure (L in series with the FC) or in parallel to the DC bus, which is a perfect voltage source with its associated capacitor CBUS. Of course, the capacitor value has to remain low in order to prevent difficult starting procedures.

Furthermore, the fuel cell system is a complex device with many auxiliary components. Hence, a significant part of the electrical power generated is used internally, which means that the real stack current is greater than iFC. More precisely, the relationship between the FC voltage vFC and the output current iFC is given by the following equations:

( ) ( ) ( ) ln exp FC Cell Stack Stack l Stack V N E R j A j j m nj = -⋅ -⋅ + - (2) ( ) ( ) 2
, and here: 1

Stack Stack Stack FC FC Cell I j I I I A α β γ = = + + + ( 3
)
Where ACell is the area of each cell, N is the stack cell number, Ecell is the reversible cell voltage, R is the membrane area specific resistance; A is the Tafel coefficient; and m and n are the two coefficients of the mass transfer equation [START_REF] He | Diagnostic tool to detect electrode flooding in proton exchange membrane fuel cells[END_REF]. Fig. 4 shows the experimental characteristics and the simulated model of the studied FC.

The experimental set-up of the hybrid Fuel Cell / Ultracapacitors power source uses a Nexa BALLARD PEM fuel cell stack of 1200 W as the main source, as shown in Figure 5. The auxiliary source is obtained with a Maxwell SC module associating two modules in series: Maxwell BMOD0250 (each module is achieved with a series connection of six individual "2.7V, 1500F" elements).

A ZS electronic load ZS1806 of 1800 W (i= 0….150A / V=0…..60V) fulfills the power load demands. The power converters DC/DC are realized with standard IGBT modules (SEMITRANS: SKM50GB123D). A real-time dSPACE DS1104 controller board is used to implement the energy management control strategy designed with Matlab/ Simulink.

B. Control description

The two-converter structure consists in associating a static converter and a control loop to each power source (Fig. 1). This control strategy generates a large number of degrees of freedom in the control design (Fig. 6). This classical hybrid structure has been validated in various works [START_REF] Thounthong | Control strategy of fuel cell/supercapacitors hybrid power sources for electric vehicle[END_REF]- [START_REF] Azib | Supercapacitors for Power Assistance in Hybrid Power Source with Fuel Cell[END_REF], [START_REF] Azib | An Innovative Control Strategy of a Single Converter for Hybrid Fuel Cell/Supercapacitors Power Source[END_REF].

The control strategy is based on a power decoupling strategy in the frequency domain of the power source. This energy-management strategy fulfills the fast energy demands of the load and respects the integrity of each source. The full details of the control strategy have been published in previous work [START_REF] Azib | An Innovative Control Strategy of a Single Converter for Hybrid Fuel Cell/Supercapacitors Power Source[END_REF].

C. Optimization algorithm

In literature, optimization of a mixed continuous and discrete model is known as hard optimization [START_REF] Floudas | Encyclopedia of Optimization[END_REF]. The four most widely known metaheuristics for hard optimization are [START_REF] Glover | Handbook of Metaheuristics[END_REF]:

-The simulated annealing method; -The tabu search; -The ant colony algorithms; -The genetic and evolutionary algorithms. Some common features clearly appear in most metaheuristics, such as the use of diversification, to force the exploration of regions of the search space, rarely visited until now, and the use of intensification, to go thoroughly into some promising regions. One common drawback for most metaheuristics still is the delicate tuning of numerous parameters; theoretical results available by now are not sufficient to really help, in practice, the user facing a new hard optimization problem.

Among them, we have chosen the Genetic Algorithms available in MATLAB® products. Because it can be used for such problems that are difficult to solve with traditional optimization techniques, including problems that are not well defined or are difficult to model mathematically. Precisely, GA can be used when computation of the objective function is discontinuous, highly nonlinear, stochastic, or has unreliable or undefined derivatives. Some of the GA advantages over the traditional numerical optimization approaches are related to:

-The optimization with continuous or discrete parameters; -The derivative information, which is not requested; -The simultaneously searches from a wide sampling of the cost surface; -The well suited implementation for parallel computers.

Nevertheless, in the field of continuous optimization, there is a significant arsenal of traditional methods used for global optimization [START_REF] Horst | Handbook of Global Optimization[END_REF], but these techniques are often ineffective if the objective function does not possess a particular structural property, such as convexity. In order to overstep this convergence issue and to reduce the GA calculation time, a Monte-Carlo analysis can be performed on the system. Indeed, it allows GA to start with initial parameter values leading close to the solution.

Besides, the definition of the optimization problem is often supplemented by the constraints. All the parameters of the adopted solutions must satisfy these constraints, or otherwise these solutions are not realizable. But including constraints in the optimization induces a larger time calculation for the GA [START_REF] Conn | A Globally Convergent Augmented Lagrangian Algorithm for Optimization with General Constraints and Simple Bounds[END_REF]- [START_REF] Miettinen | Numerical comparison of some penalty-based constraint handling techniques in genetic algorithms[END_REF]. Then the problem formulation can be reconsidered in order to avoid formulations with constraints [START_REF] Baldick | Applied Optimization, Formulation and Algorithms for Engineering Systems[END_REF]. Indeed, some transformations are critical in matching a problem to the characteristics of an algorithm.

In the studied case, constraints can be avoided using penalties method on the model [START_REF] Luenberger | Linear and Nonlinear Programming[END_REF]. However, even if stochastic methods such as genetic algorithms can be easily implemented, algorithm parameter tuning is highly important to achieve the global optimum.

The advantage of this kind of algorithm is that, by creating a large enough population, the probability of finding the global minimum point is higher. The drawback is the time needed to run such an algorithm [START_REF] Pardalos | Handbook of global optimization[END_REF]. Therefore the parameters should be well tuned. In our case, two main strategies are chosen in order to tune the GA.

Firstly, we performed a Monte Carlo (MC) simulation in order to put forward an initial population to GA. This fact will help GA in order to converge to the global optimum. Performing a MC simulation with a large population and selecting the best points is a better way for creating an initial population for GA compared with the random creation of an initial population. To perform the MC simulation, a random value is generated for each variable by using a classical generator as an extended cellular automata generator [START_REF] Kier | Cellular Automata[END_REF].

When tackling the multi-criteria optimization of a device in electrical engineering, the exhaustive sampling of the Pareto optimal front implies the use of complex and time-consuming algorithms that are unpractical from the industrial viewpoint. In several cases, however, the accurate identification of a few non-dominated solutions is often sufficient for design purposes. An evolutionary methodology of the lowest order, dealing with a small number of individuals, is proposed to obtain a cost-effective approximation of non-dominated solutions. The optimal design of a hybrid power source, based on the optimization of cost, volume, weight and performance of the device, is obtained in Table I.

IV. OPTIMIZATION RESULTS

In our case, after performing the Monte Carlo simulation with a 5000 population, 24 best points are selected among the results and they are put into GA as part of the initial population. The rest of the initial population (40 points) is created by GA using a uniform distribution. This allows GA to sweep a wider search zone.

Secondly, the other parameters of GA such as mutation rate, elite count, selection function, etc., are chosen so as to provide a maximum chance for GA to achieve the global optimum.

V. CONCLUSION

In this paper, an optimization of a parallel hybrid Fuel Cell (FC)/ Ultracapacitors (UCs) power sources for automotive applications has been described. A Model-Based Design approach has been used to accelerate the design and to reduce the size of the test bench. An optimization using a paralleled Genetic Algorithm has been performed. The aim is to find the optimized value of the system parameters within the lowest possible volume in order to fulfill a NEC-15 for a mid-size car. The final paper will include the results of the multiobjectives regarding the volume, the cost, the weight and the energy consumption of the car.

Fig. 1 .FC

 1 Fig. 1. Two-converter parallel structure for FC/UCs hybrid power system.

Fig. 3 .

 3 Fig. 3. Model and Control Scheme using Simulink.

Fig. 2 .

 2 Fig. 2. ECE-15 Driving Cycle.

Fig. 5 .

 5 Fig. 5. The components of the experimental bench.

Fig. 4 .

 4 Fig. 4. PV Characteristics.

Fig. 6 .

 6 Fig. 6. Optimization results of the Weight/Volume objective functions.

Fig. 8 .

 8 Fig. 8. Optimization results of the cost function of the number of generation.

Fig. 7 .

 7 Fig. 7. Optimization results of the Weight/Cost objective functions.

  Accepted article for publication in 2010 IEEE Vehicle Power and Propulsion Conference, Lille, France -Citation information: doi: 10.1109/VPPC.2010.5729165 Available online 10 March 2011The optimized values of the system parameters are depicted in TableII, and the obtained results in TableIII:

	22.8 22.8
	22.6 22.6
	22.4 22.4
	18.1 18.1	18.2 18.2	18.3 18.3


M Ghanes

O Bethoux

M Hilairet

J-P Barbot

Fuel Cell System Control under Converter Losses with Experimental Results

Keywords: Fuel cell, ultracapacitors, power management, energy-shaping, singular perturbation approach, converter losses

The energy management issue of a hydrogen fuel cell (FC) with ultracapacitors (UCs) for applications with high instantaneous dynamic power is considered. Singular perturbation approach is employed for the control of two converters and converter losses. The resistance load is considered unknown variable and the proposed controller uses load current and voltage measurements. The convergence of the voltage controller under converter losses is analysed by using Lyapunov theory. According to a significant FC-UCs benchmark, experimental results highlight the interest of the proposed approach and also show the difficulty due to converter losses.

PROBLEM STATEMENT

In this paper, we consider the power management issue of hydrogen Fuel Cell (FC) system associated to a reversible impulse energy source (the ultracapacitors) due to the development increasing of electric and hybrid vehicles since 2009. The FC must deliver a slowly varying current, not more than 4A/s for a 0.5kW/12.5V FC [START_REF] Thounthong | Energy management of fuel cell/battery/supercapacitor hybrid power source for vehicle applications[END_REF], and 10A/s for a 20kW/48V FC [START_REF] Corbo | Pefc stacks as power sources for hybrid propulsion systems[END_REF] as examples. That is why the FC needs to be associated with auxiliary sources. Nowadays, these auxiliary sources can either be batteries or supercapacitors (UCs). Sometimes, batteries are not able to bear high power charge and discharge conditions, whereas UCS have a high power range. Therefore, for fast power demands, supercapacitors are probably the best-suited components [START_REF] Rodatz | Optimal power management of an experimental fuel cell/supercapacitor-powered hybrid vehicle[END_REF]. This supply set leads to a design challenge requiring to conceive an adapted architecture, to choose power components, and to dene the appropriate associated control strategy [START_REF] Thounthong | Supercapacitors as an energy storage for fuel cell automotive hybrid electrical system[END_REF]. The parallel structures with only one converter [START_REF] Davat | th International Symposium on Advanced Electromechanical Motion Systems[END_REF] or two converters [START_REF] Thounthong | Supercapacitors as an energy storage for fuel cell automotive hybrid electrical system[END_REF] can fully respect the mentioned requirements. This paper is dedicated to the study of the structure with two converters as shown in Fig. 1. This architecture offers lower stress for components and an easy and a more reliable power management [START_REF] Cacciato | A critical evaluation and design of bi-directional dc/dc converters for supercapacitors interfacing in fuel cell applications[END_REF]. However it inserts extra weight and complexity (consequently a control challenge). The three major objectives of this power architecture management are to match globally the power (positive or negative) demanded by the load while satisfying FC dynamics (mainly limited by the time response of the air compressor) and monitoring the storage device (UCs) state of charge. It means that the current delivered by the FC must have smooth behaviour in order to ensure its live time, while the UCs provide the load power transient. The hybridization of the electrical supply is highly recommended to limit the dynamic effects and a special attention is also paid for the shut-down and start-up procedures to preserve the fuel cell performances and enhance its durability. Therefore, it seems clear that the DC bus regulation is managed by the ultracapacitors. (2010) or empirical control associated with a multi objective genetic algorithm optimization [START_REF] Paladini | Super-capacitors fuel-cell hybrid electric vehicle optimization and control strategy development[END_REF] that has been applied for the supervisory power train control problem in charge sustaining hybrid electric vehicles. However, these approaches are based on the a priori knowledge of the load, thus real-time control is not straightforward implementable. The frequency decoupling method using two cascaded loops of currents and voltages proposed in [START_REF] Davat | th International Symposium on Advanced Electromechanical Motion Systems[END_REF] allows to have fast transient load power provided by the UCs supply, while the FC supplies the mean power (This slow FC dynamics increases the system life time).

The control is mainly based on the change in the DC bus voltage induced by load power variations. The controller gains are tuned to ensure the closed-loop system stability, although it has not been theoretically proved.

In [START_REF] Talj | A controller tuning methodology for the air supply system of a pem fuel-cell system with guaranteed stability properties[END_REF] work, a passivity based approach is adopted experimented which enables to prove the closedloop system stability. However this control is sensitive to the load knowledge. To overcome this drawback, the authors have proposed to use either an observer estimating the unknown resistance load or an integral action.

In [START_REF] Ghanes | Singular perturbation control for coordination of converters in a fuel cell system[END_REF] a controller based on singular perturbation approach [START_REF] Khalil | Nonlinear systems[END_REF], [START_REF] Kokotovic | Singular Perturbation Methods in control: Analysis and design[END_REF]) is proposed to respect the slower FC dynamics, control of the storage device (UCs) state of charge and the power response required by the load with a DC bus regulation. This solution permits to avoid the use of a PI controller loop which is not very robust to load variation [START_REF] Davat | th International Symposium on Advanced Electromechanical Motion Systems[END_REF] and the use of an observer or an integral action [START_REF] Talj | A controller tuning methodology for the air supply system of a pem fuel-cell system with guaranteed stability properties[END_REF]. Moreover this approach is well adapted to the problem of FC-UCs control where the FC and UCs currents must be slow and fast respectively. A stability proof of the system is given. Nevertheless the converter losses are not taken account in [START_REF] Ghanes | Singular perturbation control for coordination of converters in a fuel cell system[END_REF].

In this paper, the losses of the two converters used with FC-UCs system are considered. The contribution is to analyse the problem of FC-UCs control in presence of converter losses by the proposed singular perturbation approach. Moreover, experimental results are given to highlight the validity limit of the proposed control against converter losses.

The paper is organized as follows : In section 2, the idea of the singular perturbation (SP) approach is recalled. The FC-UCs model is presented in section 3. In section 4, the proposed controller with a stability proof under converter losses is develloped. Section 5 reports the experimental results. Conclusion and on-going work are given by section 6.

BRIEF RECALLS ON SINGULAR PERTURBATION

Let us consider the following nonlinear system:

ẋ = f (x, z, ε) (1)
ε ż = g(x, z, ε) (2) with x ∈ R m , z ∈ R n ,
ε a small positive parameter and f , g two analytical vector fields of appropriate dimension. Roughly speaking, x can be see as the slow state and z as the fast variable. Nevertheless, this statement must be clarified and some assumptions and theoretical developments must be added. First of all, it is usual to decompose the system (1)-( 2) in decoupled two time scales dynamic. For this purpose, it is important to be able to compute the so-called slow manifold z = φ(x, ε), this manifold is the z behaviour when the fast transient time is finished ("outside the boundary-layer"). The slow manifold φ must verify the following equation: [START_REF] Vasil'eva | Asymptotic behavior of solutions to certain problems involving nonlinear differential equations containing a small parameter multipling the highest derivatives[END_REF]. For example the so-called frozen solution verify: 0 = g(x, α 0 (x), 0) Condition is requested for the α 0 existence:

ε φ(x, ε) = g(x, φ(x, ε), ε) (3) where φ(x, ε) = ∞ i=0 α i (x) ε i i! is computed iteratively
Assumption A.1: The Jacobian { ∂g(x,z,0) ∂z } is regular in the considered state space x ∈ D x ∈ R M and z ∈ D z ∈ R n .
This assumption is directly linked to the implicit function theorem and in nonlinear case more than one solution is possible (this particular case is outside the scoop of this short presentation). Now, it is important to know if the system (1)-( 2) converges on a slow manifold this is given by the well known Tikhonov's theorem [START_REF] Tikhonov | Ordinary differential equations[END_REF] but before recall the theorem it is necessary to analyse the fast dynamic on the boundary-layer. For this, a new state variable η = zφ is introduced and η converges rapidly to zero if the system behavior converges on the slow manifold. The η dynamics is equal to

η = 1 ε g(x, φ(x, ε) + η, ε) - ∂(φ(x, ε), ε) ∂t (4)
Setting ς = t ε , (4) may be rewritten

∂η ∂ς = g(x, φ(x, ε) + η, ε) -ε ∂φ(x, ε), ε) ∂ς (5)
Assumption A.2 The system (5) is locally in η and uniformly in x exponentially stable.

Hereafter, the Tikhonov's Theorem without considerations of time domain and existence and uniqueness of the solution (for example Lipschitz conditions are implicit). Theorem 1. Under the assumptions A.1 and A.2 and for ε ∈ R + sufficiently small, after a transient t 1 the dynamics (1)-( 2) evolves on a slow dynamic on his dynamic is equal to:

ẋ = g(x, φ(x), ε) (6)
In many applications, ( 6) is approximated at first order (the frozen solution of φ) ẋ = g(x, α 0 (x), 0) (7) As in the considered application the dynamics are too dynamically closed, some fast actuators (high-gain feedback) are used [START_REF] Marino | High-gain feedback in non-linear control systems[END_REF]. Considering, for the seek of simplicity, the following dynamical system

χ = f (χ, ζ) (8) ζ = g(χ, ζ) + β(χ)u (9) with χ ∈ R m , ζ ∈ R n , u ∈ R n and β regular for all χ. Then setting, for example u = -1 ε β(χ) -1 (ζ -α 0 (χ)) the dynamics become: χ = f (χ, ζ) (10) ε ζ = εg(χ, ζ) -(ζ -α 0 (χ)) (11)
Dynamics ( 10)-( 11) are similar to the dynamics (1)-( 2), thus, it is possible to use the theorem 1 and the slow dynamic of ( 10)-( 11) in first approximation is equal to

χ = f (χ, α O (χ)) (12)
Remark 1. In this paper, we often restrict our purpose to the first approximation of the slow manifold i.e. φ α O (x). Nevertheless, for example, when the system behaviour is too closed to a fold it is necessary to do a better φ's approximation.

FUEL CELL-ULTRACAPACITORS MODEL: TWO TIME SCALES MODEL

The "fuel cell-ultracapacitors" (FC-UCs) system (figure 2) is modelled by the following 5 th order non-linear state space model :

ẋ1 = (1 -α 1 ) x 4 + (1 -α 2 ) x 5 -x 3 -i d C ẋ2 = - x 5 C sc ẋ3 = -R l x 3 + x 1 L l ẋ4 = -(1 -α 1 ) x 1 + z L f c ẋ5 = -(1 -α 2 ) x 1 + x 2 L sc (13)
with state space

x(t) = [x 1 , x 2 , x 3 , x 4 , x 5 ] t = [v b , v sc , i l , i f c , i sc ] t , control inputs u(t) = [u 1 ; u 2 ] t = [1 -α 1 ; 1 -α 2 ] t and measures y(t) = x, z(t) = v f c .
Moreover the actual power converter has power losses which has to be taken into account using the term i d . Current i d refers to the to the losses current of FC converter (i d1 ) and UC converter (i d2 ). These current depend of the main currents (x 4 and x 5 ) and are expressed as the sum of a linear term and a quadratic term

i d = i d1 + i d2 = (V 0 + R 0 x 4 )x 4 x 2 (14)
+ (V 0 + R 0 x 5 )x 5 z
where V 0 , R 0 represent respectively the threshold voltage and the internal resistance of switches.

CONTROL DESIGN

Two fast actuators are applied to the FC-UCs model (13). These fast actuators are implemented with two PI controllers. The dynamics of both PI are chosen such that assumptions A.1-A.2 are verified and hence outside the boundary-layer the dynamic is on the slow manifold and the system behaviour is given by equation ( 7).

DC/DC

v b i sc L sc u 2 (t) T 1 T 2 v sc C sc DC/DC v f c i f c L f c u 1 (t) C R l L l i sc i l
Fig. 2. FC-UCs system.

According to (7), the fast actuators allows to get a slow reduced dynamics of ( 13) in input-state representation where i * f c is considered as the input u 1 and i * sc as the input

u 2 : ẋ1 = 1 C ( z 3 x 1 u 1 + x 2 x 1 u 2 -x 3 -i d ) ẋ2 = - u 2 C sc ẋ3 = -R l x 3 + x 1 L l (15) with x = [x 1 ; x 2 ; x 3 ] t = [v b ; v sc ; i l ] t ; control inputs u = [u 1 ; u 2 ] t = [i * f c ; i * sc ] t , measures y = [v b ; v sc ; i l ] t and z = [i f c ; i sc , v f c ] t .
The current losses i d given by ( 14) is rewritten as follows

i d = (V 0 + R 0 * u 1 )u 1 + (V 0 + R 0 u 2 )u 2 x 1 . (16)
It is important to note that the model ( 15) is valid only if all its outer closed loops dynamics are slower than the dynamics of the fast actuators, namely the PI controllers which assign the current i f c and i sc . This will be a constrain in the design of the slow controllers.

Slow control design

The desired equilibrium point are the following one

x * = [v * b ; v * sc ; v * b R l ],
with v * b and v * sc the DC bus and UCs desired voltages. As the constraint on i f c is such that its dynamic must be very slow (i.e. di f c dt < 10As -1 ), u 1 is considered as a slow input and has only a slow effect on the convergence of the dynamic (15). Consequently the gain on the u 1 loop must be very small and the variation of the component of u 1 must be filtered directly or implicitly. From these considerations u 1 is designed as

u 1 = x 1 z 3 i lm - C cs T slow e sc ( 17
)
where i lm is the load current average (in reality a load's current filtered with a low past filter). Moreover, e sc = v scv * sc and T slow are chosen in order that the feedback on e sc will be also sufficiently slow. It is important to note that e sc is not filtered because v * sc has a constant profil and v sc is proportional to the i sc integral and then implicitly filtered.

Fast control design

As u 1 is a slow input the second input u 2 can be faster than u 1 but slower that the current loops driven by both PI fast actuators, so u 2 is:

u 2 = x 1 x 2 [- C T f ast e b - z 3 x 1 ( x 1 z 3 i lm - C cs T slow e sc )+x 3 +i d ] (18)
with e b = v bv * b and T f ast is chosen in order that the feedback on e b remains slower that both PI fast actuators but faster than the loop on e sc . From ( 17) and ( 18), it is possible to set the following proposition. Proposition 1. Under the controls ( 17) and ( 18), the system (15) including converters losses (i d ) is locally exponentially stable for v * b and v * sc constants with a bias on x 2 = v sc .

Proof.

Considering the following Lyapunov candidate function

V = V 1 + V 2 + V 3 with V 1 = 0.5e i , V 2 = 0.5e 2 b and V 3 = 0.5e 2 sc
, where e i = (x 3i lm ) 2 and i lm = v * b /R l . By applying the controls ( 17)-( 18) to the reduced model ( 15), the differentiation of V gives:

V = e i -R l e i + e b L l - e 2 b T f ast -e sc ( x1 x2 [ C T f ast e b -z3 x1 ( x1 z3 i lm -Ccs T slow e sc ) + x 3 -i d ] C cs )
Choosing T f ast << 1, it follows that e b converges to zero independently of e i . Consequently e i also converges exponentially to zero and this uniformly with respect to e sc . Then

V = - R l L l e 2 i - 1 T f ast e 2 b - z 3 T slow x 2 e 2 sc + x 1 i d x 2 C cs e sc
Which is locally exponentially stable for x 2 > 0 with a bias convergence of e sc to x1i d x2Ccs depending on i d . Moreover T slow must be chosen in order that the dynamic of e sc is slower with respect to PI fast actuators. Remark 2. It is important to note that the bias may be cancelled by adding an adaptive skim. This will be done in a future work. Moreover due to practical implementation, the control is digital and some problems appear due to sampling frequency and quantification [START_REF] Monaco | Sampled-data energetic management of a fuel cell/capacitor system[END_REF].

EXPERIMENTAL RESULTS

To test and validate the performances of the proposed FC controller (17-18), a dedicated FC Control Benchmark is given by Fig. 5-(a)-(c) where the reference DC bus voltage v b is set equal to 50V and the load current i l varies between 1 and 20 A. Due to the safety constraints of academic laboratory use, this power cycle is representative of a reduced vehicle power demand, where the load requirement consists in raising and lowering power edges between 50 and 500W.

Experimental set-up.

Fig. 3 shows the experimental set-up used under which the proposed control is tested. The test bench is based on a Nexa fuel cell, Maxwell ultracapacitor bank, a programmable electronic load of 1800 W (i max = 150A/V max = 60V ) and a dSPACE DS1103 controller board. The supply sources are interconnected to the DC bus using two choppers built with standard IGBT modules The switching frequency of the PWM is set to 15 kHz. Electric characteristics of both sources are shown in table I. Concerning the load, the proposed FC controller uses only the measurement of load current and load voltage while the resistance load is an unknown variable.

Real time DS1103 card

Experiments without converter losses

Fig. 5 shows the response of the FC hybrid system according to the reference trajectories of a significant Benchmark Fig. 5-a)-b)-f), when the converter losses ( 16) are not included in controller u 2 (18), i.e. i d = 0. The parameters of u 2 (18) are chosen as follows: T f ast = 50ms and T slow = 2s. It can be remarked that the DC bus showed by Fig. 5-b) is not well regulated by the proposed fast controller (18) (because of converter losses which are not taking into account) and its value is fluctuate with respect to load current reference. To overcome this difficulty we have tried to decrease the value of T f ast but the error of DC bus regulation was near 20% which is not acceptable. Nevertheless a smooth response of the FC current1 (Fig. 5-(c)) is ensured during fast power demands of the load (Fig. 5-a)) so as to maintain the FC state of health. However, the energy balance is well achieved, characterized by a SC voltage (Fig. 5-f)) value reaching the desired set point (21V) at steady state by the proposed slow controller (17) with a small static error as stated in proposition 1 and demonstrated by proof 4.2.

I pac_max = 40A, Pente_Ipac_max = 30A/s T rapide = 50 ms, T lent = 2 s i load (t) [A] v BUS (t) [V] i FC (t) [A] v FC (t) [V] i UC (t) [A] v UC (t) [V] t (s) t (s) t (s) t (s) t (s) t (s) a) b) c) d) e) f)

Experiments with offline converter losses identification

In this section, the converter losses i d are included in the fast controller controller u 2 (18) as showed by Fig. 6. The converter losses are offline identified in steady state conditions. As the load current reference changes from zero to 20A, then the current converter losses is deduced from equation ( 16) where u 1max = 45A, V 0 = 1, 5V and R 0 = 170mΩ. The experimental results of the FC hybrid system according to the significant Benchmark Fig. 6-a)-b)-f) are displayed in Fig. 6. Compared to previous case (without converter losses) it can be noted that the DC bus showed in Fig. 5-b) is now well regulated by the proposed fast controller (18) and its error is less than 2% which is acceptable. The same conclusion is made for FC current and SC voltage (Fig. 6-f)) responses, as done for the case without converter losses.

Ipac_max = 45A, Pente_Ipac_max = 20A/s Trapide = 10 ms, Tlent = 0,25 s Vpertes = 3V, Rpertes = 0.17 ohm

CONCLUSION

Using both singular perturbation control and Lyapunov theory, this article demonstrates that it is possible to prove the stability of the FC/UCs system satisfying unknown variable load. Converter losses are taken into account in this proof. FC and UCs currents have to be slow and fast respectively which is not a major constraint since it respects the devices specifications. The experimental results, presented on a reduced scale test bench representative of a electrical vehicle application, highlight the interest of the proposed approach. However, the obtained results are linked to an off-line method identification test which enables to determine the converter losses. Consequently, our on-going work is to estimate the converter losses online using a suitable observer.

Fig. 1 .

 1 Fig. 1. Two converters parallel structure studied in this work. To reach the objectives, high-performance controllers are readily available based on the system state Jiang and Dougal (2006), fuzzy logic Kisacikoglu et al. (2009), proportional-integral controllers Thounthong et al. (2009), Davat et al. (2009), RST controller Caux et al. (2005), passivity Becherif (2006), flatness Payman et al. (2008) or model predictive control Vahidi et al. (2006). Alternative approaches exist such as optimal control Rodatz et al. (2005), dynamic programming Brahma et al.

Fig. 3 .

 3 Fig. 3. Experimental set-up.

Fig. 5 .

 5 Fig. 5. Experimental results without converter losses.Main variables (iL (t), i F C (t), i U C (t), v BU S (t), v U C (t))using the two tuning parameters (T f ast = 50ms and T slow = 2s) and the two bounds (I F CM AX = 40A and ( diF C dt M AX = 30A/s).

  results with converter losses. Main variables (i L (t), i F C (t), i U C (t), v BU S (t), v U C (t)) using the four tuning parameters (V 0 = 3V , R = 0, 170Ω, T f ast = 10ms and T slow = 0, 2s) and the two bounds (I F CM AX = 40A and ( diF C dt M AX = 30A/s).

The FC rate of change is nearly equal to 10A/s as for an FC example given by[START_REF] Thounthong | Energy management of fuel cell/battery/supercapacitor hybrid power source for vehicle applications[END_REF],[START_REF] Corbo | Pefc stacks as power sources for hybrid propulsion systems[END_REF]

Acknowledgment This research was supported by the PEPS project "GESE : Gestion Echantillonne des Systmes Energtiques", 2012.


Jianming Yang
email: jianming.yang@unistra.fr

Patrice Twardowski

Philippe Gérard

Y U Wenhui

Joël Fontaine

Chromatic analysis of harmonic Fresnel lenses by FDTD and angular spectrum methods

In this paper, we present a detailed and rigorous study of cylindrical harmonic Fresnel lenses (HFLs) using finite difference time domain method (FDTD) and angular spectrum method (ASM). The HFL is a kind of diffractive lens that can have maximum diffraction efficiency at several discrete harmonic wavelengths, which is suitable for some broadband applications. Previous studies on HFLs were investigated mainly in the domain of paraxial approximation. By using our proposed calculation method, we have determined the efficiency, focal length, maximum focus intensity and full width at half maximum (FWHM) of the focal spot for several harmonic numbers and for F-numbers of 0.5, 1 and 3. To compare with the paraxial approximation, we have presented the response to both s-polarized and p-polarized light with constant refractive index and real dispersive material BK7. Moreover, we have also analyzed the cases with oblique illumination. We have shown that the harmonic wavelengths do not change with F/# and that the diffraction efficiency and FWHM of the focus increase as F/# increases. New results on harmonic wavelengths shift and oblique angle of incidence response have been detailed.

INTRODUCTION

The diffractive Fresnel lens (DFL), also called kinoform, was firstly invented by a research group of IBM [START_REF] Lesem | The kinoform: a new wavefront reconstruction device[END_REF][START_REF] Jordan | Kinoform lenses[END_REF]. Since then, tremendous interests have been aroused due to its compactness. However, the dispersion of a DFL is large compared to optical glass (Abbe number of -3.5), excluding direct broadband application. As its dispersion sign is opposite to the normal optical glass, hybrid diffractive-refractive achromatic lenses have been studied and used [START_REF] Stone | Hybrid diffractive-refractive lenses and achromats[END_REF]. The optical power of the diffractive element is chosen much lower than the refractive one. These kinds of structures can reduce the dimension and weight of optical systems. In order to achromatize by only using diffractive optical elements (DOEs), two main paths have been proposed: stacks of DOEs with different Abbe numbers and harmonic lenses. However, the stack of DOEs can only increase the broadband efficiency. The Abbe number is still around -3.5 [START_REF] Kleemann | concepts for broadband high--efficiency DOEs[END_REF][START_REF] Mao | High diffraction efficiency of threelayer diffractive optics designed for wide temperature range and large incident angle[END_REF][START_REF] Piao | Diffraction efficiency change of multilayer diffractive optics with environmental temperature[END_REF][START_REF] Wang | Athermalization and thermal characteristics of multilayer diffractive optical elements[END_REF][START_REF] Xue | Diffraction efficiency evaluation for diamond turning of multilayer diffractive optical elements[END_REF]. A harmonic Fresnel lens (HFL) or multi-order diffractive lens can keep high diffraction efficiency and constant focal length at several discrete wavelengths [START_REF] Sweeney | Harmonic diffractive lenses[END_REF][START_REF] Faklis | Spectral properties of multiorder diffractive lenses[END_REF]. Recently, metasurfaces have been proposed to function as flat lenses, which paves a new way for optical components [START_REF] Groever | Meta-Lens Doublet in the Visible Region[END_REF][START_REF] Arbabi | Multiwavelength polarization-insensitive lenses based on dielectric metasurfaces with meta-molecules[END_REF][START_REF] Arbabi | Subwavelengththick lenses with high numerical apertures and large efficiency based on high-contrast transmitarrays[END_REF][START_REF] Aieta | Multiwavelength achromatic metasurfaces by dispersive phase compensation[END_REF][START_REF] Aieta | Aberration-free ultrathin flat lenses and axicons at telecom wavelengths based on plasmonic metasurfaces[END_REF]. However, the mass manufacturing for this kind of nanostructure is still challenging.

With broadband properties, HFLs can be applied in many situations, such as intraocular lenses [START_REF] Lane | Multifocal intraocular lenses[END_REF], contact lenses and optical detection in microflow cytometry [START_REF] Hunt | Kinoform microlenses for focusing into microfluidic channels[END_REF]. The previous works focused on HFLs have been carried out with scalar approach [START_REF] Sweeney | Harmonic diffractive lenses[END_REF][START_REF] Faklis | Spectral properties of multiorder diffractive lenses[END_REF][START_REF] Yan | Research on Dispersion Performance of Harmonic Diffractive Lenses[END_REF]. The results become inaccurate when F/# is smaller than 1 [START_REF] Bendickson | Scalar integral diffraction methods: unification, accuracy, and comparison with a rigorous boundary element method with application to diffractive cylindrical lenses[END_REF]. Based on the boundary element method (BEM) or the improved first Rayleigh-Sommerfeld method (IRSM), metallic cylindrical micro mirrors or cylindrical Fresnel lenses have been deeply studied [START_REF] Bendickson | Scalar integral diffraction methods: unification, accuracy, and comparison with a rigorous boundary element method with application to diffractive cylindrical lenses[END_REF][START_REF] Bendickson | Metallic surface-relief on-axis and off-axis focusing diffractive cylindrical mirrors[END_REF][START_REF] Ye | Applications of improved first Rayleigh-Sommerfeld method to analyze the performance of cylindrical microlenses with different f-numbers[END_REF][START_REF] Dong | Analysis of the focal characteristics of cylindrical lenses made of anisotropically dielectric material based on rigorous electromagnetic theory[END_REF][START_REF] Lin | Focusing performance of the closed-boundary cylindrical microlenses analyzed by the boundary element method[END_REF][START_REF] Ye | Rigorous electromagnetic analysis of metallic cylindrical focusing micromirrors with high diffraction efficiency, achromatic aberration and long focal depth[END_REF][START_REF] Ye | Improved first Rayleigh-Sommerfeld method applied to metallic cylindrical focusing micro mirrors[END_REF][START_REF] Ye | Improved first Rayleigh-Sommerfeld method for analysis of cylindrical microlenses with small fnumbers[END_REF]. However, to our knowledge, rigorous studies of harmonic Fresnel lenses (HFLs) has not been performed for the moment. Compared to ordinary Fresnel lenses, the broadband property of HFL is more effective.

The finite difference time domain (FDTD) is very suitable to calculate the broadband properties as only one simulation that covers the selected spectrum is needed. To reduce the computational time, the height of the FDTD region is limited to nearly the same height as the HFLs. For the field outside the HFLs, as the medium is homogeneous, we can use the angular spectrum method (ASM) to calculate the field [START_REF] Albarazanchi | Smart multifunction diffractive lens experimental validation for future PV cell applications[END_REF][START_REF] Di | Rigorous concept for the analysis of diffractive lenses with different axial resolution and high lateral resolution[END_REF][START_REF] Jiang | A rigorous unidirectional method for designing finite aperture diffractive optical elements[END_REF][START_REF] Mellin | Limits of scalar diffraction theory and an iterative angular spectrum algorithm for finite aperture diffractive optical element design[END_REF]. In this paper, we apply FDTD plus ASM to calculate the rigorous properties of HFLs. To reduce the calculation burden, we consider two-dimensional cases which correspond to cylindrical lenses. Different F/# are considered under p and s polarized illuminations. Results with on-axis and off-axis illuminations are also presented. 2)

2.

Ba

[9,

3)

the gth dinates, accordin here P is an arbit ofile in each groov ng to Fermat trary point on ve is given by:

4)

OP ILLUM

2 2 d i f r + , ( 7
)

OP ILLUM

  Fig ha dia

(

  Fig. 2. P A co surface wave on determi

a

  Th pose O is the o le, we have PF = file. Hence, the eq Fresnel lens, ea nt vortex Oi (see rtices should no ines the minimu gives the same m nimum F/# is 0.4 tical angle on the reflected if the F/PTICAL PROP MINATIONall the simulation ed numerically b solutions from L ), where n is the r ply perfect matche ur sides. All the tation of the field an example of out th ph=9, F/#=1, t he illumination w

1 .

 1 mbered i is a hyp grooves 0, 1 and image focal po L. The vortex coi 1 , Eq[START_REF] Xue | Diffraction efficiency evaluation for diamond turning of multilayer diffractive optical elements[END_REF]. As an exam origin of this F/# he lens material. n this limitation.HFL WITHgnetic near-field rcially available use a maximum of the correspon ) as the boundary re performed in the air we use A D and ASM field p the lens is 400µm o 0.5 µm. In the ( 6 )

  Fig wi λ=

  Fig an

  Fig be ca de ma ma sim

  Fig dif Th len