domain: spi phys chim

text:
V Guérineau

G Vilmart

N Dorval
email: nelly.dorval@onera.fr

Anda Julian-Jankowiak
email: aurelie.jankowiak@onera.fr

Comparison of ZrB 2 -SiC, HfB 2 -SiC and HfB 2 -SiC-Y 2 O 3 oxidation mechanisms in air using LIF of BO 2 (g)

Keywords: A) ceramic, (C) high temperature corrosion, (C) oxidation, (B) spectroscopy, (B) Laser Induced Fluorescence

The oxidation behaviour of ZrB 2 -SiC, HfB 2 -SiC and HfB 2 -SiC-Y 2 O 3 is studied using the real-time Laser-Induced Fluorescence (LIF) detection of BO 2 (g) evaporated from samples heated up to 1873 K in dry air. For each composition, the relevant steps of oxidation (silica formation, volatilization, etc.) can be precisely determined. Moreover, the influence of composition on the oxidation behaviour, and more precisely on the B 2 O 3 /SiO 2 ratio in the glassy phase, can be understood and described by monitoring the LIF signal from the BO 2 (g) radical. Thus, this technique has confirmed the potential of the HfB 2 -SiC-Y 2 O 3 composition as a Ultra High Temperature Ceramic (UHTC) material.

Introduction

Ultra-High Temperature Ceramics (UHTC) are extensively investigated as materials for applications in very high temperatures (>2073 K) and oxidising or corrosive atmospheres. Among this material family, UHTC borides, such as ZrB 2 and HfB 2 -based ceramics, are particularly interesting since the melting points of these materials are high, and the mechanical, thermal and chemical properties are excellent. Hence, they might be ideal candidates for developing thermal protection systems on hypersonic aerospace vehicles or atmospheric re-entry vehicles which require a good oxidation resistance under severe conditions.

Thus, studying and understanding the oxidation resistance of such materials is essential.

The influence of the SiC addition or of other additives (e.g. TaSi 2 , MoSi 2 ) has been largely studied [START_REF] Opila | Oxidation of ZrB 2 -and HfB 2 -based ultra-high temperature ceramics: effect of Ta additions[END_REF][START_REF] Talmy | Oxidation of ZrB 2 ceramics modified with SiC and group IV-VI transition metal diborides[END_REF][START_REF] Williams | Oxidation of ZrB 2 -SiC ultra-high temperature composites over a wide range of SiC content[END_REF][START_REF] Zhang | Processing and characterization of ZrB 2 -SiC W ultra-high temperature ceramics[END_REF][START_REF] Zhang | Oxidation of ZrB 2 and ZrB 2 -SiC ceramics with tungsten additions[END_REF] so that, the oxidation behaviour of ZrB 2 -SiC and HfB 2 -SiC is relatively well known up to 2073 K in air [START_REF] Carney | Oxidation behavior of zirconium diboride silicon carbide produced by the spark plasma sintering method[END_REF][START_REF] Fahrenholtz | Oxidation of ultra-high temperature transition metal diboride ceramics[END_REF][START_REF] Monteverde | The thermal stability in air of hot-pressed diboride matrix composites for uses at ultra-high temperatures[END_REF][START_REF] Carney | Oxidation resistance of hafnium diboride-silicon carbide from 1400 to 2000 °C[END_REF][START_REF] Shugart | SiC depletion in ZrB 2 -30 vol% SiC at ultrahigh temperatures[END_REF][START_REF] Shugart | Initial stages of ZrB 2 -30 vol% SiC oxidation at 1500°C[END_REF][START_REF] Hu | Oxidation mechanism and resistance of ZrB 2 -SiC composites[END_REF] . A few studies have also been carried out at higher temperatures and under other atmospheres than air [START_REF] Guérineau | Oxidation mechanisms under water vapour conditions of ZrB 2 -SiC and HfB 2 -SiC based materials up to 2400°C[END_REF][START_REF] Opila | Oxidation of ZrB 2 -SiC[END_REF][START_REF] Nguyen | Oxidation of ultrahigh temperature ceramics in water vapor[END_REF][START_REF] Alfano | Microstructural characterization of ZrB 2 -SiC based UHTC tested in the MESOX plasma facility[END_REF][START_REF] Marschall | Temperature jump phenomenon during plasmatron testing of ZrB 2 -SiC ultrahightemperature ceramics[END_REF][START_REF] Carney | Qualitative analysis of hafnium diboride based ultra high temperature ceramics under oxyacetylene torch testing at temperatures above 2100[END_REF][START_REF] Paul | UHTC-carbon fibre composites: preparation, oxyacetylene torch testing and characterisation[END_REF] , for instance in the presence of humidity or in an oxyacetylene flame.

Although the key role of the borosilicate layer [START_REF] Opila | Oxidation of ZrB 2 -and HfB 2 -based ultra-high temperature ceramics: effect of Ta additions[END_REF][START_REF] Tripp | Thermogravimetric study of the oxidation of ZrB 2 in the temperature range of 800° to 1500[END_REF][START_REF] Opeka | Oxidation-based materials selection for 2000°C + hypersonic aerosurfaces: Theoretical considerations and historical experience[END_REF] and finally of the ZrO 2 or HfO 2 layer [START_REF] Guérineau | Oxidation mechanisms under water vapour conditions of ZrB 2 -SiC and HfB 2 -SiC based materials up to 2400°C[END_REF][START_REF] Carney | Qualitative analysis of hafnium diboride based ultra high temperature ceramics under oxyacetylene torch testing at temperatures above 2100[END_REF][START_REF] Li | Effect of heat flux on ablation behaviour and mechanism of C/C-ZrB 2 -SiC composite under oxyacetylene torch flame[END_REF] has been demonstrated, a very few studies attempted the in situ investigation of its behaviour during an oxidation test. The thermogravimetric analysis (TGA) [START_REF] Tripp | Thermogravimetric study of the oxidation of ZrB 2 in the temperature range of 800° to 1500[END_REF] and the optical emission spectroscopy have been used elsewhere [START_REF] Sevast'yanov | Behavior of a sample of the ceramic material HfB 2 -SiC (45 vol %) in the flow of dissociated air and the analysis of the emission spectrum of the boundary layer above its surface[END_REF][START_REF] Playez | Optical emission spectroscopy during plasmatron testing of ZrB 2 -SiC ultrahightemperature ceramic composites[END_REF] . However, TGA only gives the sum of the weight gains and weight losses and only provides insights about the chemical evolution of the borosilicate layer. Investigations by emission spectroscopy imply that the chemical species are formed in sufficient quantity in an excited state to be detected through the radiative deactivation process, thus making this technique slightly sensitive. Previously, we have proposed the Laser-Induced Fluorescence (LIF) technique to detect gaseous BO 2 (g) radicals above the sample during oxidation of ZrB 2 and ZrB 2 -SiC samples [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] . These radical species come from the volatilization of B 2 O 3 from the glassy layers. The LIF technique enables the detection of BO 2 radicals with much better sensitivity, spatial and temporal resolutions than emission spectroscopy. Laser excitation of BO 2 is performed at 547.3 nm to promote the A 2  u (0,0,0) -X 2  g (0,0,0) vibronic band originating from the ground state. The subsequent fluorescence is detected via the A 2  u (0,0,0)X 2  g (1,0,0) red-shifted band at 580 nm [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] . In contrast to the TGA analysis, oxidation detection is delayed since B 2 O 3 has to be formed prior to the detection of BO 2 (g) [START_REF] Sevast'yanov | Behavior of a sample of the ceramic material HfB 2 -SiC (45 vol %) in the flow of dissociated air and the analysis of the emission spectrum of the boundary layer above its surface[END_REF] . However, the further evolutions of the glassy layer are markedly detected using the LIF technique such as the onset of the oxidation of SiC or the intense volatilization of the glassy layer at temperatures higher than 1500 K. This paper presents a comparison of the oxidation mechanisms of ZrB 2 -SiC, HfB 2 -SiC and HfB 2 -SiC-Y 2 O 3 by using the same real-time LIF detection strategy of evaporated BO 2 above the ZrB 2 , ZrB 2 -SiC samples heated in air, a technique which has been applied for the first time in our previous article [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] . The influence of the composition on the oxidation behaviour is described in details and the potential of the LIF technique for thermal oxidation detection is confirmed here in the case of the HfB 2 -based ceramics. After weighing, the powders are attrition-milled in cyclohexane using zirconia media (for ZS) or WC media (for HS and HSY), dried and sieved (50 µm). Fully dense samples (20 to 40 mm in diameter) are obtained using Spark Plasma Sintering (SPS, FCT HD25) in the MATEIS laboratory (Lyon, France). More details concerning manufacturing have already been reported elsewhere [START_REF] Guérineau | Oxidation mechanisms under water vapour conditions of ZrB 2 -SiC and HfB 2 -SiC based materials up to 2400°C[END_REF] . The bulk density and open porosity of these materials are measured by the Archimedes' method. Then, the densification level is calculated as the ratio of the apparent density to the theoretical density of the powder mixture.

Materials and apparatus

Materials

Sintering conditions and densification levels of the samples are reported in

LIF experiment

During the thermal oxidation runs in air of the samples heated with a continuous wave 2 kW CO 2 laser, BO 2 (g) is detected in the evaporation plume by monitoring its fluorescence signal at 580 nm induced by laser excitation tuned to 547.3 nm as mentioned above [START_REF] Johns | The absorption spectrum of BO 2[END_REF] . A very detailed description of the set-up, the method and the experimental protocol were reported in our previous study [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] .

The dedicated chamber is equipped with optical windows for laser-based measurements. The CO 2 laser irradiates the target at normal incidence angle from the top of the chamber (figure 1 a). The laser beam is focused on the target in order to irradiate the whole surface. The sample is mounted on an alumina holder in the chamber centre. The support can be translated upwards in order to adjust the height of the sample surface with respect to the probe laser beam axis. Surface temperature is measured with a bichromatic pyrometer in the 1273-2773 K range (Modline 6R-2565, Ircon). Synthetic air (Alphagaz 1) is injected in the chamber by means of a mass flow regulator, and is evacuated via a pumping system to maintain the pressure constant during the oxidation run.

The probe laser source is a tunable dye laser (Quantel, TDL50) pumped by a pulsed Nd:YAG laser (Quantel, YG781, 10 Hz repletion rate). The pulse duration is 6 ns and the linewidth is 0.07 cm -1 . The output energy is 15 mJ/pulse and is attenuated to 1 mJ using neutral density filters. The probe laser beam is crossing the chamber horizontally through anti-reflection coated windows in the visible range. The beam is focused with a spherical lens (f = 600 mm) to obtain a 1 mm laser spot at the sample location in the chamber. During a thermal oxidation run, the laser energy is recorded at the chamber exit by a photodetector (PD10, Ophir) in order to monitor the beam attenuation through the evaporation plume (absorption pathway in figure 1 b)

Fluorescence light is collected at right angle from the probe laser beam axis (fluorescence pathway in figure 1 b) by a 100 mm focal-length lens coupled to an optical fibre (1 mm in diameter). The fluorescence is imaged with a magnification of 1 onto a monochromator (Jobin-Yvon, H-20) equipped with a 1 mm wide exit slit (4 nm bandpass centred at 580 nm). The signal is then detected, amplified, time-integrated (34 ns temporal gate) and averaged (over 10 laser shots) by means of a photomultiplier tube (Photonis, XP2017B) and a Boxcar system (SRS, SR250) [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] . A Labview  interface is used to acquire integrated fluorescence signals via an acquisition PCI card (1.25 MHz, 16 bit) at 10 Hz.

Experimental protocol

Once the sample is positioned on the holder, the chamber is evacuated prior to synthetic air injection with a 1.3 l/min flow rate. The pressure is regulated to 0.1 MPa, then the heating of the sample and the recording of both the LIF signal and the laser transmission can start. The CO 2 -laser power is monitored manually, thus allowing to hold temperature steady to perform timely several adjustments and to ensure that the LIF signal is properly monitored [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] .

The nature and the thickness of each oxidised layer were identified through ex situ examinations in the Scanning Electron Microscope (SEM, DSM 962) and by Energy Dispersive Spectroscopy (EDS) of both the surface and the polished cross-sections of the oxidised materials.

In addition, TGA experiments (SETSYS Evolution 16/18, SETARAM) are carried out to measure the mass variations during oxidation of each composition. Samples (2 × 2 × 4 mm 3 in size), are slowly heated (2 K/min) from room temperature to 1823 K under a 50 ml/min air flow.

Results

The variation of amplitude of LIF and laser transmission signals as a function of time for ZS, HS and HSY samples is shown in figure 2. As for ZS samples [START_REF] Guérineau | In situ study of the oxidation of ZrB 2 and ZrB 2 -SiC materials by monitoring the LIF signal of BO 2 radicals[END_REF] , four key steps can be identified for HS and HSY samples and are numbered from ① to ④ in figure 2 but with noticeable differences in the amplitude of signal or in the temperature range associated with each step, as compared in Table 3.

Step ① corresponds to the appearance of the LIF signal indicating that a detectable quantity of BO 2 (g) has been formed. This gaseous species comes from the volatilization of the B 2 O 3 (l) glassy layer formed at the surface of the sample ( Reactions 1 and 2).

MeB 2 + 5/2 O 2  MeO 2 (s) +B 2 O 3 (l) (1)
B 2 O 3 (l)  B 2 O 3 (g) (2)
where Me = Zr or Hf.

For all samples, weight gains are observed at lower temperature with TGA (1013 K for ZS and 948 K for HS and HSY) as shown in figure 3. Thus, as already observed, a threshold in boria amount and thus in B 2 O 3 (g) has to be overstepped before the appearance of the LIF signal. At this temperature (<1473 K), SiC is not yet oxidised (figure 2) [START_REF] Williams | Oxidation of ZrB 2 -SiC ultra-high temperature composites over a wide range of SiC content[END_REF] . Thus, the following steady rise in the LIF signal between steps ① and ② is attributed to the enhanced volatilization of B 2 O 3 , leading to an increased quantity of BO 2 (g). This is confirmed by the slower increase in weight gain as measured by TGA (figure 3).

An inversion of trend is observed at step ②: the LIF signal slowly decreases despite the continuous rise of sample heating, up to step ③, indicating that less BO 2 (g) is formed in this temperature range (Table 3). This striking observation is associated with the beginning of the passive oxidation of SiC through reaction 3. Therefore, the enrichment of the glassy phase with silica increases its stability (higher viscosity) compared to the previous boria glassy layer (reaction 4).

SiC + 4O 2  SiO 2 (s/l) + CO 2 (g) (3) SiO 2 (l)+B 2 O 3 (l)  {SiO 2 :B 2 O 3 }(l) (4)
Several studies have already shown that a B 2 O 3 concentration gradient is formed in the borosilicate glassy layer as preferential evaporation of B 2 O 3 is observed at the surface [START_REF] Karlsdottir | Oxidation of ZrB 2 -SiC: influence of SiC content on solid and liquid oxide phase formation[END_REF][START_REF] Karlsdottir | Convection patterns in liquid oxide films on ZrB 2 -SiC composites oxidized at a high temperature[END_REF][START_REF] Karlsdottir | Formation of oxide scales on zirconium diboridesilicon carbide composites during oxidation: relation of subscale recession to liquid oxide flow[END_REF] .

This lower concentration, and thus the lower activity of B 2 O 3 at the surface explains the LIF signal decrease [START_REF] Shugart | Initial stages of ZrB 2 -30 vol% SiC oxidation at 1500°C[END_REF][START_REF] Sevast'yanov | Behavior of a sample of the ceramic material HfB 2 -SiC (45 vol %) in the flow of dissociated air and the analysis of the emission spectrum of the boundary layer above its surface[END_REF][START_REF] Karlsdottir | Formation of oxide scales on zirconium diboridesilicon carbide composites during oxidation: relation of subscale recession to liquid oxide flow[END_REF] . Passive oxidation of SiC is also detected using TGA, but at a higher temperature.

Sample Detection of LIF signal (T 1 , ①)

Step ② onset (T 2 ).

Step ③ onset (T 3 ) Between steps ③ and ④ small bumps in the LIF signal, immediately followed by a decrease in the signal are observed. This is called the "dynamic phase" and can be explained by a power peak which makes the oxide layer less stable. Thus, the B 2 O 3 -rich glassy layer, where the activity of B 2 O 3 is higher, is exposed and, therefore, the volatilization of B 2 O 3 is favoured. The glassy layer is thinner and less protective against oxidation, so that the oxidation of the material underneath is favoured and, as a consequence, quickly reforms a stable borosilicate glassy layer, and the LIF signal decreases. The subsequent signal bumps have the same origin (power peak) and the same consequences. Between each bump, the average intensity of the LIF signal increases as the temperature is continuously increasing.

ZS 1263 K 1553 K 1643 K HS 1223 K 1433 K 1573 K HSY 1223 K* 1523 K 1588 K
Compared to TGA, between 1300 and 1673 K, the weight gain increases strongly for HSY, and more smoothly for HS. ZS has a different TGA behaviour, especially at 1473 K. Instead of a decrease in the mass gain slope, which is observed for HS and HSY, there is an increase in this slope. This indicates a continuous oxidation of the materials.

Finally, the high intensity LIF signal peak coupled with a dramatic drop in the laser transmission signal is observed in step ④. At that moment, the probe laser beam is strongly attenuated because of the formation of a dense smoke of oxide particles in the plume, and the LIF signal is no longer detected. This step is also triggered by a power peak, but at this high temperature, the volatilization of the glassy layer is very intense, and the active oxidation of SiC is also favoured via the Reaction 5.

SiC + O 2 (g)  SiO(g) + CO(g).

(

) 5
The drop in the transmission signal is linked to the catastrophic volatilization of the glassy layer leading to a high amount of B 2 O 3 and SiO particles in the atmosphere just above the sample. However, at 1673 K, the mass gain measured by TGA increases noticeably for all the compositions (figure 3), which confirms that the volatilization of the borosilicate glassy phase promotes the oxidation of the material underneath.

These general trends are observed for all compositions. However, some differences are noticeable, such as the temperatures (T1, T2 and T3) of the onset of the three steps ① to ③, the LIF signals between steps ③ and ④ and the LIF and transmission signals at step ④.

These differences are discussed in the next section.

Discussion

Formation of the boria glassy phase

Considering Table 3, it can be noticed that the temperature at which the LIF signal is first detected (T 1 ) is slightly higher for ZS than for HS and HSY. Interestingly, weight gain is detected at lower temperatures in the TGA signal (around 1113 K for ZS and around 948 K for HS and HSY) indicating that the LIF signal BO 2 (g) is detected once a boria glassy layer is formed at the surface of the sample (figure 3). The later apparition of the B 2 O 3 boria glassy phase and thus of the BO 2 (g) detection for ZS samples is observed using both experiments.

This can be explained considering the volatility diagrams of ZrB 2 and HfB 2 [START_REF] Fahrenholtz | The ZrB 2 volatility diagram[END_REF][START_REF] Li | Thermodynamic calculation of HfB 2 volatility diagram[END_REF] and the diagrams of Ellingham (figure 4). It is clearly shown that oxidation of HfB 2 is favoured at a lower temperature than for ZrB 2 .

Formation of the borosilicate layer

The temperature at which the LIF signal stops increasing (T 2 ) is higher for ZS and HSY than for HS (~100 K in difference). This phenomenon is worth noting since T 2 marks the moment when SiO 2 (l) effectively integrates the boria network to form a borosilicate layer, making it more stable and less subject to volatilization. SEM micrographs of the three materials after oxidation are shown in figures 5 and 6.

The oxidised samples all develop a classic 3-layer oxide scale above the bulk material, from top to bottom: a SiO 2 -rich glassy layer, a MeO 2 layer (Me = Zr or Hf) and a SiC-depleted MeB 2 layer. The SiC-depleted layer is formed by active oxidation of SiC (reaction 5), and the reaction product SiO(g) can diffuse towards the surface [START_REF] Fahrenholtz | Thermodynamic analysis of ZrB 2 -SiC oxidation: formation of a SiCdepleted region[END_REF] . Once SiO(g) reaches the borosilicate layer, it can react with dissolved O 2 (g) to form SiO 2 (l), thus participating in the replenishment of the glassy layer [START_REF] Guérineau | Oxidation mechanisms under water vapour conditions of ZrB 2 -SiC and HfB 2 -SiC based materials up to 2400°C[END_REF][START_REF] Parthasarathy | Modeling oxidation kinetics of SiC-containing refractory diborides[END_REF][START_REF] Simonenko | Promising ultra-high-temperature ceramic materials for aerospace applications[END_REF] (reaction 6).

SiO(g) + ½ O 2 (g)  SiO 2 (l) (6)
From the SEM observations, it should be noted that the SiC-depleted layer has a very different relative thickness depending on the nature of the sample. HS sample develops a very thick SiC-depleted HfB 2 layer (~500 µm, 90% of the total oxide scale), HSY develops a thick SiC-depleted HfB 2 layer (~250 µm, 60% of the total oxide scale) and ZS develops a comparatively thin SiC-depleted ZrB 2 layer (~60 µm, 25% of the total oxide scale) [START_REF] Shugart | SiC depletion in ZrB 2 -30 vol% SiC at ultrahigh temperatures[END_REF][START_REF] Karlsdottir | Formation of oxide scales on zirconium diboridesilicon carbide composites during oxidation: relation of subscale recession to liquid oxide flow[END_REF][START_REF] Fahrenholtz | Thermodynamic analysis of ZrB 2 -SiC oxidation: formation of a SiCdepleted region[END_REF] , as reported in Table 4. 5 a andb). Another aspect that should be taken into consideration is that, at the experiment temperature (T<1900 K), tetragonal-ZrO 2 is a better oxygen conductor than monoclinic-HfO 2 [START_REF] Trubelja | Ionic conductivity of the fluorite-type hafnia-RE 2 O 3 solid solutions[END_REF][START_REF] Corman | Phase equilibria and ionic conductivity in the system ZrO 2 -Yb 2 O 3 -Y 2 O 3[END_REF] . Thus, this point accentuates the fact that the active oxidation of SiC is not favoured in the ZS composition.

SiC

The dynamic phase

The beginning of step ③ always corresponds to a peak in the LIF signal (figure 2) and a strong increase in the weight gain (figure 3). In this temperature range (>1573 K) the volatilization of B 2 O 3 becomes very intense [START_REF] Tripp | Thermogravimetric study of the oxidation of ZrB 2 in the temperature range of 800° to 1500[END_REF] For ZS, during the dynamic phase, the average level of the LIF signal increases after each peak (figure 2). For this sample, the stabilisation provided by SiO 2 is less effective since the glassy layer is B 2 O 3 -rich, and the volatilization of B 2 O 3 is, therefore, thermo-activated throughout the dynamic phase. The laser transmission signal also decreases steadily;

consequently, there is a continuous formation of B 2 O 3 droplets. Compared to HSY, the LIF signal takes more time to reach the low level after each signal peak (~180 s compared to ~60 s for HSY). This means that a longer time is required for the sample to regenerate a glassy layer with a good stability. This is mostly due to the fact that the glassy layer of the oxidised ZS sample is B 2 O 3 -rich, which means that it is more subjected to volatilization, and that a larger quantity is required in order to get the same protection, due to the lower viscosity of B 2 O 3 -rich borosilicate glass. On the contrary, the glassy layer developed by HSY contains relatively more SiO 2 , so that it is more stable, and the material developed a thick SiC-depleted layer, which means that the replenishment of the layer is made with a SiO 2 -richer glass, explaining why the glassy layer more rapidly attains stability.

The LIF signal behaviour of HS during oxidation is similar to that of ZS. Each LIF signal peak is followed by an increase in the average LIF signal level. But, contrary to the other samples, the laser transmission signal was almost constant throughout the experiment.

The LIF peaks are also very sharp. To explain these observations, it should be remembered that the borosilicate glass formed by HS is likely very SiO 2 -rich, since the SiC-depleted zone is 9 times thicker than the HfO 2 layer. The deficiency in B 2 O 3 explains that there is no or not enough formation of B 2 O 3 droplets and, therefore, the transmission level is unchanged. It can be noticed that the LIF signal peak at 2100 s (1783 K) is accompanied by a short drop in the transmission level, thus meaning that, at this temperature, some B 2 O 3 droplets are formed in a short time. As evidenced in figure 5, the borosilicate layer is likely thin and SiO 2 -rich during the test. Contrary to the HSY sample, where the height of LIF signal peaks was about a half of the average LIF signal level, the LIF signal peaks observed for HS are much more prominent, and can be ten times higher than the average LIF signal level. Due to the B 2 O 3 deficiency and, therefore, the thinness of the HS borosilicate layer, this glassy phase is much more sensitive to a power peak. After the power peak, the borosilicate layer is even thinner, and this triggers the oxidation of the underlying HfB 2 , which will refill the borosilicate glass at the surface.

Intense volatilization of B 2 O 3

Step ④ also illustrated some differences in the oxidation behaviour of the different materials. For the ZS sample, the final LIF signal peak was accompanied by a drastic drop in transmission. For the HS material, this step is characterised by a moderate LIF signal peak and a decrease in the transmission signal. And, finally, for HSY, there is a high intensity LIF signal peak, but there is only a small laser transmission decrease (from 0.7 to 0.6).

Regarding

Conclusion

The oxidation behaviours of ZrB 2 -SiC, HfB 2 -SiC and HfB 2 -SiC-Y 2 O 3 materials were studied using the real-time LIF detection of evaporated BO 2 (g), above the samples heated up to 1873 K in a dry air atmosphere. For each composition, the starting point of the influence of SiC on the stability of the glassy layer and the catastrophic volatilization of the glassy layer can be precisely determined. Moreover, the influence of the composition on the oxidation behaviour can be understood and described using this technique.

In order to form a stable glassy layer, the B 2 O 3 /SiO 2 ratio must be finely tuned. If this ratio is too high, as it seems to be the case for ZS, then the glassy layer is not stabilised enough by SiO 2 , and the glassy layer will completely volatilize at ~1873 K. If it is too low, for instance for the HS material, then there is not enough "matter" to create a thick glassy layer that will completely cover the surface of the sample. On the contrary, a good compromise is achieved for HSY: this sample is not only able to generate a thick borosilicate layer, but the preferential volatilization of B 2 O 3 will trigger the precipitation of HfO 2 particles, that will in turn increase the stability of the glassy layer in the considered temperature range.

As it was already shown, LIF is a promising method, which can be quantitative, to study in situ, and thus continuously, the oxidation behaviour of UHTC. Moreover, the influence of the composition and its ability to form a protective layer can be studied and determined. Thus, the interest of the LIF technique in the detection of thermal oxidation is clearly confirmed. Tables captions:

Figures captions

  and explains the constant increase in the average level of LIF signal observed in figure2. However, in the present cases, the volatilization of boria is attenuated by the presence of silica in the glassy layer and the composition of the sample induces some differences in the LIF and transmission signals. For HSY, once the dynamic phase has begun, the average level of the LIF signal remains almost constant, despite the continuous increase in temperature. In parallel, the transmission slowly decreases from 1473 to 1823 K. This means that there is a continuous formation of B 2 O 3 droplets and, therefore, a steady volatilization of the borosilicate phase. But, the HSY material developed a very thick borosilicate layer containing a large amount of recrystallised HfO 2 particles (figure6 b). The thickness of the glassy layer explains its resilience, i.e. its capacity to recover an average LIF signal level after the power peak. Moreover, Talmy et al.[START_REF] Talmy | Oxidation of ZrB 2 ceramics modified with SiC and group IV-VI transition metal diborides[END_REF] have shown that the presence of immiscible phases in a glassy layer could increase its viscosity and, therefore, limit the diffusion of oxygen. During its convection, the B 2 O 3 -rich glass is able to solubilise Y 2 O 3 and HfO 2 . When reaching the surface, B 2 O 3 preferentially evaporates, so that Y 2 O 3 and HfO 2 attain their saturation concentration and precipitate. Therefore, the continuous volatilization of B 2 O 3 throughout the oxidation test favours the continuous precipitation of the solubilised components, which, in turn, increases the stability of the glassy layer.

  ZS sample, the LIF and transmission signals are associated with a catastrophic volatilization of the B 2 O 3 -rich borosilicate layer. This is confirmed by the final microstructure of the oxidised ZS sample (figure 5 a), where the glassy layer is very thin or non-existent. However, SEM micrographs of the surface of the oxidised ZS sample (figure 6 a and b) clearly show two features: the first is that, despite the thinness of the glassy layer, itstill almost covers the whole surface of the sample, and the second is that the glassy layer is filled with recrystallized ZrO 2 (small particles and large dendritic structures). The crystallites in the glassy layer constitute an indirect proof of the intense volatilization of B 2 O 3 . The fact that the borosilicate layer covers the entire surface reveals the low viscosity of this layer, probably due to its poverty in SiO 2 . The step ④ of HSY is also characterised by an intense LIF signal peak, but not a drastic transmission drop. As explained earlier, the thickness of the glassy layer developed by HSY allows the borosilicate glass to be more resilient towards the power peak. Even if the LIF signal of BO 2 increases, there are not enough B 2 O 3 droplets formed to drastically lower the laser transmission and, therefore, the volatilization cannot be qualified as catastrophic. This is confirmed by the SEM observations; the glassy layer is thick (figure5 e and f) and covers the whole surface (figure6 e and f). As in the case of ZS, the overwhelming presence of crystallites and dendritic structures of Y 2 O 3 -doped HfO 2 at the surface constitutes an evidence of the constant volatilization of B 2 O 3 during the oxidation test.The case of HS is very different. As explained above, the step ④ of HS is not characterised by an intense LIF signal peak. But a closer examination of the transmission signal reveals that there is a very short but intense drop at 2400 s (1833 K). SEM observations of the HS microstructure showed that the glassy layer was thin (figure5 c and d) and did not cover the whole surface of the sample (figure6 c and d). The fact that there is simply not enough B 2 O 3 in the remaining borosilicate layer to generate enough detectable BO 2 can explain that there is no intense LIF signal at the end. The lack of B 2 O 3 means that the viscosity of the borosilicate layer increases, thus, explaining that the glassy layer does not cover the entire surface of the sample.

Figure 1 :

 1 Figure 1: (a) a schematic view of the chamber and (b) LIF set-up.

Figure 2 :

 2 Figure 2: Surface temperature (left axis), normalised LIF and transmission signals (right axis) vs. time of (a) ZrB 2 -SiC, (b) HfB 2 -SiC and (c) HfB 2 -SiC-Y 2 O 3 materials during an oxidation test in dry air under 0.1 MPa. Key steps (numbered 1 to 4) are described in the text.

Figure 3 :

 3 Figure 3: Mass variation vs. time of HS, ZS and HSY. Temperature ramp is 2 K /min. The TGA is performed under a 50 ml/min flux of air. Inset: a zoom on TGA signals between 800 and 1200 K.

Figure 4 :

 4 Figure 4: Diagram of Ellingham for Si, Zr and Hf elements.

Figure 5 :

 5 Figure 5: SEM micrographs of polished cross-sections of oxidised samples: (a, b) ZS, (c, d) HS and (e, f) HSY. The different layers are (1) borosilicate glass, (2) MeO 2 (Me= Zr or Hf), (3) SiC-depleted MeB 2 and (4) bulk material.

Figure 6 :

 6 Figure 6: SEM micrographs of the surfaces of oxidised samples: (a, b) ZS, (c, d) HS and (e, f) HSY materials. In grey, the borosilicate layer and in white, zirconia or hafnia inclusions and dendritic structures. (f) Recrystallised HfO 2 is omnipresent (in white) as dendrites and crystallites in the borosilicate layer (grey).

Figure 1 :

 1 Figure 1: (a) a schematic view of the chamber and (b) LIF set-up.

Figure 2 :Figure 3 :

 23 Figure 2: Surface temperature (left axis), normalised LIF and transmission signals (right axis) vs. time of (a) ZrB 2 -SiC, (b) HfB 2 -SiC and (c) HfB 2 -SiC-Y 2 O 3 materials during an oxidation test in dry air under 0.1 MPa. Key steps (numbered 1 to 4) are described in the text.

Figure 4 :Figure 5 :Figure 6 :

 456 Figure 4: Diagram of Ellingham for Si, Zr and Hf elements.

Table 1 : Grade, purity and particle size of starting powders.

 1 Three UHTC compositions are studied: ZrB 2 + 20 vol% SiC, labelled ZS, HfB 2 + 20 vol% SiC, labelled HS and HfB 2 + 20 vol% SiC + 3vol% Y 2 O 3 , labelled HSY. Properties of the starting powders are listed in Table 1.

	Powder	Provider	Grade	d 50 (µm)	Purity
	ZrB 2	H.C. Starck	A	2.8	>97.8
	HfB 2	H.C. Starck	A	7.6	>97
	SiC	H.C. Starck	BF12	0.6	>98.5
	Y 2 O 3	Ampere Industrie -	<5	99.99

Table 2 .

 2

	Composition	Sintering process	Applied pressure	Thermal treatment	Open porosity	Densification rate / th
	ZrB 2 + 20 vol% SiC	SPS	28 MPa	2323 K, 2 min	0.1 %	>99 %
	HfB 2 + 20 vol% SiC	SPS	16 MPa	2223 K, 5 min	1.8 %	97 %
	HfB 2 + 20 vol% SiC + 3vol% Y 2 O 3	SPS	7 MPa	2153 K, 5 min	0.1 %	98.5 %

Table 2 : Sintering conditions and open porosity of the studied UHTC materials.

 2

Table 3 : Temperatures of the onset of the three key steps T1, T2 and T3. *extrapolated temperature.

 3

Table 4 : Thickness of each layer in the oxidised layers for all the tested compositions (average values from 10 measurements).

 4 O2 should be low enough to oxidise SiC without oxidising HfB 2[START_REF] Fahrenholtz | Thermodynamic analysis of ZrB 2 -SiC oxidation: formation of a SiCdepleted region[END_REF][START_REF] Shugart | SiC depletion in ZrB 2 -30 vol% SiC at ultrahigh temperatures[END_REF][START_REF] Opila | SiC recession caused by SiO 2 scale volatility under combustion conditions: II, thermodynamics and gaseous-diffusion model[END_REF] . If the glassy layer is less protective towards oxygen diffusion, the formation of a SiC-depleted layer is hindered in favour of the oxidation of both SiC and HfB 2 . This explains why, when Y 2 O 3 is added to a HfB 2 -SiC composition, the stabilisation provided by SiO 2 (l) is triggered at a higher temperature. Then, comparing ZS and the two other compositions, ZS exhibits the thinnest SiC-depleted layer and the thickest glassy layer, indicating that the glassy layer formed at the top surface of ZS samples is B 2 O 3 -rich compared to the glassy phases developed by HS and HSY. Thus, the glassy layer of ZS is more permeable towards oxygen, thus promoting the oxidation of both SiC and ZrB 2 (figure

		-depleted layer	Glassy phase layer	Total oxidised layer
		(µm)	(µm)	(µm)
	ZS	60	180	240
	HS	500	55	555
	HSY	250	170	420

Moreover, HS develops a very thin HfO 2 layer compared to the SiC-depleted HfB 2 layer or to other compositions. The composition of the glassy layer surface is directly dependent on the different reaction products that replenish it during the oxidation test. B 2 O 3 (l) only comes from the oxidation of HfB 2 whereas SiO 2 (l) only comes from (directly or via active oxidation) the oxidation of SiC. Thus, the borosilicate layer developed by HS will likely be enriched in SiO 2 (l) compared to HSY and ZS. This enrichment leads to an early stabilisation of the glassy layer, and that is why step ② happens at a lower temperature for HS. For HSY, step ② happens at 1523 K, closer to the step ② of ZS (1553 K). First, HSY develops a thinner SiCcapability to act as an oxygen diffusion barrier, since a low viscosity favours the diffusion of oxygen through the borosilicate layer. To form a SiC-depleted layer, the oxygen partial pressure P

Table 1 : Grade, purity and particle size of starting powders.Table 2 : Sintering conditions and open porosity of the studied UHTC materials.

 12

Table 3 : Temperatures of the onset of the three key steps T1, T2 and T3. *extrapolated temperature.Table 4 : Thickness of each layer in the oxidised layers for all the tested compositions (average values from 10 measurements).

 34

Acknowledgements

The authors would like to thank M. Bejet and T. Schmid (ONERA) for modifying the oxidation chamber and, finally, Dr P. Beauchêne and Dr. R. Valle (ONERA) for their kind and efficient assistance.

Data availability

The raw/processed data required to reproduce these findings cannot be shared at this time due to legal or ethical reasons.


Alexis Godart
email: alexis.godart@utbm.fr

Hervé Manier
email: herve.manier@utbm.fr

Christelle Bloch
email: christelle.bloch@univ-fcomte.fr

Marie-Ange Manier
email: marie-ange.manier@utbm.fr

Marie Manier Métaheuristique

Métaheuristique hybride pour un problème de tournées de véhicules multimodales et multi-visites

Keywords: Pickup & Delivery Problem, logistique urbaine, mobilité urbaine, synchronisation, mutualisation, collaboration

teaching and research institutions in France or abroad, or from public or private research centers.

Introduction

Nos travaux de recherche explorent une variante du problème de tournées de véhicules avec fenêtres de temps, qui s'inspire du Vehicle Routing Problem with Pickup & Delivery [START_REF] Parragh | A survey on pickup and delivery problems. Part I: Transportation between customers and depot[END_REF]- [START_REF] Ropke | An Adaptive Large Neighborhood Search Heuristic for the Pickup and Delivery Problem with Time Windows[END_REF] et du Dial-A-Ride Problem [START_REF] Masson | The dial-a-ride problem with transfers[END_REF], et s'inscrivent dans le cadre du projet de recherche Mobility In Smart Cities (MISC). La modélisation proposée tient compte des contraintes et caractéristiques intrinsèques au milieu urbain, parfois communes au transport de personnes et de marchandises.

L'originalité de nos travaux repose sur la possibilité d'effectuer des transferts avec stockage, avec les contraintes de précédence entre les deux modes de transport, ou bien d'effectuer des transferts sans stockage et en synchronisant donc ces deux modes de transport, le tout dans le cadre du développement du transport multimodal. Aussi, nous offrons la possibilité pour chaque mode de transport de visiter plusieurs fois un même site (multi-visites) afin, par exemple, de différer les opérations, et de proposer un choix plus large pour ordonnancer la prise en charge des demandes.

Méthodes

Une formulation mathématique de ce problème a été proposée à l'aide d'un programme linéaire en nombres mixtes [START_REF] Godart | Un modèle d'optimisation linéaire pour résoudre un problème riche de tournées de véhicules en milieu urbain[END_REF]. Les résultats obtenus permettaient d'identifier les gains sur la distance totale parcourue, et parfois même en termes de faisabilité, mais mettaient aussi en évidence l'explosion combinatoire empêchant le passage à l'échelle ; on propose alors une résolution par métaheuristique hybride, à base d'algorithme génétique (1), d'algorithme glouton (2) et de programmation linéaire en nombres mixtes (3).

Le problème initial considère simultanément : a. l'affectation des demandes aux véhicules ; b. la création des tournées de véhicules, qui comprend l'ordonnancement des opérations de ramassage/dépose, la proposition de transfert des demandes avec/sans stockage en vérifiant les relations d'interdépendances temporelles et spatiales entre les véhicules, les contraintes temporelles de visite sur chaque site, sur chaque demande, ainsi que les contraintes de capacité des véhicules et des points de transfert. Dans notre modèle, les demandes sont exprimées par une origine et une destination, avec une quantité associée.

Notre proposition consiste à résoudre d'abord le sous-problème d'affectation des demandes aux véhicules via un codage génétique (1) en variables binaires (figure). Le problème étant non-sélectif, chaque demande doit être affectée au moins à un véhicule. Cependant, dans le cadre du transport multimodal, il est possible d'associer une demande à plusieurs véhicules, auquel cas nous devons répartir les opérations de ramassage, de dépose, et de transfert (localisation et date) entre ces véhicules.

Figure -Exemple d'individu pour l'affectation des demandes aux véhicules

Une fois les affectations obtenues, un algorithme glouton (2) de liste est utilisé afin d'obtenir, pour chaque demande, un classement des véhicules en fonction de la proximité physique et temporelle du site origine et destination de la demande. Puis, dans le cas où plusieurs véhicules sont affectés à une demande, on se propose de transférer la demande de telle sorte que chaque véhicule effectue une partie de l'itinéraire de la demande. Enfin, on détermine le point de transfert à utiliser parmi ceux disponibles, pour chaque couple de véhicules, en utilisant la même stratégie.

Connaissant désormais les opérations à effectuer pour chaque véhicule, il ne nous reste plus qu'à les ordonnancer tout en respectant les contraintes temporelles, de capacité et de précédence. Pour ce faire, nous proposons de réutiliser le programme linéaire (3) développé au début de nos travaux pour résoudre le problème d'ordonnancement des opérations à l'optimalité, en minimisant la distance parcourue par les véhicules. Les affectations des demandes et des opérations aux véhicules sont soumises au modèle linéaire sous la forme de contraintes sur les variables de décision. Nous ajoutons les contraintes temporelles (dates d'arrivée et de départ au plus tôt et au plus tard) qui conditionnent la tournée de chaque véhicule et l'itinéraire de chaque demande en utilisant la propagation de contraintes sur les fenêtres de temps.

Expérimentation & Résultats

Notre approche a été testée sur ce problème dans sa variante statique, avec un ensemble de données déterministes. Nous comparons les résultats obtenus dans notre approche avec ceux obtenus par méthode exacte grâce au programme linéaire sur 5 instances de taille variable. Puis nous estimons la capacité de notre modèle à passer à l'échelle, en testant sur des instances de plus grosse taille, inspiré de la littérature. Le modèle linéaire est implémenté sur IBM CPlex® en version 12.7.1 et résolu sur un processeur Intel Core i7 4 coeurs cadencés à 3,5 GHz. D'autres instances pourront illustrer l'efficacité de la méthode sur des instances de tailles plus conséquentes. L'étape suivante de nos travaux sera d'améliorer le temps de calcul en explorant d'autres méthodes approchées afin de construire d'autres métaheuristiques, sans dégrader la qualité des solutions proposées. L'une des perspectives de nos travaux est l'exploration de ce problème dans sa variante dynamique.


On the Collision Property of Chaotic Iterations Based Post-Treatments over Cryptographic Pseudorandom Number Generators

Luigi Marangio, Christophe Guyeux and Jacques M. Bahi

Femto-ST Institute, UMR 6174 CNRS Université de Bourgogne Franche-Comté France

Abstract-There is not a proper mathematical definition of chaos, we have instead a quite big amount of definitions, each of one describes chaos in a more or less general context. Taking in account this, it is clear why it is hard to design an algorithm that produce random numbers, a kind of algorithm that could have plenty of concrete appliceautifat (anul)d bions. However we must use a finite state machine (e.g. a laptop) to produce such a sequence of random numbers, thus it is convenient, for obvious reasons, to redefine those aimed sequences as pseudorandom; also problems arise with floating point arithmetic if one wants to recover some real chaotic property (i.e. properties from functions defined on the real numbers). All this considerations are synthesized in the problem of the Pseudorandom number generators (PRNGs). A solution to these obstacles may be to post-operate on existing PRNGs to improve their performances, using the so-called chaotic iterations, i.e., specific iterations of a boolean function and a shift operator that use the inputted generator. This approach leads to a mathematical description of such PRNGs as discrete dynamical systems, on which chaos properties can be investigated using mathematical topology and measure theory. Such properties are well-formulated, and they allow us to characterize which functions improves the sensitivity to the seed, the expansivity, the ergodicity, or the topological mixing of the generator resulting from such a post-processing. Experience shows that choosing relevant boolean functions in these chaotic iterations improves the randomness of the inputted generator, for instance when considering the number of statistical tests of randomness passed successfully. If we focus on the cryptographical application of PRNGs, there are two main classical notions to be considered, namely collision and avalanche effect. In this article, we recall the chaotic properties of the proposed post-treatment and we study the collision property in families of pseudorandom sequences produced by this process.

Index Terms-Pseudorandom numbers generator, cryptography.

I. INTRODUCTION

To simulate some well-known chaotic real functions, such as a logistic map or the Arnold's cat map, is the main idea of many algorithms developed until now which aspire to be good PRNGs, programs with a lot of applications for instance in cryptography. Indeed, it is reasonable to think that elements of chaos can improve the random-like quality produced by such algorithms. A theoretical result that establish a link between a chaotic map on R and its floating-point counterpart, it has not yet been stated. Conversely, there are some results in the opposite direction, indicating that the numerical truncation may change drastically the statistical properties of orbits (see, e.g., [START_REF] Galatolo | Probability, statistics and computation in dynamical systems[END_REF], [START_REF] Guiheneuf | Dynamical properties of spatial discretizations of a generic homeomorphism[END_REF]). A first attempt to avoid the use of a chaotic real map, is to consider the socalled parallel asynchronous linear iterations (PALI). Even if the domain of definition is the floating point one, the effect of round-off errors in this situations was analyzed (for instance) in [START_REF] Bahi | Asynchronous iterative algorithms for nonexpansive linear systems[END_REF], with numerical analysis techniques, in particular for linear fixed point systems. A second attempt to design such kinds of PRNGs is to avoid the use of floating point arithmetic, by considering Boolean functions such as the logical negation, and then an asynchronous iteration scheme that includes this function coupled with a shift. In [START_REF] Guyeux | A topological study of chaotic iterations. application to hash functions[END_REF], the authors proved that such iterations, viewed as an operator on a suitable discrete dynamical system, satisfy various topological properties of chaos, like mixing.

Due to the finite structure of the problem it can be an hard task formalizing a notion of pseudorandomness for binary sequences, however Mauduit and Sárközy [START_REF] Maudit | On finite pseudoradom binary sequences i: The measure of pseudorandomness, the legendre symbol[END_REF] made an attempt in doing so. In this new framework, it is relatively easy to redefine some of the classical notions that arise from cryptography, such as collision and avalanche effect, as pointed out by Tóth in [START_REF] Tóth | Collision and avalanche effect in families of pseudorandom binary sequences[END_REF].

In this article we recall how to design a PRNG based on chaotic iterations, how to describe it as a dynamical system and, for the sake of completeness, we recall some of the topological properties of chaos previously obtained. Finally, we provide new investigations of such a dynamical system using the theoretical set proposed by Mauduit and Sárközy.

This research work is organized as follows. In the first two sections we recall all the basic facts, such as how to post treat an existing PRNG by using chaotic iterations. Then, in the last section, we show that the collision free property of a family of binary sequences is preserved if we post treat this family using chaotic iterations. This research work ends by a conclusion section, in which our contributions are summarized and intended future work is outlined.

II. BASIC RECALLS

We first explain how to properly design chaotic PRNGs by using discrete dynamical systems. To accomplish this task we must introduce some notations and terminologies.

In the remainder of this article, N is the set of natural (nonnegative) numbers, while N * stands for the positive integers 1, 2, 3, . . . s n denotes the n th term of a sequence s while X N is the set of all sequences whose elements belong to a given set X; also we denote by (x n ) n∈N a sequence of those sequences. Instead of working with sequences of length N in {0, 1} N , it will be useful considering sequences in {-1, +1} N . The set {-1, 1} is denoted 2.

Definition 1 Let N ∈ N * , f : 2 N -→ 2 N be a function, and s ∈ [1, N] N be a sequence of integers between 1 and N. The so-called chaotic iterations are defined by x 0 ∈ 2 N and

∀n ∈ N * , ∀i ∈ [1; N], x n i = x n-1 i if s n = i f (x n-1 ) sn if s n = i.
In other words, at the n th iteration, only the s n -th component of the vector x n is updated; let us explain it with an example.

Example 1 For the sake of concreteness, let us consider that N = 3. Let s be the sequence s = (123123123123...) and let f be the Boolean negation (where in this case the boolean negation of -1 is 1 and viceversa). If we start with input (111), the chaotic iterations will produce the following output

(1, 1, 1) → (-1, 1, 1) → (-1, -1, 1) → (-1, -1, -1) → → (1, -1, -1) → . . .
At each iteration step, we look at the correspondent element of the sequence s (at the n-th step, we look at the n-th element of s), and then we apply f to the element of the input sequence suggested by the element of s that we are considering. In the example, for instance, at step number 5 we need to change the second bit of the input, since the fifth element of s is 2.

We shall use the term "chaotic" (and similar) with various meaning throughout this paper and each of these terms has no link a priori with the other ones. In the definition of chaotic iterations, there is a sequence s of components to update as input, and a sequence of binary vectors as output. In other words, given a function f , these chaotic iterations transform a sequence s = (s n ) n∈N in another sequence (x n ) n∈N . If s is provided by a pseudorandom number generator, we thus have defined a way to modify the produced sequence, leading to a post-treatment on this generator.

Let N ∈ N * and f :

2 N → 2 N . Any sequence u ∈ [1, N]
N , provided by an inputted pseudorandom number generator, defines a "chaotic iterations based" PRNG, which is denoted by CIPRNG 1 f (u). It is defined by [START_REF] Bahi | A novel pseudo-random generator based on discrete chaotic iterations[END_REF], [START_REF] Bahi | Improving random number generators by chaotic iterations. Application in data hiding[END_REF]:

x 0 ∈ 2 N , ∀n ∈ N, ∀i ∈ [1, N], x n+1 i = f (x n ) i if i = u n , x n i else.
The outputted sequence produced by this

CIPRNG 1 f (u) gen- erator is (x n ) n∈N .
However for several reasons to be discuss later, it will be useful to consider another version of this generator, in which is allow to update more than one digit at each iteration step. If P ⊂ N * is a non empty and finite set of integers, any couple

(u, v) ∈ [1, N] N × P N defines another chaotic iterations based PRNG, which is denoted by CIPRNG 2 f (u, v) (defined in [9]). It is defined as follows:        x 0 ∈ 2 N ∀n ∈ N, ∀i ∈ [1, N], x n+1 i = f (x n ) i if i = u n x n i else ∀n ∈ N, y n = x v n .
The outputted sequence produced by this generator is (y n ) n∈N . In other words, the first inputted generator u is the sequence of components to update, while the second inputted generator v provides the number of iterates between two outputs of the CIPRNG 2 f (u, v) generator. Note that

CIPRNG 1 f (u) is equal to CIPRNG 2 f u, (1
) n∈N , where (1) n∈N is the sequence that is uniformly equal to 1. Pseudorandom number sequences have plenty of applications in computer science and in particular in cryptography; in [START_REF] Bahi | Efficient and cryptographically secure generation of chaotic pseudorandom numbers on gpu[END_REF], it has been shown how to use CIPRNG to produce a cryptographycally secure PRNG. In [START_REF] Guyeux | A topological study of chaotic iterations. application to hash functions[END_REF] the chaotic property of CIPRNG has been analyzed, and some of them are listed in the next section; also simulations show the good behavior of such PRNGs. An interesting particular situation has been analyzed: if we consider the boolean negation as function involved in CIPRNG (a good PRNG should be fast and the boolean negation is really cheap from a computational point of view), simulations shown that CIPRNG 2 f (u, v) really improved the chaotic behavior of the inputted generator, while CIPRNG 1 f (u) does not work as good as its "big brother". In an article that will appear, authors provided a theoretical explanation to this phenomenon by using ergodic theory; in the next sections (after some recalls) it will be shown that actually CIPRNG 2 f (u, v), where f is the boolean negation and under some hypothesis, satisfies a classical cryptography property, namely to be collision free, while there is no reason to conclude that such a property is satisfy also by CIPRNG 1 f (u). Thus we conclude that CIPRNG 2 f (u, v) can be successfully used to improve the cryptography properties of a family of pseudorandom binary sequences.

III. CHAOTIC BEHAVIOUR OF CIPRNG

In this section are listed results concerning the topological randomness of the presented class of generators. They appear, for instance, in [START_REF] Guyeux | A topological study of chaotic iterations. application to hash functions[END_REF]. Firstly we recall the theoretical framework that can emphasize such a behavior and how to adapt the definition of CIPRNG to fit this framework (i.e., CIPRNG can be described as a dynamical system on a suitable metric space). Even if the results listed here are about CIPRNG 1 f , they can be easily extended to CIPRNG 2 f .

A. A notion of randomness: Devaney's theory of chaos

Consider a topological space (X , τ ) and a continuous function f : X → X on (X , τ ). Definition 2 Let U, V be any pair of opens subsets of X . If there exists an integer k > 0 such that f k (U ) ∩ V = ∅, then the function f is said to be topologically transitive Definition 3 If the set of periodic points for f is dense in X (i.e. for any point x in X , any neighborhood of x contains at least one periodic point), then f is regular on (X , τ ); let us recall that an element x is a periodic point for f of period

n ∈ N, n > 1, if f n (x) = x.
Definition 4 (Devaney's formulation of chaos [START_REF] Devaney | An Introduction to Chaotic Dynamical Systems[END_REF])

The function f is chaotic on (X , τ ) if f is regular and topologically transitive.

Banks et al. have proven in [START_REF] Banks | On Devaney's definition of chaos[END_REF] that, when the topological space is a metric one (X , d), chaos implies sensitivity, defined below:

Definition 5
The function f has sensitive dependence on initial conditions if there exists δ > 0 such that, for any x ∈ X and any neighborhood V of x, there exist y ∈ V and n > 0 such that d (f n (x), f n (y)) > δ. δ is called the constant of sensitivity of f .

To show that CIPRNG 1 f (u) has a chaotic dependence regarding modifications on u is equivalent to prove that "chaotic iterations" are indeed chaotic according to the definition of Devaney recalled above.

We first recall how to define a suitable metric space where chaotic iterations are continuous. For further explanations, see, e.g., [START_REF] Bahi | Topological chaos and chaotic iterations, application to hash functions[END_REF]. Let δ be the discrete Boolean metric, δ(x, y) = 0 ⇔ x = y. Given a function f , define the function

F f : [1; N] × 2 N -→ 2 N by: (k, E) -→ E j .δ(k, j) + f (E) k .δ(k, j) j∈[1;N]
where + and . are the Boolean addition and product operations. Consider the phase space: X = [1; N] N × 2 N , and the map defined on X by:

G f (S, E) = (σ(S), F f (i(S), E)) , (1)
where σ is the shift function defined by σ :

(S n ) n∈N ∈ [1, N] N -→ (S n+1 ) n∈N ∈ [1, N] N and i is the initial function i : (S n ) n∈N ∈ [1, N] N -→ S 0 ∈ [1; N].
Then the chaotic iterations proposed in Definition 1 can be described by the following discrete dynamical system, whose topological chaos can now be studied:

X 0 ∈ X X k+1 = G f (X k ). (2)
To do so, a relevant distance between two points X = (S, E), Y = ( Š, Ě) ∈ X has been introduced in [START_REF] Bahi | Topological chaos and chaotic iterations, application to hash functions[END_REF] as follows: d(X, Y ) = d e (E, Ě) + d s (S, Š), where

           d e (E, Ě) = N k=1 δ(E k , Ěk ), d s (S, Š) = 9 N ∞ k=1 |S k -Šk | 10 k . (3)
It has been established in [START_REF] Bahi | Topological chaos and chaotic iterations, application to hash functions[END_REF] that, Proposition 1 G f is continuous in the metric space (X , d).

Let f be a map from 2 N to itself. The asynchronous iteration graph associated with f is the directed graph Γ(f ) defined by: the set of vertices is 2 N ; for all x ∈ 2 N and i ∈ [1; N], the graph Γ(f ) contains an arc from x to F f (i, x).

It has been proven in [START_REF] Bahi | On the link between strongly connected iteration graphs and chaotic boolean discrete-time dynamical systems[END_REF] that, Additionally to the Devaney's chaos, a discrete dynamical system can be intrinsically complicated for various other understanding of this wish, that are not equivalent one another, like:

Theorem 1 Let f : 2 N → 2 N . G f is
• Undecomposable: it is not the union of two nonempty closed subsets that are positively invariant (f (A) ⊂ A). • Total transitivity: ∀n 1, the composition function

f n = f • f • . . . • f is transitive. • Strong transitivity: ∀x, y ∈ X , ∀r > 0, ∃z ∈ B(x, r),
∃n ∈ N, f n (z) = y. • Topological mixing: for all pairs of disjoint open nonempty sets U and V , there exists n 0 ∈ N such that ∀n n 0 , f n (U ) ∩ V = ∅. These varieties of definitions lead to various notions of chaos. For instance, a dynamical system is chaotic according to Wiggins if it is transitive and sensible to the initial conditions. According to Knudsen if a dynamical system has a dense orbit while being sensible then we are able to use again the adjective chaotic for such a system. Also, when the properties of transitivity, regularity, and expansiveness are satisfied we fit the so called expansive chaos.

Let us denote by C the set of f : 2 N → 2 N such that Γ(f ) is strongly connected. We have previously established that ( [START_REF] Bahi | Evaluating quality of chaotic pseudo-random generators. application to information hiding[END_REF], [START_REF] Guyeux | A topological study of chaotic iterations. application to hash functions[END_REF]) ∀f ∈ C, G f is strongly transitive undecomposable, unstable, and chaotic as defined by Wiggins.

Further results have been obtained when considering the boolean negation f 0 , which are summarized below [START_REF] Guyeux | Chaotic iterations versus spreadspectrum: chaos and stego security[END_REF].

Theorem 2 (X , G f0 ) is topologically mixing and expansive with a constant equal to 1.

Definition 6 Let (X , d) a metric space and f : X -→ X a continuous map. A scrambled couple of points is a pair (x, y) ∈ X 2 such that lim inf n→∞ d(f n (x), f n (y)) = 0 and lim sup n→∞ d(f n (x), f n (y)) > 0, or in other words, the two orbits oscillate.

A scrambled set is a set in which any couple of points are a scrambled couple, whereas a Li-Yorke chaotic system is a system possessing an uncountable scrambled set.

We have previously stated that [START_REF] Guyeux | Le désordre des itérations chaotiques et leur utilité en sécurité informatique[END_REF], Theorem 3 Chaotic iterations are chaotic as defined by Li and Yorke.

Theorem 4 Chaotic iterations using the boolean negation have a topological entropy equal to ln(N).

Theorem 5 Chaotic iterations using the boolean negation have an exponent of Lyapunov equal to ln(N).

Finally, it has been proven in [START_REF] Bahi | On the link between strongly connected iteration graphs and chaotic boolean discrete-time dynamical systems[END_REF] that, Theorem 6 Let f : 2 n → 2 n , Γ(f ) its iteration graph, M its adjacency matrix and M a n × n matrix defined by

M ij = 1 n Mij if i = j and M ii = 1 -1 n n j=1,j =i Mij otherwise.
If Γ(f ) is strongly connected, then the output of CIPRNG 1 f (u) follows a law that tends to the uniform distribution if and only if M is a double stochastic matrix.

These results of topological chaos and uniform distribution have initially led us to study the possibility of building a pseudorandom number generator (PRNG) based on chaotic iterations. As G f , defined on the domain [1; N] N × 2 N , is built from Boolean networks f : 2 N → 2 N , we can preserve the theoretical properties on G f during implementations (due to the discrete nature of f ).

IV. COLLISION PROPERTY

A. Basic recalls

It should be clear from the last section that there are several ways to define the notion of pseudorandomness. Recently a constructive approach to this notion was initiated by Mauduit and Sárközy [START_REF] Maudit | On finite pseudoradom binary sequences i: The measure of pseudorandomness, the legendre symbol[END_REF], and during these years it has been extended to more general situations. V.Toth in [START_REF] Tóth | Collision and avalanche effect in families of pseudorandom binary sequences[END_REF] adapted the classical cryptographical notions of collision and avalanche effect to this new framework. A good survey of this new approach can be found in [START_REF] Sárk | On pseudorandomness of families of binary sequences[END_REF]; in this section we will recall the basic definition and in the next section we will discuss the behaviour of CIPRNG in this new framework, with particular focus on the collision property. Let N ∈ N, let S be a given set (e.g., a set of certain polynomials or the set of all binary sequences of a given length much less than N ), to each s ∈ S we assign a unique binary sequence E N (s) = (e 1 , . . . , e N ) ∈ {-1, +1} N , and let F = F(S) denote the family of the binary sequences obtained in this way:

F = F(S) = {E N (s) : s ∈ S}. Definition 7 We have a collision in F if there are two different elements s = s ∈ S such that E N (s) = E N (s );
if there is no collision in F then F is said to be collision free.

An ideally good family of pseudorandom binary sequence should be collision free or, at least, the number of collision should be limited. We can reformulate this notion. Definition 8 Let N ∈ N and let E N = (e 1 , . . . , e N ), E N = (e 1 , . . . , e N ) be two sequences of {-1, +1} N , then the distance between E N and E N is defined by

d(E N , E N ) = |{n : 1 ≤ n ≤ N, e n = e n }|.
Moreover the distance minimum m(F) of F can be defined as m(F) = min

s =s , s,s ∈S d(E N (s), E N (s )).
With the definitions above easily follows that a family F is collision free if and only if m(F) > 0.

B. Collisions in CIPRNG 1

Let F(S) be a family of binary sequences generated, for instance, by a PRNG (or several PRNGs); we would like to post process this family with CIP RN G and see when we are able to improve the distance minimum of the inputted family.

Firstly we investigate a very simple case, i.e. the effects of CIP RN G 1 f where f is the boolean negation; with abuse of notation (not really) let S = F(S) be a set of strategies of length k, let u = 1 N be the input of CIP RN G 1 f and let F = {CIP RN G 1 f (s) : s ∈ S} the outputted family of binary sequences. We want to study collisions in F , i.e. we want to see what kind of conditions we need in order to get m(F ) > 0.

For better notations we shall write C(s) instead of CIP RN G 1 f,N (u, s), since the input string u and the length N are fixed. Now observe that if C(s) = (c 1 (s), . . . , c N (s)), then for n = 1, 2, . . . , N we have

c n (s)c n (t) = +1 if c n (s) = c n (t), -1 if c n (s) = c n (t) thus 1 2 (1 -c n (s)c n (t)) = 0 if c n (s) = c n (t), +1 if c n (s) = c n (t) It follows that d(C(s), C(t)) = N n=1 1 2 (1 -c n (s)c n (t)) = 1 2 (N - N n=1 c n (s)c n (t)) = 1 2 (N - N n=1 (-1) |n∈s|+|n∈t|(2) ),
where |n ∈ s| is the number of occurrences of the symbols n in the string s. To understand the last equality in the expression above observe that the n-th digit of c(s) (i.e. c n (s)) is equal to 1 if n appear an even numbers of times in the strategy (since the boolean negation is such that f • f = id and we started from input u = 1 N ), otherwise is equal to -1.

So far from now we obtained that if we want to estimate

m(F ) = min s =t d(C(s), C(t)),
we should be able to estimate the following quantity

max s =t N n=1
(-1) |n∈s|+|n∈t| (2) .

However, since the strategies are computed in some way by some PRNG (or several PRNGs), is not clear why assuming the collision free property of the set of the strategies should led to the collision free property for the outputted family, that now can be expressed by the following formula

max s =t N n=1 (-1) |n∈s|+|n∈t|(2) < N.

C. Collisions in CIPRNG 2

In the last section it has been shown that there is no reasons to conclude that CIPRNG 1 preserves the collision free property; in this section it will be show that CIPRNG 2 is more suitable for this task. As in the previous section, let S = F(S) be a set of strategies generated in some way and let P = F (P) be a set of iterations strings (pseudo) randomly generated (by some PRNGs). Recall that the strategies are sequences in {1, . . . , N } k while the iteration strings are sequences in {2, 3} k1 ; the set {2, 3} can be actually any finite set of integers with k 1 < k.

We fix the input u = 1 N and we want to post process the family S × P with CIP RN G 2 f (s, p) = C(s, p), where f is the boolean negation and (s, p) ∈ S × P .

As before we want an estimate of this quantity

d(C(s, p), C(t, q)) = 1 2 (N - N n=1 c n (s, p)c n (t, q));
recalling the structure of CIP RN G 2 , we can observe that an iteration string p induces a partition of the strategy s (in block of size 2 or 3 in this case), thus we can write s = s 1 s 2 . . . s r , for some r, where |s i | = p i for i = 1, . . . , r -1 and p r ≤ |s r | ≤ p r +max{2, 3}-1. To avoid the fact that the last block of a strategy can have size different from 2 or 3, we will just consider r -1 iterations instead of r. Due to the fact that we are using the boolean negation and we started with input 1 N , we have

d(C(s, p), C(t, q)) = N n=1 1 2 (1 -c n (s, p)c n (t, q)) = 1 2 (N - N n=1 c n (s, p)c n (t, q)) = 1 2 (N - N n=1 (-1) |n∈s|+|n∈t|(2) ), = 1 2 (N - N n=1 (-1) r-1 j=1 |n∈sj |+|n∈tj |(2) )
So now our aim is to prove the following inequality

max (s,p) =(t,q) N n=1 (-1) r-1 j=1 |n∈sj |+|n∈tj |(2) < N.
Suppose N prime (in this way things works much better, but probably this is not necessary) and suppose for an absurdum that there are (s, p) = (t, q) such that We will show now that, under the assumption that the family of the iteration strings P is collision free, then there is some j * such that N n=1 (-1) 2) is neither 1 nor N , thus we reach an absurdum since N is prime.

|n∈s j * |+|n∈t j * |(
In fact, if the family P is collision free, then for any two elements, there is at least one index with different digits in those two elements. Thus, they induced a different partition of the strategies, i.e., there is a j * such that (without loss of generality) |s j * | = 2 and |t j * | = 3. There are two cases, either s j * = AA or s j * = AB, where A, B ∈ {1, . . . , N }; in each of these cases, we can compute all the possible configurations of t j * and the correspondent values. case 1 In this case the possible configurations of t j * are t j * = BCD in this case A appears with the same parity in both s and t, while B, C, D appear with different parity, and then there are other N -4 letters that appear with the same parity, thus the total is 1 + (-3) + N -4 = N -6; t j * = ABC (-3) + (N -3) = N -6; t j * = ABB (-1) + (1) + N -2 = N -2; t j * = AAB (1) + (-1) + N -2 = N -2; t j * = AAA (-1) + (N -1) = N -2.

and this prove the first case. case 2 For the second case, we will first consider all the possible configuration of t j * without any letter in common with s j * . Then, all the possible configuration with one letter in common, and finally the case with two letters in common. t j * = CDE (-5) + (N -5) = N -10; t j * = ACD observe that we can assume that the letter in common is A; here we have (1)+(-1)+(-2)+ N -4 = N -6; t j * = ACC (1) + (-1) + (N -2) = N -2; t j * = AAC (-1) + (-1) + (-1) + N -3 = N -6; t j * = ABC (+2) + (-1) + N -3 = N -2; t j * = ABA (1) + (-1) + N -2 = N -2.

V. CONCLUSION AND FUTURE WORK

Designing a pseudorandom number generator is a complicated task. In this article we presented a way to post treat a family of binary sequences in order to improve their randomness. After having recalled the structure of such a post treatment and its properties, we proved in the new framework introduced by Mauduit and Sárközy, that the post treatment preserves the collision free property of the family, a classical notion in cryptography.

In future work we intend to extend our investigation to the avalanche effect property, and more in general we would like to proof that such a post treatment actually increase the distance minimum of a family of binary sequences. Seems reasonable to think that the structure of the chaotic iterations strictly depend from the boolean function involved in the process, we would like to explore a possible characterization of such a property.

ACKNOWLEDGEMENT

I would like to say thank you to professor Stefano Galatolo for his guide through this new world (and for saying me where to put his name in this article).


Zaher Al Chami
email: zaher.al-chami@utbm.fr

Hervé Manier
email: herve.manier@utbm.fr

Marie-Ange Manier
email: marie-ange.manier@utbm.fr

Différentes approches pour la résolution d'un problème multipériodes de collecte et de livraison

Keywords: PDPTW, Multi-périodes, Méthodes exactes, Méthodes approchées

Introduction :

Le problème de collecte et de livraison PDP (Pickup and Delivery Problem) consiste à planifier un ensemble des tournées de véhicules, initialement basés au dépôt, afin de satisfaire toutes les demandes de transport des marchandises entre les points de collecte et ceux de livraison.

Dans cette étude, nous nous intéressons à la résolution d'une variante avec fenêtres temporelles dans lesquelles les visites doivent être accomplies. De plus, les demandes sont appairées ce qui donne la possibilité de transporter simultanément différents types de produits. Ainsi, chaque demande lie une origine à sa destination précise.

Chaque origine peut être un site fournisseur ou le dépôt et chaque destination peut être associé à un site client ou encore au dépôt. Dans un contexte urbain, cela génère trois types de flux : le premier relie le dépôt et un client (flux entrant), le deuxième relie un fournisseur et le dépôt (flux sortant) et le dernier fait le lien entre un fournisseur et un client (flux intra-urbain). La satisfaction de toutes les demandes de transport sur une période donnée est parfois impossible à cause de la capacité limitée des véhicules, de l'existence de fenêtres temporelles et des contraintes de précédence. Dans ce cas, il y aura des sites nonvisités pendant cette période. L'enjeu de notre travail est de résoudre une variante multi-périodes dans laquelle chaque site peut accepter d'être servi pendant une ou plusieurs périodes. L'objectif est de satisfaire toutes les demandes en minimisant la distance totale parcourue.

Nous avons développé quatre approches différentes pour résoudre la variante « Mu-PDPTWPD » (Multi-periods PDP with Time Windows and Paired Demands). La première méthode est de traiter la globalité du problème afin d'obtenir des solutions optimales. Par contre, les trois autres approches consistent à résoudre le problème multi-périodes en tant qu'une séquence de Mo-SPDPTWPDs (Monoperiod Selective PDPTWPDs). L'aspect sélectif est traduit par le fait que la visite de tous les sites n'est pas obligatoire.

Méthodes de résolution : 2.1 MILP global :

Dans une étude précédente, nous avons élaboré le premier modèle mathématique qui résout la variante « Mu-PDPTWPD » [START_REF] Manier | Shippers' collaboration in city logistics[END_REF]. Cette méthode exacte nous donne les solutions optimales mais son inconvénient est le temps de résolution très long.

MILP séquentiel :

Cette deuxième méthode consiste à utiliser une méthode exacte pour la résolution du Mo-SPDPTWPD introduite dans [START_REF] Al Chami | New model for a variant of pick up and delivery problem[END_REF] basée sur un MILP. Cette méthode sera lancée autant de fois que le nombre de périodes. L'objectif de ce MILP est de maximiser le profit collecté puisque dans un problème sélectif, un profit est associé à chaque site. Du coup, nous pouvons dire que cette méthode peut nous donner une solution faisable rapidement mais un peu loin de l'optimale puisque la minimisation de la distance n'est pas un critère dans le choix de la solution finale.

Approche lexicographique séquentielle :

Pour pallier le point négatif de l'approche précédente, nous avons décidé d'utiliser une approche lexicographique d'une manière séquentielle. Cette approche consiste à maximiser le profit puis à minimiser la distance en fixant la valeur du profit à sa valeur optimale [START_REF] Al Chami | A lexicographic approach for the bi-objective selective pickup and delivery problem with time windows and paired demands[END_REF]. Donc, le MILP sera lancé deux fois par période. Cette méthode peut générer de très bonnes solutions mais le temps de résolution nécessaire pour obtenir de telles solutions représente un inconvénient.

HGA séquentiel :

Compte tenu des inconvénients de deux premières méthodes séquentielles, nous avons développé une troisième approche basée sur l'algorithme génétique hybride introduit dans [START_REF] Al Chami | A hybrid genetic algorithm to solve a multi-objective Pickup and Delivery Problem[END_REF]. Cette méthode nous a donné de bonnes solutions en un temps très court.

Résultats :

A notre connaissance, il n'y a pas d'instances dans la littérature qui traitent exactement notre variante. Pour cela, nous avons généré nos propres instances afin de tester nos méthodes de résolution. Le tableau 1 présente les résultats des différentes approches sur quelques instances. Nous pouvons remarquer que « MILP global » fournit les solutions optimales. En revanche, cette méthode n'arrive à trouver les solutions optimales en un temps raisonnable (2 heures) que pour 36 instances sur 48 générées. L'approche « MILP séquentiel » génère des solutions avec un gap moyen de 21.60% par rapport à la distance optimale et un gain moyen de 58.19% au niveau du temps de résolution. En plus, la méthode « Lexicographique séquentiel » donne des bonnes solutions avec un gap moyen de 13.4% par rapport à la distance mais avec un gain de temps de résolution moyen de 26.77% seulement. Ces deux premières méthodes séquentielles ont permis de trouver des solutions pour 44 instances en un temps raisonnable. Avec la dernière approche « HGA séquentiel », nous trouvons des solutions pour les 48 instances avec un gain moyen en temps de 66.12% et avec un gap moyen en distance de 13.74% par rapport à la méthode globale.

Conclusion et perspectives :

Dans cette étude, nous avons élaboré une nouvelle façon d'aborder le Mu-PDPTWPD en le considérant comme une séquence de Mo-SPDPTWPDs. De plus, nous avons généré les premières instances de cette variante dans la littérature. Plusieurs tests ont été menés afin de valider les avantages et les inconvénients d'une telle décomposition. Pour la suite, nous allons travailler à développer une méthode approchée qui peut résoudre notre variante dans sa globalité dans un temps de résolution acceptable.

Ce travail est soutenu par l'ANR (Agence Nationale de la Recherche) dans le cadre du projet TCDU (Transport Collaboratif dans la Distribution Urbaine). Ce projet ANR-14-CE22-0017 est labellisé par le