Imade Koutiri
email: imade.koutiri@ensam.eu

Etienne Pessard

Patrice Peyre

Ouafae Amlou

Thibaut De Terris

Influence of SLM process parameters on the surface finish, porosity rate and fatigue behavior of as-built Inconel 625 parts

Keywords:

This paper is dedicated to understanding fatigue crack initiation for an Inconel 625 manufactured by SLM, using a hatching + contour procedure. In the first part of the paper, an optimum set of parameters was found to deliver the best surface roughness combined with low porosity. This process optimization, mostly focused on adjusting the volume energy density aimed at finding a compromise between an optimum densification state and a minimum number of contaminating spatters. Secondly, a fatigue test campaign has been conducted on as-built SLM samples or polished samples. The analysis of failure surfaces allowed identifying different heterogeneities at the origin of the fatigue damage for each configuration. The embedded particles on the surface of as-build specimens have been shown to play an important role in fatigue and need to be optimized or taken into account in the fatigue strength design of SLM components.

Introduction

Additive Layer Manufacturing (ALM) currently impacts several industrial sectors, ranging from medical to aerospace and energy. Since 2003, there has been a 42% increase of parts produced by additive manufacturing. ALM processes, generating mechanical parts layer by layer, allow realizing complex shapes which are not achievable by traditional processes at significantly reduced cost. A limitation for the development of this process is still considered to be the mechanical properties obtained of as-built ALM parts, and the possible occurrence of defects, allowing the direct use-or not-of these parts in an industrial environment.

Static properties of ALM materials have been the subject of a tremendous number of works in the past 3 years, and it has been shown that ALM alloys could exhibit correct static mechanical resistance compared to wrought or cast alloys. For instance, [START_REF] Chastanda | Fatigue characterisation of titanium Ti-6Al-4V samples produced by additive manufacturing[END_REF] show that fatigue properties of Titanium Ti6Al4V parts built by SLM are comparable to casting after stress relieving and comparable to wrought processes after an HIP heat treatment.

The durability of ALM parts has become an important topic to address for many applications, fatigue failure of materials being one of the most common damage modes in cyclically loaded metallic materials. Many works have resulted in a better knowledge of the phenomena involved. These studies contributed to a better understanding of the mechanisms of crack initiation involved on ALM materials and proposed new models to predict the material fatigue behavior.

In many papers, it appears that high cycle fatigue cracks, for part with as-built surface, usually initiate on:

-porosities on the surface or sub-surface, generated by ALM process for different materials like Ti6AL4V [START_REF] Ellyson | Characterisation of bendig vibration fatigue of SLM fabricated Ti6Al4V[END_REF], AlSi10Mg [START_REF] Brandl | Additive manufactured AlSi10Mg samples using selective laser melting (SLM): microstructure, high cycle fatigue and fracture behavior[END_REF], AlSI12 for very high cycle fatigue [START_REF] Siddique | Very high cycle fatigue and fatigue crack propagation behavior of selective laser melted AlSi12 alloy[END_REF] and 630 stainless steel [START_REF] Akita | Defect-dominated fatigue behavior in type 630 stainless steel fabricated by selective laser melting[END_REF]. Such subsurface porosities are mostly located in the transition zone between contour path and hatching path. -Elevated surface roughness, especially on inclined parts, that tends to generate stress concentration factors and favor crack's initiation. -Un-melted particles located close to the surface of specimens. Such particles have been observed already on a Ti6AL4V [START_REF] Fatemi | Torsional fatigue behavior of wrought and additive manufactured Ti6Al4V by powder bed fusion including surface finish effect[END_REF] and on a AlSi10Mg [START_REF] Brandl | Additive manufactured AlSi10Mg samples using selective laser melting (SLM): microstructure, high cycle fatigue and fracture behavior[END_REF].

Following this last statement, one of the well-known issues in SLM is the contamination of the powder bed surface by large spatters [START_REF] Khairallah | Laser powder bed fusion aditive manufacturing : physics of complex melt flow and formation mechanisms of pores, spatter and denudation zones[END_REF] generated during laserpowder bedmeltpool interaction. The generation of such spatters is usually favored by high volume energy densities VED = P/ (S. V) (high power P, low scan velocity V, small laser-irradiated surface S) that tend to destabilize melt-pools and surrounding powder beds. However, it is also well known that high VED values are required to obtain an optimum densification state. For this reason, and because a compromise is hard to find, most of SLM parts are systematically submitted to post-machining steps (sand blasting, electro-polishing, grinding) to improve surface finish. Such machining steps are not easy to carry out on complex shapes, and in some specific cases, as-built SLM samples can be envisaged for parts submitted to reduced static or cyclic loading. Consequently, mechanical investigations on as-built samples are still of importance.

The current study focuses mainly on the effects of material's heterogeneities (surface roughness and porosity) generated in an Inconel 625 part by the SLM process. In a first step, process parameters such as laser power P, beam diameter D, scanning speed V and orientation angle θ of the piece versus the building platform will be varied in order to optimize surface roughness and porosity rate (considered as the percentage of void spaces in a material over the total volume). In a second step, as-built specimens manufactured with optimum SLM conditions will be submitted to cyclic loading in order to evaluate fatigue limits and to detect the origin of the crack initiation.

Experimental procedure: Inconel 625 and SLM conditions

Selective Laser Melting (SLM) is one of the two powder bed method (with EBM) for the additive manufacturing of metallic alloys. It allows the processing of a large range of metals including Ni-based alloys, stainless steels, cobalt-chromium, aluminum and titanium alloys.

Inconel 625

A wide range of Ni-based superalloys, like Inconel 718, Nimonic 623 has been manufactured already with additive layer process, including Inconel 625. This nickel based superalloy is widely used in applications like supercritical water reactors or engine with high temperatures. It combines high mechanical properties with good corrosion resistance due to its high Cr content. Nickel is the principal element of this alloy and the principal constituent of the austenitic matrix γ. Ni participates also in the composition of secondary phases (γ′, γ″ and δ) and gives to the material a good resistance on a large range of temperature (from cryogenic temperatures to 900 °C).

The following table summarizes literatures data for mechanical properties of Inconel 625 in a wrought or SLM state (for horizontal and vertical direction). Our results and literature data show a variability in the mechanical properties (especially in elongation rate) depending on the set of process parameters and the building direction.

The Inconel 625 powder, obtained by inert gas atomization, was supplied by PRAXAIR with a 10-45 μm particle range. This particle distribution was confirmed by a laser granulometer device (CILAS) as shown in Fig. 1. The particle shape has also been observed by Scanning Electron Microscopy (SEM). Particles present a mostly spherical shape with no satellites (Fig. 1), and no agglomeration.

Prior to SLM, the powder was systematically heated for 24 h at 50 °C in order to remove the residual moisture adsorbed on particles. This heating stage allows improving the powder spreading in order to have homogeneous powder layers.

SLM conditions

All the specimens investigated in this study are manufactured by a SLM 125HL machine provided by the SLM Solutions GmbH company. In the first part of the study, we will vary the irradiation parameters: laser power P (between 100 and 400 W), scanning velocity V (0.5-1.8 m/s), beam diameter D (70-160 μm), for both the contour and hatching parameters (Table 2). A specific focus will also be given to the influence of tilt angle θ (0°to 35°) between the platform and the part. The other process parameters (hatch distance, hatching strategy = stripes …) will be fixed and will not be considered. The objective of this optimization step is to obtain parts with the best surface roughness and the lower porosity ratio, which are both known to promote fatigue crack initiation. In this first stage, dedicated to the parametric study, only blocks of 10 × 18 × 4 mm were manufactured. For each set of parameters, a volume energy density (VED) has been calculated using Eq. (1). Such an equation allows considering the amount of energy induced by the laser scanning on the powder bed as an unified process parameter. This specific analytical formulation of laser energy input VED can be directly correlated to penetration depths of beads as shown analytically by [START_REF] Fabbro | Analysis and possible estimation of keyhole depths evolution using laser operating parameters and materials properties[END_REF] and confirmed experimentally in [START_REF] Gunenthiram | Experimental analysis of spatter generation and melt-pool behavior during the powder bed laser beam melting process[END_REF].

= VED J m P W V m s A m ( / ) ( ) / ( / ) . ( ²) 3 (1)
With A = spot area (π. D²/4) All the process parameters (Contour and hatching parameters) are reported in Table 2.

The second part of the study is dedicated to the production of tensile and fatigue samples (Fig. 2), using optimized process parameters (P, V, D). The geometry of tensile and fatigue specimens are shown in Fig. 3. A rounded edge was carried out to avoid local fatigue crack initiation.

3. Influence of SLM parameters on surface finish 3.1. Influence of laser power (contour stage) and building angle All the samples were characterized using a VEECO Dektak 150 profilometer, using a 12 μm curvature radius stylus and a 1 mg weight to scan the surfaces. The measurement technique is based on contact stylus profilometry. The non-filtered surface finish values Sa presented here below (Table 2) are the mean values of three scan surfaces of 1 × 1 mm.

We first focus here on the effect of the tilt angle and contour laser power on the surface roughness, at constant beam diameter (70 μm) (samples 1-6). Even if the hatching stage is also known to affect surface roughness by varying the average temperature of parts, only contour parameters have been addressed whereas hatching parameters have been kept constant (P = 175 W and V = 0.6 m/s). Two tilt angles θ (10°a nd 35°) and three laser powers (100, 140 and 180 W) have been tested. 35°was selected as the threshold angle above which supporting is usually required in SLM, due to heat accumulation in bottom sides (because of the low conductivity of powder bed). The scanning velocity is fixed and equal to 0.6 m/s. Finally, six blocks (specimen 1-6), presented in Table 2 have been obtained. For specimen 1-6, the VED values are related only to contour parameters contrary to the other specimens.

Results are presented in Figs. 4 and 5 considering Sa average roughness, which measures the arithmetic average deviation of the nonfiltered measured profile from the centerline of the measured profile. Globally, S a values are comprised between 6 and 11 μm and tend to decrease with heat input (=laser power), indicating a smoothing effect, except for the 35°angle, for which specific boundary conditions (lower heat dissipation) tend to provoke a roughness increase between 140 W and 180 W (Fig. 4). This result confirms afterwards that supporting is mandatory above θ = 35°to favour heat dissipation and limit deleterious surface finish.

As shown in Fig. 5, even for low building angles (θ = 10°) a roughness difference is systematically found between top and bottom sides of the specimen. [START_REF] Leary | Surface roughness optimisation for selective laser melting (SLM): accommodating relevant and irrelevant surfaces[END_REF] and [START_REF] Triantaphyllou | Surface texture measurement for additive manufacturing[END_REF] observed the same behaviour in a Ti6Al4V. This indicates different temperature ranges between up and bottom sides, whatever the building angle. These first results confirm that both the angle between the platform and the part and the laser power influence the surface roughness, which can be modified by a factor 2 depending on process parameters. To explain these results, it can be assumed that the influence of laser power is two-fold: on the one hand, higher powers tend to favour layer re-melting and provoke a smoothing effect, and on the other hand excessive laser powers, combined with large angles θ favour sticking of surrounding particles on overhanging parts in the powder bed (Fig. 6).

3D scans corresponding to six tests presented in Fig. 4 confirm the results already presented.

-For θ= 10°, a decrease of the surface ratio of large "red" agglomerates on the surface with laser power is obtained. Such agglomerates disappear for P = 180 W. -For θ = 35°, large agglomerates decrease between 100 and 140 W then increase-again for P = 180 W.

This confirms that an increase of laser power tends to increase the size of the melted area, allowing a better re-melting of previous layers and a better spreading of large particles stuck and incrusted onto the surface (θ = 10°). Therefore, the quantity of agglomerated powder decreases, which improves the surface roughness. Such macro-agglomerates can mostly be attributed to the combination of two factors:

(1) liquid spatters formed during laser-powder-melt-pool interaction, which tend to increase with VED values and contaminate powder beds and manufactured parts, or (2) to the surrounding powder bed where powder grains located in the close vicinity of the fused walls can coalesce, grow and stick onto the surface. Those incrusted particles are susceptible to play a dominant role during the cyclic tests and favor fatigue crack initiations as they are fully embedded in the surface as shown in Fig. 7.

Influence of contour scan speed and beam diameter

In a second stage, we addressed the effect of scanning velocity and focus diameter of the contour phases on the resulting surface roughness considering a unique building angle θ = 0°. The objective was to check if the use of a single and unified energy parameter, including P, V, and D was sufficient to represent roughness variations. For these specific trials (samples 7-16 in Table 2), the laser power was increased and set constant to 400 W, the maximum available laser power of our SLM machine. Three beam diameters (70, 160 and 240 μm), plus four scanning velocities (1, 1.4, 1.8 and 2 m/s) were tested. With such conditions, the VED range was larger than previously used with a 100-180 W laser power, and similar contour / hatching process parameters were used.

At P = 400 W, one can observe (Fig. 8): Even if a clear smoothing effect is globally obtained with an increase of VED, if we consider similar building angles (≈0-10°), beam diameters (70 μm) and VED values (75-100 J/mm 3 ), better surface finishes are obtained for P = 140-180 W than for P = 400 W (samples 2,3 with Ra 6-7 μm versus samples 7,8 with Ra = 10-12 μm in Table 2). This is assumed to be due to an increase of spattering effect for P = 400 W.

The clear difference evidenced between surface finishes obtained with similar VED (for specimens 11-16 with VED between 6-20 J/ mm3) but different beam diameters mostly indicate that VED values are not sufficient to represent correctly the influence of SLM parameters on surface finish. Sa results presented here above mostly traduce the influence of scan speed at constant laser power: lower scan speed promote large but turbulent melt-pools which deteriorate surface finish, by increasing the amount of agglomerated spatters particles on the built surfaces [START_REF] Khairallah | Laser powder bed fusion aditive manufacturing : physics of complex melt flow and formation mechanisms of pores, spatter and denudation zones[END_REF].

Finally, such results clearly show that the correlation between surface roughness, which depends mostly on surface temperatures achieved on SLM part, and process parameters, is not a simple task, and that the use of VED as a representative factor is not fully satisfactory, when considering a large range of laser spot diameters. This could be attributed to a balance between the density of spatters contaminating SLM surfaces, which tend to increase with smaller scan speeds (at constant beam diameters) and to the global surface smoothing effect with smaller spot-diameters.

Influence of hatching process parameters on porosity rate

The second part of this study addresses the variation of porosity rate with hatching parameters, with the objective to minimize the void fraction, porosity being one of the microstructural heterogeneity which can initiate a fatigue crack. Various hatching parameters were used (Laser power, scanning velocity and angle of building), while keeping constant contour parameters (P = 140 W, V = 0.5 m/s).

According to previous works [START_REF] Spierings | Comparison of density measurement techniques for additive manufactured metallic parts[END_REF], different techniques can be used for measuring densities on additive manufactured metallic parts: Archimedes' method, microscopic analysis of cross-sections and X-ray scanning, with up to 1% gap at elevated porosity rate.

Provided that wetting of metallic surfaces in the selected liquid is optimized, the Archimedes method is supposed to deliver more reliable and more global results as the whole sample volume is taken into account instead of some single cross sections. Moreover, this method is very easy to use, fast and cheap. Based on these results, density measurements have been carried out with the Archimedes' double weighing method to measure the density of parts. The basic principle consists of measuring the mass in normal air and in the fluid (water in the current paper). In the present configuration, the measuring accuracy of the Mettler Toledo MS104TS balance is 0.0001 g.

The following equation allows the calculation of the density ρ p of the part:

= - -+ ρ m m m ρ ρ ρ . ( ) p a a f l fl a a (2)
With m a = mass in the air, m fl mass in the fluid, ρ a = air density, ρ fl = fluid density In the current study, we have investigated the possible link between porosity rate and volume energy density. Eight conditions have been tested including various hatching parameters (Laser power, scan velocity) (Table 3).

Two building angles (10°and 35°) have also been tested as a possible contributor to densification rate. A "stripe" scanning strategy has been chosen for the building of the different samples, with a constant hatch distance of 140 μm. For each specimen, three measurements have been performed. Results are presented in Table 3. For a constant building angle (0°), the results (Fig. 9) corresponding to specimens 7-21 indicate that the mean porosity rate value globally decreases with VED values (J/mm 3 ) above a VED threshold (≈7 J/ mm 3 ). The main explanation is that high VED promote larger meltpools susceptible to better re-melt inter-layer porosities. For 3 different beam diameters, and various (P,V) values a similar evolution is shown with VED. Moreover, for all the (P, V, D) conditions considered, a correct densification (to less than 0.5% porosity rate) is shown above a similar VED threshold (close to 30-40 J/mm 3 ).

The inclination of parts (from 10 to 35°) has a clear effect on densification rates because higher inclinations increase average temperatures due to heat concentration in the lower sides of SLM walls (Fig. 10) for specimens 22-29. This explains why, in % porosities versus VED curves the densification starts for 2 times lower VED on 35°inclined parts (Fig. 10).

All these data globally tend to confirm the possibility to use VED as a rather reliable process parameter combining P, V and D, for a constant building angle θ and layer height Δh.

However, the variation of laser beam spatial distribution with zdefocusing may explain why the use of VED is not fully sufficient to predict porosity rate with a good accuracy for a large range of diameters (variations of the porosity rates are shown).

Finally, for optimized conditions, a 0.1-0.25% minimum porosity rate is achievable with SLM, to be compared with the values obtained by [START_REF] Spierings | Comparison of density measurement techniques for additive manufactured metallic parts[END_REF] with SLM on a 316L stainless steel, which are comprised between 0.3 and 0.5%. [START_REF] Tang | Prediction of lack-of fusion porosity for powder bed fusion[END_REF] has summarized different values of porosity rate for different material, which are comprised between 0.8 and 3%, in order to propose a geometry based simulation to predict porosity caused by insufficent overlap of melt pools in powder bed fusion.

An interesting point to notice on cross-sections micrographs (Fig. 11) is the change of size and geometry of pores in the volume and close to the surface. Two different geometries of pores are clearly shown on cross-sections:

-In the sample volume, the geometry of pores is circular and their size is between 10 to 40 μm (0.3% of porosity). -Close to the surface of the specimen, their sphericity decreases and their mean diameter increases up to 150 and 200 μm (2.2% of porosity). This last family of porosity can be explained by the overlap between the two strategies, (contour + hatching), or by decelerating-accelerating issues of the scan head near the contour, which favor deeper beads, and possible key-hole closure effects, similar to welding issues [START_REF] Sames | The metallurgy and processing science of metal additive manufacturing[END_REF]. [START_REF] Kasperovich | Correlation between porosity and processing parameters in TiAl6V4 produced by selective laser melting[END_REF] has done a similar study in a Ti6AlV4 produced by SLM with the effect of more parameters. The scan-strategy used is the zigzag scan vectors. They observed that the velocity exerts the most dominant influential parameter (with porosity between 0,1% at 125 J/mm² and 3,2% at 300 J/mm²) as obtained in our study.

Fatigue behavior of as-built and re-polished SLM samples

The aim of this part is to better understand the fatigue mechanisms for parts manufactured by the SLM process. The originality of this study is to compare specimens with skin-process (as built SLM surface) and polished specimen (without skin-process), using a 220 SiC paper resulting in a smoothing effect (Sa ≈ 0.8 μm for polished samples compared with 7 μm for the as-built condition).

Selection of optimum SLM parameters

Following process investigations carried out in the previous two parts, SLM parameters were chosen for manufacturing fatigue samples, mostly considering the reduction of porosity rate and surface finish which are both expected to be major contributors to fatigue crack initiation. For this reason, SLM parameters, corresponding to high VED values (93 J/mm 3 for the contour and 87 J/mm 3 for the hatching), and a low building angle (10°) were used (Table 4).

Tensile tests

Prior to fatigue tests, two tensile tests were carried out under quasistatic load conditions (10 -3 s -1 ) directly on as-built SLM samples with the orientation shown in Fig. 2. Results presented here below, are at least as good than previous ALM works on Inconel 625 (Table 1), despite the presence of ≈0.3% porosity rate (Table 5).

Fatigue tests

Seven fatigue tests were conducted to determine the number of cycles to failure for as-built and polished specimens, at different stress amplitudes (with a 25 MPa increment). Even if the number of specimens is limited, the results and the failure surface obtained give important information and ways of improving SLM materials for fatigue applications.

All tests were carried out at a test frequency of 20 Hz with a load Rratio equal to -1 using the specimen geometry presented in Fig. 3. Specimens that survived 2 × 10 6 cycles have been retested at a higher stress level, increased by 25 MPa, to estimate the fatigue strength and to observe the fatigue damage mechanism on the failure surfaces. The results are summarized in Table 6 and plotted as a S-N Wöhler curve (Fig. 12).

The estimated fatigue limits are between 175 and 200 MPa (i.e. approximately 20% of the UTS). This is much lower than the fatigue strengths obtained with other manufacturing processes [START_REF] Koutiri | A probabilistic model for the high cycle fatigue behaviour of cast aluminium alloys subject to complex loads[END_REF]. Considering previous "SLM and fatigue" studies on other classes of materials, like an AlSI10Mg [START_REF] Brandl | Additive manufactured AlSi10Mg samples using selective laser melting (SLM): microstructure, high cycle fatigue and fracture behavior[END_REF] with a fatigue limit equal to approximately 30% of the UTS for the same building orientation, the fatigue strengths determined here can also be considered to be very low.

The fatigue results clearly indicate an increase in the number of cycles to failure (i.e. an average of approximately +2 × 10 5 cycles for stress amplitudes of 225, 250 and 275 MPa) for the polished specimens, even if a higher difference could have been expected between the two conditions, due to the difference in surface roughness. However, for the same number of cycles, the difference in stress amplitude, observed between both conditions remains low. Despite the high surface roughness of as-built specimen, the difference in the fatigue strength between the polished and the as-built specimens is not pronounced (about 10%). This low difference can be explained by the similar defect size observed at the initiation sites for the two studied batches. This point is developed in next section.

Analysis of failure mechanisms on SLM + re-polished specimens

The failure surface of polished specimens has been observed using a Scanning Electron Microscopy (SEM) in order to address the fatigue mechanisms.

Observations of specimens P 1 , P 6 and P 4 (Fig. 13) show that the initiation of the fatigue crack occurs on sub-surface without any defects like inclusions or porosity at the initiation site. For these specimens, initiation is localized in the material matrix and is controlled by local plasticity, according to the less favorably crystal orientation (Dang [START_REF] Van | Introduction to fatigue analysis in mechanical design by the multiscale[END_REF]. [START_REF] Nicoletto | Anisotropic high cycle fatigue behavior of Ti-6Al-4V obtained by powder bed laser fusion[END_REF] has observed the same mechanisms on smooth specimens in a Ti6Al4V obtained by powder bed laser fusion. Other fatigue initiation mechanisms are observed in samples P 2 , P 3 , P 5 and P 7 where cracks initiate from pores created during the SLM process. In Fig. 13 (sample P 3 ), a 60 μm diameter pore is located at 150 microns below the surface. It is believed that the associated porosity (which is not fully spherical) can be correlated to the sub-surface pores shown in Fig. 13, and attributed to a key-hole like welding mode or a lack of inter-layer dilution. [START_REF] Konecna | Microstructure and directional fatigue behavior of Inconel 718 produced by selective laser melting[END_REF] observed the same mechanism in an Inconel 718 on smooth specimens.

It is therefore concluded that fatigue damage in this material is controlled by two different fatigue mechanisms: one is associated with relatively large pores and the other is controlled by local plasticity within the material matrix.

Failure mechanisms of as-built SLM specimens

In order to characterize the fatigue mechanisms on specimens carried out with a skin ("contour") process, failure surfaces have also been observed on as-built SLM samples to explain the difference obtained with polished specimens.

The initiation site in specimen B 2 is at a large pore in which some un-melted particle powder can be observed (see Fig. 14). This initiation mechanism is the same as that observed for the polished specimens P 2 , P 3 , P 5 and P 7 . For the specimen B 5 which has a particular high fatigue strength there are no pores observed at the initiation site, as per the polished specimens P 1 , P 4 and P 6 .

The observations show that for specimens B 1 , B 3 , and B 4 a similar crack initiation mode occurs, controlled by the presence of large embedded particles, and localized on the surface of the specimen (Fig. 15 for specimen B1 and Fig. 16 for specimen B3). These particles have a diameter of about 150 μm which seems large compare to the grain size distribution (Fig. 1). Spatter ejection, generated during the powder bed laser beam melting process, discussed in recent papers [START_REF] Wang | Mechanisms and characteristics of spatter generation in SLM processing and its effect on the properties[END_REF][START_REF] Gunenthiram | Experimental analysis of spatter generation and melt-pool behavior during the powder bed laser beam melting process[END_REF] could explain the origin of this kind of defect.

Specimen B 1 shows an unusual failure surface, in which a competition between two different heterogeneities can be observed (Fig. 15B 1b) and B 1 -c)). This promotes the formation of two fatigue crack initiation sites (one at a sub-surface pore and the other at a 150 μm embedded particle).

The specimen B 3 was analyzed in greater detail prior to and after fatigue testing. Prior to fatigue tests, a 3D surface scan was carried out to detect and visualize the position and the size of embedded particles. It can be seen in Fig. 15B 3 -C that approximately twenty particles are present on the specimen surface. The critical defect can be identified on the surface scan using the fracture surface analysis. Fig. 15 shows the surface scan and a close-up of the critical defect of the as-built specimen B 3 .

On the failure surface, it can be observed that the interface between the particle and the material matrix is not completely smooth. A partially embedded particle seems to create an initial defect which could be very critical in fatigue. It can be observed in Fig. 16 that the 3D profilometer can be used to detect this kind of heterogeneity, as shown in Fig. 14B 3 -c/d.

Finally, three different fatigue damage mechanisms are observed. The most damaging is probably due to particle ejected during the SLM process and embedded on the surface. The second one is controlled by the degree of porosity, created during the process close to the surface. The third and least damaging fatigue crack initiation mechanism observed is crack initiation due to local plasticity which appears in the material matrix.

These observations tend to show that the usual roughness parameters like Sa or St, which gives a more global description of the surface are probably not suitable to characterize the surface effect on fatigue. It seems to be more efficient to quantify and optimize the particle population in term of density, size or geometry in order to better assess the fatigue strength of a SLM material. New surface parameters must also be developed to propose a link between the surface state and the fatigue strength of this type of material.

Conclusion and future work

The SLM manufacturing of Inconel 625 samples has been carried out with a double objective: (1) determining optimum process parameters leading to reduced surface finish and porosity rate and understanding underlying physics, and (2) investigating the fatigue behavior and failure mechanisms on as-built and (as-built + re-polished) samples.

The first part was focused on parameters to deliver the best surface finish. It appears that the correlation between surface roughness, which depends mostly on surface temperatures achieved on SLM part, and process parameters, is not a simple task, and that the use of VED as a surfaces showing the crack initiation location for the as-Built specimen B1. representative factor is not fully satisfactory, especially when considering a large range of laser spot diameters. Globally, Sa was shown to increase with lower scan speeds, to decrease with higher powers, and to increase severely on down-skin sides for large building angles.

In a second part, hatching parameters have been optimized to reach an optimum densification state. For a constant building angle (10°), the porosity rate globally decreases with VED values (J/mm 3 ) above a VED threshold. Logically, high enough VED promote large melt-pools susceptible to better re-melt inter-layer porosities. The inclination of parts (from 10 to 35°) also has a clear effect on densification rates because higher inclinations increase average temperatures due to heat concentration in the lower sides of SLM walls. Contrary to the surface roughness, the possibility to use VED as a rather reliable process parameter combining P, V and D, for a constant building angle θ and layer height Δh is confirmed.

Finally, this process optimization was found to be a compromise between an optimum densification state and a minimum number of contaminating spatters which degrade the surface roughness.

An experimental fatigue campaign has been conducted on as-built specimens and polished specimens. Despite the poor surface roughness of as-built specimen, the polished and the as-built fatigue batches exhibit similar fatigue strength. This can be explained by the defects associated with crack initiation sites which are similar in size for the studied specimens.

Three different fatigue damage mechanisms have been observed and described:

-Particles probably ejected during the SLM process, which stick to and become embedded on the surface. -Porosity created during the SLM process, located on the surface or sub-surface of the material. -Local plasticity which appears in the matrix material, depending on material microstructure.

This study provides an initial step for future modeling work to describe the fatigue behavior of as-built SLM specimens. For instance, in the future, a probabilistic approach using the weakest link concept [START_REF] Pessard | Modelling the role of non metallic inclusions on the anisotropic fatigue behaviour of forged steel[END_REF] could be introduced to model the competition between the three damage mode that have been observed.

Fig. 1 .

 1 Fig. 1. (a) SEM analysis of powder particles, (b) Inconel 625 Grain distribution.

Fig. 2 .

 2 Fig. 2. Workflow of the SLM pre-processing.

Fig. 3 .

 3 Fig. 3. Geometry of the specimens tested in fatigue.

Fig. 4 .

 4 Fig. 4. Sa Roughness vs Laser Power for a tilt angle equal to 10°and 35°(upper surface) and a scan speed of 0.5 m/s.

Fig. 5 .

 5 Fig. 5. (a) Average roughness Sa measured on bottom and top surfaces for aθ = 10°building angle, (b) SEM image for P = 140W, V = 0.5 m/s (sample 2, upper side), (c) SEM image for sample 2 (lower side).

Fig. 6 .

 6 Fig. 6. 3D maps (2.5 mm × 3 mm) of the SLM surfaces for P = 100, 140, 180 W and θ = 10, 35°(top surfaces) (V = 0.6 m/s).

Fig. 7 .

 7 Fig. 7. Non-melted particle (D ≈ 100 μm) on a SLM surface (a) 3D scan and (b) SEM image.

Fig. 8 .

 8 Fig.8. Effect of beam diameter on the Sa parameter (P = 400 W, contour strategy, 0°a ngle).

Fig. 9 .

 9 Fig. 9. Dependence of porosity content versus process parameters in In625 for P < 200 W and P = 400 W (10°angle).

Fig. 10 .

 10 Fig. 10. Dependence of porosity content versus process parameters in In625 for two distinct building angles (10°and 35°).

Fig. 11 .

 11 Fig. 11. 2D quantification of the voids fraction a) in the cross-section (0.3% of spherical porosities) b) close to the surface (2.2% of porosity) (sample 24: 150 W-0.6 m/s-10°).

Fig. 12 .

 12 Fig. 12. S-N curves for polished and as-built batches.

Fig. 13 .

 13 Fig. 13. SEM images of fracture surfaces showing the crack initiation location for the polished samples.

Fig. 14 .

 14 Fig. 14. SEM images of fracture surfaces showing the crack initiation location for the as-Built specimen B1.

Fig. 16 .

 16 Fig. 16. B3-a and B3-b: SEM images of fracture surfaces showing the crack initiation located at a surface particle; B3-c and B3-d: Surface scan at different scale of the as-built specimen showing the critical particle at the origin of the failure initiation.

Table 1

 1 Mechanical properties of Inconel 625 (Literature data).

	Ref.	Process	UTS (MPa)	Yield strength (MPa)	Elongation (%)	P (W)	V (m/s)	Powder thickness (μm)
	(Yadroitsev et al., 2007)	Wrought Inconel 625	930	517	40
	(Yadroitsev et al., 2007)	SLM Horizontal samples	1030	800	8-10	50	0.13	50
	(Yadroitsev et al., 2007)	SLM Vertical samples	1070	720	8-10	50	0.13	50
	(Wang et al., 2016)	SLM Horizontal samples	878	641	30	360	0.4	50

Table 2

 2 Process parameters (laser power, scan speed, spot diameter and tilt angle) and resulting surface finish.

	Number of specimen		θ (°)	Contour parameters		VED (contour))	Hatching parameters		Roughness*
				P (W)	V (m/s)	D (μm)	J/mm 3	P (W)	V (m/s)	D (μm)	Sa (μm)
	Tests for Surface roughness	1	10	100	0.5	70	52	175	0.6		11/14
		2	10	140	0.5	70	72.8	175	0.6		7.5/11.8
		3	10	180	0.5	70	93.6	175	0.6		5.9/6.2
		4	35	100	0.5	70	52	175	0.6		11/16.5
		5	35	140	0.5	70	72.8	175	0.6		7.5/20.5
		6	35	180	0.5	70	93.6	175	0.6		10.6/23
		7	0	400	1	70	104	400	1		10.5
		8	0	400	1.4	70	74.3	400	1.4		12.3
		9	0	400	1.8	70	57.8	400	1.8		12.4
		10	0	400	2	70	52	400	2		9.5
		11	0	400	1	160	19.9	400	1		21
		12	0	400	1.4	160	14.2	400	1.4		18.9
		13	0	400	1.8	160	11.1	400	1.8		17.1
		14	0	400	2	160	10	400	2		15.2
		15	0	400	1	240	8.8	400	1		21
		16	0	400	1.4	240	6.3	400	1.4		20.3

(* bold = upper surface, normal = lower surface).

Table 3

 3 Set of SLM Parameters and porosity ratios.

		Test	θ (°)	Contour parameters		Hatching parameters		VED (hatching)	Porosity
				P (W)	V (m/s)	D (μm)	P (W)	V (m/s)	D (μm)	J/mm 3	(%)
	Tests for porosity	7	0	400	1	70	400	1	70	104	0,51 ± 0.05
		8	0	400	1.4	70	400	1.4	70	74.3	1,17 ± 0.4
		9	0	400	1.8	70	400	1.8	70	57.8	1,85 ± 0.34
		10	0	400	2	70	400	2	70	52	2,43 ± 0.9
		11	0	400	1	160	400	1	160	19.9	0,36 ± 0.09
		12	0	400	1.4	160	400	1.4	160	14.2	0,91 ± 0.12
		13	0	400	1.8	160	400	1.8	160	11.1	2,96 ± 0.8
		14	0	400	2	160	400	2	160	10	2,7 ± 0.4
		15	0	400	1	240	400	1	240	8.8	1,16 ± 0.09
		16	0	400	1.4	240	400	1.4	240	6.3	1,30 ± 0.6
		17	0	400	1.8	240	400	1.8	240	5,1	5,06 ± 0.2
		18	0	400	2	240	400	2	240	4,6	7,20 ± 0.3
		19	0	180	0,5	70	175	0,6	70	79,3	1,44 ± 0,4
		20	0	140	0,5	70	175	0,6	70	79,3	0,7 ± 0,16
		21	0	140	0,5	70	175	0,7	70	68	1,29 ± 0,22
		22	10	140	0.5	70	150	0.6	70	68	0.2 ± 0.08
		23	10	140	0.5	70	200	0.6	70	90,7	0.37 ± 0.18
		24	10	140	0.5	70	175	0.7	70	68	2.1 ± 0.51
		25	10	140	0.5	70	175	0.8	70	59,5	6.32 ± 1.26
		26	35	140	0.5	70	150	0.6	140	17	0.36 ± 0.04
		27	35	140	0.5	70	200	0.6	140	22,6	0.12 ± 0.02
		28	35	140	0.5	70	175	0.7	140	17	0.66 ± 0.06
		29	35	140	0.5	70	175	0.8	140	14,8	1.61 ± 0.30

Table 4

 4 SLM Process parameters carried out on static tensile and fatigue samples.

	θ (°)	Contour parameters		Hatching parameters
		P (W)	V (m/s)	D (μm)	P (W)	V (m/s)	D (μm)
	10	180	0.5	70	175	0.6	70

Table 5

 5 Results of tensile tests on as-built SLM samples.

		E (GPa)	σ Y (MPa)	A (%)	UTS (MPa)
	Sample 1	177	711	35	976
	Sample 2	179	727	35	971

Table 6

 6 Fatigue data, fully reversed tensile loads (R = -1).

	Polished specimens			As-built specimens
	Specimen Name	σ a (MPa)	N (cycles)	Specimen Name	σ a (MPa)	N (cycles)	Sa (μm)
	P 1	300	678 663	B 1	275	403 887 e	8.56
	P 2	275	662 916 p	B 2	250	573 095 p	9.85
	P 3	250	779 395 p	B 3	225	901 787 e	10.69
	P 4	225	1 641 598	B 4	200	2 000 000	5.46
	P 5	200	1 373 450 p	B 4	225	1 069 045 e
	P 6	175	2 000 000	B 5	200	2 000 000	4.82
	P 6	200	1 190 074	B 5	225	2 000 000
	P 7	175	2 000 000	B 5	250	35 5349
	P 7	200	261 270 p
	p initiation from a porosity.

e initiation from an embedded particle.

Acknowledgement

This work was carried out under the framework of the PSPC project FAIR. Authors thank BPI for financial support.


Adrien Scheuer
email: adrien.scheuer@ec-nantes.fr

Emmanuelle Abisset-Chavanne
email: emmanuelle.abisset-chavanne@ec-nantes.fr

Francisco Chinesta
email: francisco.chinesta@ensam.edu

Roland Keunings
email: roland.keunings@uclouvain.be

Microscopic modelling of orientation kinematics of non-spherical particles suspended in confined flows using unilateral mechanics

Keywords:

The properties of reinforced polymers strongly depend on the microstructural state, that is, the orientation state of the fibres suspended in the polymeric matrix, induced by the forming process. Understanding flow-induced anisotropy is thus a key element to optimize both materials and process. Despite the important progresses accomplished in the modelling and simulation of suspensions, few works addressed the fact that usual processing flows evolve in confined configurations, where particles characteristic lengths may be greater than the thickness of the narrow gaps in which the flow takes place. In those circumstances, orientation kinematics models proposed for unconfined flows must be extended to the confined case. In this short communication, we propose an alternative modelling framework based on the use of unilateral mechanics, consequently exhibiting a clear analogy with plasticity and contact mechanics. This framework allows us to revisit the motion of confined particles in Newtonian and non-Newtonian matrices. We also prove that the confined kinematics provided by this model are identical to those derived from microstructural approaches (Perez et al. (2016) [1]).

Introduction

Over the last decades, composite materials, made of a suspending matrix and a reinforcement composed of fibres used to fortify the matrix in terms of strength and stiffness, were successfully introduced in the aerospace and automotive industries and proved to be a lightweight alternative to produce structural and functional parts. The mechanical properties of such reinforced polymers, however, strongly depend on the orientation state of their microstructure, which is established during the forming process [START_REF] Advani | Flow and Rheology in Polymer Composites Manufacturing[END_REF]. Predicting the evolution of this orientation state is thus a key yet complex task, since the motion of the reinforcing fibres is impacted by the flowing matrix and the interactions with the neighbouring fibres and cavity walls. Thus, flow-induced anisotropy must be understood and modelled in order to optimize both materials and processes.

There is a long history of studies of the motion of slender bodies suspended in a viscous fluid, starting from the seminal work of Jeffery back in 1922 [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF]. A vast literature dedicated to fibre and non-spherical particle suspensions is available, studying extensively different modelling scales and exploring the impact of the concentration regime and the nature of the suspending matrix. Schematically, the three main modelling scales involved when addressing the orientation kinematics of suspended particles can be summarized as follows: (i) the microscopic scale, the scale of a single particle; (ii) the mesoscopic scale, the scale of a population of particle, whose conformation is usually represented by a probability density function (pdf) of orientation, giving an unambiguous and complete description of the orientation state; and (iii) the macroscopic scale, the scale of the part, whose conformation is often given by the first moments of the aforementioned pdf, providing a coarse yet concise description of the orientation state in the part. Depending on the level of detail and accuracy required for a given application, a specific scale, or a combination of them might be chosen. For further detail on that subject, including the so-called multiscale approach, we refer the reader to the review [START_REF] Petrie | The rheology of fibre suspensions[END_REF] and the monograph [START_REF] Binetruy | Flows in Polymers, Reinforced Polymers and Composites. A Multiscale Approach[END_REF], as well as to the references therein. Only a succinct overview of the microscopic modelling is proposed thereafter.

In [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF], Jeffery derived the expression of the hydrodynamic torque exerted on an ellipsoidal particle immersed in an unbounded creeping flow of a Newtonian fluid. He then obtained an equation of motion by assuming that the particle rotates so as to achieve instantaneous zero torque, resulting in the so-called Jeffery equation. Considering a spheroid (axisymetric ellipsoid) and defining the orientation of a particle by the unit vector p along its principal axis, Jeffery's equation reads

ṗ J = • p + κ(D • p -(D : (p ⊗ p)p)) (1)
where and D are, respectively, the skew-symmetric and symmetric parts of the unperturbed velocity gradient ∇v of the flow, and κ is the shape factor of the spheroid, given by κ = r 2 -1 r 2 +1 with r the aspect ratio of the particle. Slender bodies like fibres and rods can be assimilated as infinite aspect ratio ellipsoids (κ ≈ 1).

Jeffery's equation was experimentally verified by Taylor [START_REF] Taylor | The motion of ellipsoidal particles in a viscous fluid[END_REF] and Mason [START_REF] Trevelyan | Particle motions in sheared suspensions I. Rotations[END_REF]. Bretherton [START_REF] Bretherton | The motion of rigid particles in a shear flow at low Reynolds number[END_REF] showed that the equation is also valid for any axisymmetric particle provided that an effective aspect ratio is determined. Hinch and Leal [START_REF] Hinch | The effect of Brownian motion on the rheological properties of a suspension of non-spherical particles[END_REF][START_REF] Hinch | Constitutive equations in suspension mechanics. Part I[END_REF][START_REF] Hinch | Constitutive equations in suspension mechanics. Part II[END_REF][START_REF] Hinch | Rotation of small non-axisymmetric particles in a simple shear flow[END_REF] also studied Jeffery's model, addressing the impact of Brownian motion and proposing constitutive equations for the behaviour of suspensions. However, only a few works, either experimental [START_REF] Moses | Investigation of the fiber motion of fiber motion near solid boundaries in simple shear flow[END_REF][START_REF] Stover | The motion of rodlike particles in the pressure-driven flow between two flat plates[END_REF] or numerical and theoretical [START_REF] Jayageeth | Dynamics of short fiber suspensions in bounded shear flow[END_REF][START_REF] Ozolins | Simple models for wall effect in fiber suspension flows[END_REF], address the fact that flows of industrial interest take place in narrow gaps, whose thickness is of the same order of magnitude or smaller than the length of the reinforcing fibres, thus constraining the space of possible orientation, and as a consequence the kinematics.

In [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF], we proposed a multiscale model to describe the orientation development of a dilute suspension of fibres confined in a narrow gap. The microscopic model was based on a dumbbell representation [START_REF] Bird | Dynamic of Polymeric Liquid[END_REF] of the confined rod, with hydrodynamic and contact forces (normal to the gap wall) acting on it. This confinement model was later extended in [START_REF] Scheuer | Second-gradient modelling of orientation development and rheology of dilute confined suspensions[END_REF] to include unilateral contacts and non-uniform strain rates at the scale of the rod. In any case, the resulting kinematics are a combination of the unconfined Jeffery kinematics and a correction term that prevents the fibre from leaving the flow domain, that is,

ṗ = ṗ J + ṗC (2)
The equation of motion of such a confined rod, derived in [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF] and summarized in Eq. ( 2), presents thus significant similarities with equations of elastoplaticity. Indeed, we could draw a parallel between, on the one hand, the classical unconfined kinematics and the elastic deformation, and, on the other hand, the confined motion of the particle and elastoplatic deformation.

Hence, the purpose of this short communication is to explore the alternative modelling framework based on unilateral mechanics to revisit the motion of confined particles in Newtonian and non-Newtonian matrices.

The paper is organised as follows. In Section 2, we derive the model for the confined kinematics of suspended particles using unilateral mechanics. Then, in Section 3, we discuss how this general framework allows us to build the confined kinematics of fibres and spheroids immersed in a Newtonian (based on Jeffery's model [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF]) or second-order (based on Brunn's model [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF]) viscoelastic fluid. Finally, in Section 4, we draw the main conclusions and present some perspectives of this approach.

Remark. In this paper, we consider the following tensor products, assuming Einstein's summation convention:

-if a and b are first-order tensors, then the single contraction "•" reads (a • b) = a j b j ; -if a and b are first-order tensors, then the dyadic product "⊗" reads (a ⊗ b) jk = a j b k ; -if a and b are respectively second and first-order tensors, then the single contraction "•" reads (a • b) j = a jk b k ; -if a and b are second-order tensors, then the double contraction ":" reads (a : b) = a jk b kj .

Confined orientation kinematics using unilateral mechanics

In this section, we derive step by step the kinematics of a confined suspended particle using the framework of unilateral mechanics.

(i) Additive decomposition. We assume that the particle kinematics (particle rotary velocity) can be decomposed into an unconfined (U) and confined (C) contribution, according to

ṗ = ṗU + ṗC (3)
By definition, the orientation vector is subject to the normalization condition

p • p = 1 (4)
(ii) Unconfined contribution. The unconfined equation of motion ṗU is either given by a model from the literature, e.g., Jeffery's equation [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF] for ellipsoidal particles immersed in a Newtonian fluid or Brunn's model [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF] in the case of a second-order (non-Newtonian) matrix, or could be estimated from experimental observations. In order to satisfy the normalization condition Eq. ( 4), the unconfined kinematics should verify ṗU • p = 0. (iii) Allowed domain and confinement condition. We define a function f : S → R called the confinement condition (the equivalent of the yield condition in plasticity mechanics) and constrain the admissible orientation states p ∈ S (with S the unit sphere) to lie in the flow domain D, ensuring the gap walls' impenetrability, defined as

D = p ∈ S | f (p) = p • n - H L ≤ 0 ( 5
)
where L is the semi-length of the particle, H the gap semi-width, and n denotes a unit vector normal to the gap wall. We assume, without loss of generality, that p points towards the upper hemisphere.

We refer to the interior of D, denoted by int(D) and defined as intD) = {p ∈ S | f (p) < 0}, as the unconfined domain, whereas the boundary ∂D, given by ∂D = {p ∈ S | f (p) = 0}, is called the confinement surface. (iv) Confined contribution and consistency requirement. The confined kinematics ṗC is obtained from the gradient of the con- finement condition introduced above, i.e. ṗC = -γ

d f dp ( 6
)
where γ is the consistency parameter. The derivative of f (p) with respect to p, enforcing the normalization condition,

reads d f dp = n -(p • n)p (7)
which is obtained by subtracting from the derivative of f (p) with respect to p its projection onto the p direction in order to ensure that ṗC • p = 0.

The consistency parameter γ is assumed, on the one hand, to obey the Kuhn-Tucker complementary conditions γ ≥ 0,

f (p) ≤ 0,
and

γ f (p) = 0 ( 8 )
and, on the other hand, to satisfy the consistency requirement

γ ḟ (p) = 0 (9)
To obtain the derivative of f with respect to time, we proceed as follows:

ḟ = d f dp • ṗ (10) = (n -(p • n)p) • ( ṗU -γ (n -(p • n)p)) (11)
= ṗU • n -γ (1 -(p • n) 2 ) ( 12
)
Thus the value of γ is given by

γ = ( ṗU • n) (1 -(p • n) 2 ) ( 13
)
(v) Summary. The resulting kinematics can be summarized as follows:

f < 0 ⇐⇒ p ∈ int(D), γ = 0 =⇒ ṗ = ṗU unconfined f = 0 ⇐⇒ p ∈ ∂D, ⎧ ⎪ ⎨ ⎪ ⎩ ḟ = 0 and γ > 0 =⇒ ṗ = ṗU + ṗC ḟ = 0 and γ = 0 =⇒ ṗ = ṗU ḟ < 0 =⇒ γ = 0 =⇒ ṗ = ṗU confinement force-free contact unconfined detachment (14)
The force-free contact actually corresponds to the case where the confined particle is touching the gap wall, but is not trying to leave the flow domain (to "push" on the wall). Remark 2.1. For the sake of clarity, we assumed in the present section and throughout the remainder of this article that under confinement, both extremities of the suspended particle are in contact with the gap walls and the centre of gravity of the particle is fixed in the mid-plane of the flow channel. These assumptions are, however, relaxed in Appendix A, where we address, using the same framework just introduced, the general case by introducing the position of the particle centre of gravity in the confinement condition, allowing in the meantime contact with only one gap wall.

Discussion

The confined kinematics obtained in the previous section read ṗ = ṗU + ṗC = ṗU -

( ṗU • n) (1 -(p • n) 2 ) (n -(p • n)p) (15)
which coincides exactly with the expression obtained in [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF] following a microstructural approach for a confined rod immersed in a Newtonian fluid ( ṗU was thus given by Jeffery's kinematics).

However, the unilateral mechanics approach developed in the previous section does not assume anything on the shape of the suspended particles or the nature of the matrix fluid. Only an expression of the unconfined kinematics is actually necessary, which allows us to extend straightforwardly our model describing the motion of a confined particle to situations where the microstructural approach from [1] might be tedious.

In this section, we thus discuss three scenarios: (i) a rod immersed in a Newtonian fluid (developed in [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF]), (ii) a spheroid in a Newtonian fluid, and (iii) a spheroid in a second-order fluid. We also provide some numerical illustrations in the case 15), we recover the kinematics derived in [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF] using a microstructural approach (hydrodynamic and contact forces acting on the dumbbell representation of a rod). Fig. 1 depicts the evolution of the orientation of a rigid fibre immersed in a shear flow of a Newtonian fluid. The solid line shows the trajectory of the confined rod and the dotted line the trajectory of an hypothetical unconfined particle starting from the same initial orientation. Fig. 1a shows the evolution of the rod orientation as a trajectory on the unit sphere. Starting from the same initial condition, the confined (solid red) and unconfined (dotted blue) rods both follow the same Jeffery orbit. When it touches the wall, the confined rod is constrained to slide along the gap wall, and finally catches an unconfined Jeffery orbit tangent to the wall and aligns in the flow. This abrupt change in the trajectory as it touches the gap wall can also be observed in Fig. 1b, where the components of the unit vector of orientation p are represented.

Confined ellipsoid suspended in a Newtonian fluid

Considering spheroidal particles (axisymmetric ellipsoids), we insert now the general Jeffery equation [Eq. ( 1)] in Eq. ( 15). The equivalence with the kinematics obtained from a microstructural approach on a tri-dumbbell is detailed in Appendix B. Newtonian fluid. The unconfined particle (dotted blue) undergoes its classical kayaking motion, whereas the confined one (solid red) first slides along the wall, and is then constrained on the largest kayaking orbit possible in the narrow gap, that is, the orbit tangent to both gap walls. Note that the orbit points tangent to the channel walls thus correspond to force-free contacts.

Confined ellipsoid suspended in a second-order fluid

Leal [START_REF] Leal | The slow motion of slender rod-like particles in a second-order fluid[END_REF] and Brunn [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF] published some important theoretical works to describe, respectively, the motion of a rod and an ellipsoid in a second-order fluid, in the limit of low Weissenberg numbers, which constitutes the counterpart of Jeffery's equation in the case of a viscoelastic suspending matrix. The constitutive equation for the second-order fluid is given by Giesekus [START_REF] Giesekus | Die simultane translations-und rotations bewegung einer kugel in einer elastovisken flussigkeit[END_REF] σ = -pI + 2ηD + 2η k (11)

0 D • D + k (2) 0 ∂D ∂t + v • ∇D + • D -D • ( 16
)
where k [START_REF] Hinch | Constitutive equations in suspension mechanics. Part II[END_REF] 0 and k

(2) 0 are material parameters. Brunn [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF] derived the equation of evolution of the orientation of a particle, and the resulting kinematics reads

ṗB = • p + κ(D • p -(D : (p ⊗ p)p)) -(I -p ⊗ p) • D • (H 2 D • p + H 1 (D : (p ⊗ p))p) (17)
where H 1 and H 2 depend on the material parameters and the ellipsoid aspect ratio as given in Brunn [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF].

In Appendix C, we briefly show how to write Brunn's kinematics [Eq. [START_REF] Bird | Dynamic of Polymeric Liquid[END_REF]] in the form of Jeffery's kinematics Eq. ( 1)

with an effective velocity gradient ∇v. Consequently, the microstructural validation developed in Appendix B can also be used in this section.

In the case of a second-order fluid, it is well known that the kayaking motion of the particle drifts towards the shear plane for oblate spheroids (r < 1) or towards the vorticity axis for prolate spheroids (r > 1) [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF][START_REF] Borzacchiello | Orientation kinematics of short fibres in a second-order viscoelastic fluid[END_REF][START_REF] D'avino | Particle dynamics in viscoelastic liquids[END_REF]. Fig. 3 depicts the evolution of the orientation of a spheroid of aspect ratio r = 4 immersed in the same shear flow, but now the suspending matrix is a second-order fluid (in this illustration, we have k = -0.09, and thus, H 1 = 0.084 and H 2 = 0.0056). Again, the confined spheroid (solid red) is constrained to exert its kayaking and drifting motion towards the vorticity axis of the flow in the narrow flow domain.

Conclusion and perspectives

This work proposes an alternative route for deriving the confined orientation kinematics of dilute suspensions of fibres and ellipsoidal particles. This approach, based on unilateral mechanics, allows us to extend directly our previous microstructural model developed for confined rods in a Newtonian fluid [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF] to ellipsoidal particles and viscoelastic matrices.

The same strategy might then be applied at the macroscopic scale, working on the so-called second-order orientation tensor a [START_REF] Advani | The use of tensors to describe and predict fibre orientation in short fibre composites[END_REF], a = (p ⊗p) ψ(p) dp, to derive a macroscopic model for confined suspensions. However, defining the adequate confinement condition f (a) ≤ 0 in that case is a delicate task that will be addressed in a future work.

Appendix B. Confined kinematics of an ellipsoid using the dumbbell approach

In order to address the confined kinematics of an ellipsoid immersed in a Newtonian fluid flow, we consider its corresponding tri-dumbbell representation. For the sake of simplicity, we consider in this appendix the 2D case, that is the bi-dumbbell represented in Fig. 4, which will provide the orientation kinematics of a spheroid (axisymmetric ellipsoid). The orientation of the two principle axes are given by the unit vectors p 1 and p 2 .

On each bead acts a hydrodynamic force (Stokes drag), scaling with the difference of velocity between the fluid at the bead position and the bead itself. For the bead located at p L 1 , that force reads

F H (p 1 L 1 ) = ξ(v 0 + ∇v • p 1 L 1 -v G -ṗ1 L 1 ) (26)
where ξ is a friction coefficient and v G and v 0 denote respectively the velocity of the centre of gravity G and the velocity of the fluid at G. Without any loss of generality, we assume that L 1 is the spheroid longest axis, on which contact with the gap wall will occur. The contact forces on the upper and lower beads read

F C (p 1 L 1 ) = μn and F C (-p 1 L 1 ) = -μn (27)
with n = 0 0 1 T (since the contact force is orthogonal to the wall as soon as friction is neglected -it could be easily introduced), and its intensity μ is a priori unknown.

When both beads of the dumbbell aligned along direction p 1 are in contact with the gap walls (more general situations were analysed in [START_REF] Scheuer | Second-gradient modelling of orientation development and rheology of dilute confined suspensions[END_REF]), balance of forces yields 2ξ(v 0v G ) = 0 (28) that is, v 0 = v G , the particle's centre of gravity is moving with the fluid velocity at that position.

Balance of torques thus

yields 2p 1 L 1 × ∇v • p 1 L 1 -ṗ1 L 1 + μ ξ n + 2p 2 L 2 × (∇v • p 2 L 2 -ṗ2 L 2 ) = 0 (29)
Introducing the fact that ṗi = ω × p i , i = 1, 2 (with ω the angular velocity vector), and taking into account that

p i × ω × p i = ω, i = 1, 2, we have ω = L 1 L 2 1 + L 2 2 (p 1 × (∇v • p 1 L 1 + μ ξ n)) + L 2 L 2 1 + L 2 2 (p 2 × (∇v • p 2 L 2 )) (30)
Thus, the spheroid rotary velocity is given by

ṗ1 = ω × p 1 = L 1 L 2 1 + L 2 2 ((p 1 × (∇v • p 1 L 1 + μ ξ n)) × p 1 ) + L 2 L 2 1 + L 2 2 ((p 2 × (∇v • p 2 L 2 )) × p 1 ) (31)
Applying the triple vector product formula (a

× b) × c = (a • c)b -(b • c)a, the previous equation reads: ṗ1 = L 1 L 2 1 + L 2 2 ((∇v • p 1 L 1 + μ ξ n) -((∇v • p 1 L 1 ) • p 1 + μ ξ (n • p 1 ))p 1 ) - L 2 L 2 1 + L 2 2 (((∇v • p 2 L 2 ) • p 1 )p 2 ) ( 32
)
We now develop the last term of this equation in order to obtain an expression that only depends on p 1 . First, we decompose the velocity gradient according to ∇v = D + , ((∇v

• p 2 ) • p 1 )p 2 = ((D • p 2 ) • p 1 )p 2 + (( • p 2 ) • p 1 )p 2 (33)
Then we develop each of these terms. Using the fact that D is symmetric, we have ((∇v Imposing the impenetrability condition [START_REF] Perez | A multi-scale description of orientation in simple shear flows of confined rod suspensions[END_REF] (v G + ṗ1 L 1 )

((D • p 2 ) • p 1 )p 2 = p 2 (p T 1 • D • p 2 ) = p 2 (p T 2 • D • p 1 ) = (p 2 ⊗ p 2 ) • D • p 1 (34)
• p 1 L 1 + μ ξ n) -((∇v • p 1 L 1 ) • p 1 + μ ξ (n • p 1 ))p 1 ) - L 2 2 L 2 1 + L 2 2 (D • p 1 -(p 1 ⊗ p 1 ) • D • p 1 -• p 1 + (p 1 ⊗ p 1 ) • • p 1 ) (39)
• n = 0 ( 42
)
allows us to obtain the intensity μ of the contact force. Since the particle centre of gravity is in the mid-plane of the shear flow, v G = 0 and thus

ṗ1 L 1 • n = ṗ1 J L 1 • n + L 2 1 L 2 1 + L 2 2 ( μ ξ (1 -(n • p 1 ) 2 )) = 0 (43)
The value of μ is then given by μ = -

L 2 1 + L 2 2 L 2 1 ξ L 1 (1 -(n • p 1 ) 2 ) ( ṗ1 J • n) (44)
Eventually, the confined kinematics of a rigid spheroid are given, as expected, by ṗ1 = ṗJ 1 -

( ṗ1 J • n) (1 -(n • p 1 ) 2 ) (n -(n • p 1 )p 1 ) (45)

Fig. 1 .

 1 Fig. 1. Confined rod suspended in a Newtonian fluid (solid line: confined particle; dotted line: unconfined particle).

of a linear shear flow v = γ z 0 0 T 3 . 1 .

 031 , with γ = 1 in a channel of height H = 0.25 L Confined rod suspended in a Newtonian fluid Inserting Jeffery's equation for rods ṗJ rod = ∇v • p -(∇v : (p ⊗ p)p) as unconfined kinematics in Eq. (

Fig. 2 .

 2 Fig. 2. Confined spheroid (aspect ratio r = 4) suspended in a Newtonian fluid (solid line: confined particle; dotted line: unconfined particle).

Fig. 2

 2 Fig. 2 depicts the evolution of the orientation of a spheroid of aspect ratio r = 4 immersed in the same shear flow of a

Fig. 3 .

 3 Fig. 3. Confined spheroid (aspect ratio r = 4) suspended in a second-order fluid (solid line: confined particle -dotted line: unconfined particle).

Fig. 4 .

 4 Fig. 4. Hydrodynamic and contact forces acting on a confined suspended spheroid.

or ṗ1 = • p 1 + L 2 1 -L 2 2 L 2 1 + L 2 2 D • p 1 - L 2 1 -L 2 2 L 2 1 + L 2 2 (p 1 ⊗ p 1 ) • D • p 1 2 1 -L 2 2 L 2 1 +L 2 2 = r 2 - 1 r 2

 222222221112222212 (n • p 1 )p 1 )) (40)The first part of Eq. (40) is actually the classical Jeffery equation for a spheroid, since L +1 = κ is the spheroid shape factor and the second part is thus the confinement contribution, (n • p 1 )p 1 ))(41)

  Since p 1 and p 2 are mutually(p 1 ⊗ p 1 ) + (p 2 ⊗ p 2 ) = I • p 2 ) • p 1 )p 2 = (I -(p 1 ⊗ p 1 )) • D • p 1 (36)Similarly, using the fact that is skew-symmetric, we have(( • p 2 ) • p 1 )p 2 = p 2 (p T 1 • • p 2 ) = -p 2 (p T 2 • • p 1 ) = -(p 2 ⊗ p 2 ) • • p 1 • p 2 ) • p 1 )p 2 = -(I -(p 1 ⊗ p 1 )) • • p 1

		(35)
	and thus Eq. (34) now reads
	((D (37)
	and thus
	(( (38)
	Finally, coming back to Eq. (32), we have
	ṗ1 =	L 1 1 + L 2 L 2 2

Acknowledgements

A. Scheuer is a Research Fellow of the "Fonds de la recherche scientifique de Belgique" -F.R.S.-FNRS.

Appendix A. Confined orientation kinematics using unilateral mechanics -interaction with only one gap wall

In this appendix, we extend the framework introduced in Section 2 in situations where the suspended particle's centre of gravity does not necessarily lies in the channel mid-plane, allowing interaction with only one gap wall. In that case, the confinement condition now reads

with x G the position of the particle centre of gravity G. The allowed domain D is thus now given by

Note that in this case, we consider the case where G lies in the upper-half of the channel.

Using these new definitions, the derivation of the confined kinematics is obtained following the same rationale as that presented in the step (iv) in Section 2. The derivation of f with respect to p [Eq. [START_REF] Trevelyan | Particle motions in sheared suspensions I. Rotations[END_REF]] is left unchanged, whereas the derivation of f with respect to time [Eq. [START_REF] Hinch | Constitutive equations in suspension mechanics. Part I[END_REF]] now reads

with v G = ẋG the velocity of the particle's centre of gravity. Finally, the value of the consistency parameter γ is given by

leading to the following expression for the confined kinematics

This expression is exactly the one that we derived in [START_REF] Scheuer | Second-gradient modelling of orientation development and rheology of dilute confined suspensions[END_REF] when addressing the kinematics of a confined particle interacting with only one gap wall using a microstructural approach (dumbbell).

Appendix C. Rewriting Brunn's kinematics with an effective velocity gradient

Brunn's orientation kinematics for a spheroid particle immersed in a second-order fluid reads [START_REF] Brunn | The slow motion of a rigid particle in a second-order fluid[END_REF]:

where the spheroid shape factor is given by κ = r 2 -1 r 2 +1 , with r the aspect-ratio of the particle, and H 1 and H 2 are given by Brunn,

0

. Developing Eq. ( 46), we have

We now define an effective velocity gradient:

and the effective vorticity and strain rate tensors ˜ =


M Batteux

T Prosvirnova

A Rauzy Mtp

Enhancement of the AltaRica 3.0 stepwise simulator by introducing an abstract notion of time

AltaRica 3.0 is an event-based, object-oriented modeling language dedicated to (probabilistic) safety analyses of complex systems. It makes it possible to design models at higher level than done with formalisms traditionally used for safety analyses (fault trees, Markov Chains, stochastic Petri nets, etc.), without increasing the complexity of calculations of risk indicators. Several assessment tools have been developed for AltaRica 3.0, including a stepwise simulator. This tool is of a great help for the design and the validation of AltaRica 3.0 models. It is the analog for modeling of debuggers for programming. In this article, we show how the AltaRica 3.0 stepwise simulator has been greatly enhanced by the introduction of an abstract notion of time. The key mathematical property is that abstract and concrete simulation are bisimilar: any concrete (timed, stochastic) execution can be simulated by an abstract execution and reciprocally any abstract execution corresponds to at least one concrete execution. This important result paves the way to the design of efficient model-checking algorithms, e.g. generators of sequences of events leading to a failure state.

INTRODUCTION

AltaRica 3.0 is an event-based, object-oriented modeling language dedicated to (probabilistic) safety analyses of complex systems [START_REF] Prosvirnova | The altarica 3.0 project for model-based safety assessment[END_REF]. It makes it possible to design models at higher level than done with formalisms traditionally used for safety analyses (fault trees, Markov Chains, stochastic Petri nets, etc.), without increasing the complexity of calculations of risk indicators.

The semantics of AltaRica 3.0 is defined in terms of stochastic guarded transition systems [START_REF] Rauzy | Guarded transition systems: a new states/events formalism for reliability studies[END_REF]). AltaRica 3.0 executions are similar to those of other discrete event modeling formulisms: each time a transition gets fireable, it is scheduled and possibly fired after a certain real-valued delay, see e.g. [START_REF] Cassandras | Introduction to Discrete Event Systems[END_REF][START_REF] Zimmermann | Stochastic Discrete Event Systems[END_REF]) for introductions to (stochastic) discrete event systems. Events labeling transitions may be either deterministic or stochastic. In the later case, AltaRica 3.0 provides both built-in distributions (exponential, Weibull, . . . ) and empirical distributions.

Several assessment tools have been developed for AltaRica 3.0, see e.g. [START_REF] Prosvirnova | Automated generation of minimal cutsets from altarica 3.0 models[END_REF][START_REF] Brameret | Automated generation of partial markov chain from high level descriptions[END_REF][START_REF] Aupetit | Improving performance of the AltaRica 3.0 stochastic simulator[END_REF], including a stepwise simulator. This tool is of a great help for the design and the validation of AltaRica 3.0 models. It makes it possible to perform interactive step by step simulations, i.e. to go forth and back in sequences of events, enabling in this way to track modeling errors, unexpected behaviors and so on. With that respect, stepwise simulators play a similar role for discrete event modeling as debuggers like GDB or DDD for programming, see e.g. [START_REF] Matloff | The Art of Debugging with GDB, DDD, and Eclipse[END_REF] for an introduction to the latter.

In this article, we show how the AltaRica 3.0 stepwise simulator has been greatly enhanced by the introduction of an abstract notion of time. So far, it did not consider the time at all. The reason was that it would have been much too tedious for the analyst to enter by hand the delay associated with a stochastic transition each time this transition gets fireable. Moreover, infinitely many real-valued delays can be chosen, letting the analyst pondering which one is the most suitable for her or his purpose. However, ignoring delays had a major drawback: the stepwise simulator allowed the firing of sequences of events with no counterpart in stochastic simulation and more generally that did not obey the timed semantics of AltaRica 3.0.

The idea is therefore to abstract away the time in stepwise simulation: each transition is now associated with a time interval. Firing a transition may modify the time intervals associated with already scheduled transitions. This idea is by no means new: it enters into the general framework of Cousot's abstract interpretation [START_REF] Cousot | Abstract interpretation: a unified lattice model for static analysis of programs by construction or approximation of fixpoints[END_REF]. The problem at stake was to make it work for the particular case of stochastic discrete event simulations. The key mathematical property here is that abstract and concrete simulations are bisimilar, see e.g. [START_REF] Milner | Communication and Concurrency[END_REF] for an introduction to this important notion: any concrete (timed, stochastic) execution can be simulated by an abstract execution and reciprocally any abstract execution corresponds to at least one concrete execution. In a word, abstract executions are in agreement with AltaRica 3.0 semantics. This important result paves the way to the design of efficient model-checking algorithms, see e.g. [START_REF] Clarke | Model Checking[END_REF] for an introduction. In particular, it makes it possible the design of generators of sequences of events leading to a failure state.

The remainder of this article is organized as follows. Section 2 presents an illustrative example. Section 3 recalls fundamental notions about timed and stochastic guarded transition systems. Section 4 introduces the abstract notion of time. Section 5 concludes this article and gives some perspectives.

ILLUSTRATIVE EXAMPLE

As an illustrative example, we shall consider a system made of two identical, periodically tested, components that evolve independently one another.

The behavior of such a component is represented by the state diagram pictured Figure 1.

The component alternates operation and test phases. It is initially working and starts with a first operation phase that lasts a constant time θ. All subsequent operation phases last a constant time π.

The component may fail when in operation, with a failure rate λ.

If the component is working when it enters a maintenance phase, this maintenance phase lasts a constant time τ . If, on the contrary, it is failed when it enters a maintenance phase, the duration of its repair is uniformly distributed between two values µ and ν.

Finally, the component is as-good-as-new after a repair. In the state diagram of Figure 1, transitions failure and repair are thus stochastic while transitions startTest1, startTest and completeTest are deterministic.

Table 1 shows a possible evolution of such a component.

Transition

Firing date

startTest1 d 1 = θ completeTest d 2 = d 1 + τ failure d 3 = d 2 + δ 1 , 0 < δ 1 ≤ π startTest d 4 = d 2 + π repair d 5 = d 4 + δ 2 , µ ≤ δ 2 ≤ ν startTest d 6 = d 5 + π completeTest d 7 = d 6 + τ . . . . . .
We can assume that θ and π are relatively big compared to τ , µ and ν and that τ is smaller than µ (itself smaller than µ).

In order not to have the two components A and B out of service due to test or maintenance at the same time, it is reasonable to take different values of θ for A and for B (the values of the other parameters being identical for the two components). For instance, if an operation phase lasts normally six months (π = 4380h), the component A can be tested after three months (A.θ = 2190h) while the component B is tested after six months (B.θ = 4380h). It this way, tests/maintenances of A and B are shifted by three months which improves the overall availability of the system. Table 2 gives some typical values of the parameters for components A and B that we shall use throughout the article. With this values of parameters in mind, the reader sees immediately the problem of using a stepwise simulator that does not consider delays associated to transitions. Many executions that would be impossible with a timed semantics become possible with a non-timed semantics, e.g.

. B.startTest1 -------→ . . . . A.startTest1 -------→ . A.completeTest ---------→ . A.startTest -------→ . A.completeTest ---------→ . . .
On other hand, asking the analyst to introduce interactively delays of stochastic transitions is not a practical solution. Not only it would be tedious, but it would let the analyst facing the choice of suitable delays, which gets quickly puzzling. Hence the need of an abstract notion of time which makes it possible to take into account delays of transitions without asking the analyst to enter them interactively.

To fulfill this need, the idea is to reason in terms of time intervals rather than in terms of dates. To explain how this idea works, we shall first recall the regular semantics of AltaRica in the next section. Then, we shall introduce its abstract semantics (in terms of time intervals) in Section 4.

TIMED/STOCHASTIC GUARDED TRANSITIONS SYSTEMS

The semantics of AltaRica 3.0 is defined in terms of stochastic guarded transitions systems [START_REF] Rauzy | Guarded transition systems: a new states/events formalism for reliability studies[END_REF], [START_REF] Batteux | Altarica 3.0 assertions: the why and the wherefore[END_REF]. We shall recall here only the notions that are important for the purpose of this article. The reader should refer to the cited articles for in depth presentations.

Definition

A guarded transitions system is a quintuple V, E, T, A, ι , where:

-V is a set of variables. V is the disjoint union of the set S of state variables and the set F of flow variables: V = S F . Each variable v of V takes its value into a finite or infinite set of constants called the domain of v and denoted as dom(v).

The global state of the system is thus a variable valuation, i.e. a member of the Cartesian product

v∈V dom(v).
-E is a set of events. Each event e of E is associated with a function delay(e) that returns a non-negative real number. delay(e) may be deterministic, in which case it returns always the same value, or stochastic, in which case it returns a value according to a certain cumulative probability distribution.

-T is a set of transitions, i.e. of triples e, G, P , where e is an event of E, G is a Boolean expression built on variables of V (called the guard of the transition) and P is an instruction that modifies the value of state variables (called the action of the transition). For the sake of the clarity, we shall write a transition e, G, P as G e -→ P .

-A is an assertion, i.e. an instruction that modifies the values of flow variables.

ι is a valuation of the variables of V , called the initial state.

Example The GTS encoding the periodically tested components described in the previous section is as follows.

The state of the component is represented by means of two state variables: state that takes its value in {WORKING, FAILED} and phase that takes its value in {OPERATION1, TEST, OPERATION}.

The events are startTest1, startTest, completeTest, failure and repair. They are associated with the delays described in the previous section.

The transitions are as follows.

state = WORKING ∧ phase = TEST failure ----→ state ← FAILED phase = OPERATION1 startTest1 ------→ phase ← TEST phase = OPERATION startTest -----→ phase ← TEST state = WORKING ∧ phase = TEST completeTest -------→ phase ← OPERATION state = FAILED ∧ phase = TEST repair ----→ state ← WORKING, phase ← OPERATION
Finally, the initial state is defined by the variable valuation: state = WORKING, phase = OPERATION1.

Composition

One of the advantages of guarded transitions systems over some other similar formalisms is that they are highly compositional.

Formally, let

M 1 : V 1 , E 1 , T 1 , A 1 , ι 1 and M 2 : V 2 , E 2 , T 2 , A 2 ,
ι 2 be two guarded transitions systems. Then the composition of M 1 and M 2 , denoted as M 1 ⊗ M 2 , is simply the guarded transitions system V, E, T, A, ι such that

V = V 1 ∪ V 2 , E = E 1 ∪ E 2 , T = T 1 ∪ T 2 , A = A 2 • A 1 and ι = ι 2 • ι 1 .
The above principle extends to any number of guarded transitions systems.

Example To represent the system of discussed example in the previous section, it suffices to create two copies of the above guarded transitions system and to compose them (which is done automatically by the AltaRica compiler).

In our example, it may be worth to introduce a flow Boolean variable failed to tell when the system is failed. The assertion defining this variable could be as follows. failed ← A.state = FAILED ∧ B.state = FAILED

Semantics

The semantics of a guarded transitions system S : V, E, T, A, ι is defined as the set of its possible executions.

To define formally the executions, we need to introduce the notion of schedule. A schedule of S : V, E, T, A, ι is a function from T to R + ∪ {+∞}. A schedule Γ is compatible with a state σ of the guarded transitions system and a date d if the following conditions hold for all transitions t : G

e - → P of T . -d ≤ Γ(t) < +∞ if G(σ) = true. -Γ(t) = +∞ if G(σ) = f alse.
Intuitively, an execution of S is a sequence:

σ 0 , d 0 , Γ 0 t 1 -→ σ 1 , d 1 , Γ 1 t 2 -→ . . . tn -→ σ n , d n , Γ n
where n ≥ 0, the σ i 's are states of S, the d i 's are dates i.e. non negative real numbers verifying 0 = d 0 ≤ d 1 ≤ . . . ≤ d n , each Γ i is a schedule compatible with σ i and d i and finally the t i 's are transitions of S.

The set of valid executions is defined recursively as follows.

The empty execution ι, 0, Γ 0 is a valid execution if the schedule Γ 0 is such that for all transitions t : G e -→ P of T :

-Γ 0 (t) = delay(e) if G(ι) = true. -Γ 0 (t) = +∞ if G(ι) = f alse. Now, if Λ = σ 0 , d 0 , Γ 0 t 1 -→ . . . tn -→ σ n , d n , Γ n , n ≥ 0, is a valid execution, then so is the execution Λ t n+1 --→ σ n+1 , d n+1 , Γ n+1 if the following conditions hold, assuming t n+1 = G n+1 e n+1 --→ P n+1 . -G n+1 (σ n ) = true.
σ n+1 = A(P n+1 (σ n )), i.e. the firing of the transition t n+1 is performed in two steps: first, state variables are updated by means of the action P n+1 of the transition, then flow variables are updated by means of the assertion A.

-

d n+1 = Γ n (t n+1 ) and there is no transition t of T such that Γ n (t) < Γ n (t n+1 ).
-Γ n+1 is obtained from Γ n by applying the following rules to all transitions t : G e -→ P of T .

-If G(σ n+1 ) = true, then:

-If G(σ n ) = true and t = t n+1 , then

Γ n+1 (t) = Γ n (t) -Otherwise, Γ n+1 (t) = d n+1 + delay(e)
-If G(σ n+1 ) = f alse, then: This sequence shows how deterministic and stochastic transitions can be intricated. In particular, dates of tests are not decided once for all. They depend on times to failure and to repair of the component.

Γ n+1 (t) = +∞

ABSTRACT SEMANTICS

Principle

The first idea to abstract the executions consists in associating an abstract delay delay with each event of the model. delay (e) is simply the image of the function delay, i.e. an interval of non-negative real numbers.

We have to be a bit careful because some intervals that are the images of distributions are closed while some others are open (to the left and/or to the right) and that we have to consider infinite bounds. A solution consists in working only with closed intervals, but in a non-standard arithmetic built over the set R + = R + ∪ { , ∞}, where and ∞ are respectively infinitely small and infinitely big numbers verifying: 3 gives the abstract delays associated with the most widely used distributions in Al-taRica 3.0. Moreover, a transition whose guard is not satisfied in the current state is scheduled in the interval [∞, ∞].

+ = and ∞ + x = ∞ for all x ∈ R + . In this way, the interval ]a, b[, a, b, ∈ R + , can be encoded as [a + , b -]. Table
(l, h) [l, h] Exponential(λ) [0 + , ∞] Weibull(α, β) [0 + , ∞] Empirical distribution [0 + , ∞]
The second idea is to consider not the date at which transitions are fired, but an interval of time within which they are fired.

Assume that we are building the sequence under study step by step and that the last transition we considered must be fired in the time interval [l, h]. Assume moreover that transitions t 1 , t 2 , . . . t n are scheduled in time intervals

[l 1 , h 1 ], [l 2 , h 2 ], . . . , [l n , h n ].
Then, we can make the following remarks.

1. We must have l ≤ l i for all i = 1, . . . , n, because the next transition cannot be scheduled in the past.

2. We must have also h ≤ h i for all i = 1, . . . , n, because if h i < h for some i, it means that the transition t i must be fired before h i , therefore the last transition must also be fired before h i .

3. For the same reason, we can choose t i as the next transition to be fired only if there is no other transition t j such that h j < l i .

4. Again for the same reason, if the transition t i is fired, it is necessarily fired in the interval [l i , h min ], where h min is the smallest of the h j 's.

5. If the transition t i is fired and the transition t j is such that l j < l i , then l j must be changed to l i so to obey our first remark.

6. Finally, if the transition t i is fired and a transition t associated with the interval [l t , h t ] becomes fireable (t can be the transition t i itself), then t must be scheduled in the interval

[l i , h min ] + [l t , h t ] = [l i + l t , h min + h t ].
We are now able to define formally the abstract semantics of guarded transitions systems (and therefore for AltaRica 3.0).

Formal Definition

The abstract semantics of a guarded transitions system S : V, E, T, A, ι is defined as the set of its possible abstract executions.

To define formally the abstract executions, we need to introduce the notion of abstract schedule. An abstract schedule of S : V, E, T, A, ι is a function from T to closed intervals over R + .

A schedule Γ is compatible with a state σ of the guarded transitions system and the abstract date [l, h] if the following conditions hold for all transitions t : G

e - → P of T , with Γ (t) = [l t , h t ]. -l ≤ l t < ∞ and h ≤ h t if G(σ) = true. -l t = h t = ∞ if G(σ) = f alse.
An abstract execution of S is a sequence:

σ 0 , d 0 , Γ 0 t 1 -→ σ 1 , d 1 , Γ 1 t 2 -→ . . . tn -→ σ n , d n , Γ n
where n ≥ 0, the σ i 's are states of S, the d i 's are abstract dates, i.e. time intervals, each Γ i is an abstract schedule compatible with σ i and d i and finally the t i 's are transitions of S.

The set of valid executions is defined recursively as follows.

The empty abstract execution σ 0 , [0, 0], Γ 0 is a valid abstract execution if the abstract schedule Γ 0 is such that for all transitions t : G e -→ P of T :

-Γ 0 (t) = delay (e) if G(ι) = true. -Γ 0 (t) = [∞, ∞] if G(ι) = f alse. If Λ = σ 0 , [0, 0], Γ 0

Bisimulation

The key mathematical property is that abstract and concrete executions are bisimilar: any concrete (timed, stochastic) execution can be simulated by an abstract execution and reciprocally any abstract execution corresponds to at least one concrete execution. -∀n ≥ 0 ∀t ∈ T Γ n (t) ∈ Γ n (t).

CONCLUSION AND PERSPECTIVES

In this article, we show how the AltaRica 3.0 stepwise simulator has been greatly improved by the introduction of an abstract notion of time. The abstract notion of time enables to reconcile both stochastic and stepwise simulations of AltaRica 3.0 models.

We show that abstract and concrete simulations are bisimilar: any concrete (timed, stochastic) execution can be simulated by an abstract execution and reciprocally any abstract execution corresponds to at least one concrete execution.

We illustrate our purpose using a motivating example that mix both stochastic and deterministic transitions.

The introduction of the abstract notion of time to the stepwise simulator paves the way to the design of efficient model-checking algorithms, and in particular to the design of generators of sequences of events leading to a failure state.

The next step of our work is the application of the presented results for the development of an efficient sequence generator for AltaRica 3.0 models.

Figure 1 :

 1 Figure 1: State diagram for a periodically tested component.

Theorem 1 .

 1 For any concrete execution Λ = σ 0 , d 0 , Γ 0 t 1 -→ . . . tn -→ σ n , d n , Γ n , of the Timed GT S = V, E, T, ι, A, delay it exists an abstract execution Λ a = σ 0 , d 0 , n , d n , Γ n , such that the following properties hold:-∀n ≥ 0 d n ∈ d n ; -∀n ≥ 0 ∀t ∈ T Γ n (t) ∈ Γ n (t). Theorem 2. Any abstract execution Λ a = σ 0 , d 0 , Γ 0 t 1 -→ σ 1 , d 1 , Γ 1 t 2 -→ . . . tn -→ σ n , d n , Γ n of the Timed GT S = V, E, T, ι, A, delay corresponds to at least one concrete execution Λ = σ 0 , d 0 , Γ 0 t 1 -→ . . . tn -→ σ n , d n , Γ n ,such that the following properties hold:-∀n ≥ 0 d n ∈ d n ;

Table 1 :

 1 A possible evolution of component for a periodically tested component.

Table 2 :

 2 Typical values of the parameters

	parameter	A	B
	θ	2190 4380
	π	4380 4380
	τ	0	0
	µ	12	12
	ν	24	24

Table 3 :

 3 Intervals associated with delay functions

	Concrete delay	Abstract delay
	Dirac(t)	[t, t]
	UniformDeviate

t 1 -→ . . . tn -→ σ n , [l n , h n ], Γ n , n ≥ 0, is a valid abstract execution, then so is the abstract execution Λ

h.

-G n+1 (σ n ) = true.

-

-There is no transition t of T such that Γ n (t) = [l, h] and h < l .

-

-Γ n+1 is obtained from Γ n by applying the following rules to all transitions t :

Example We shall consider the abstract version of the execution given in the previous section.


Hayat Zaoui

Pier Luca Palla

Stefano Giordano

Fabrizio Cleri

Maxime Verdier

David Lacroix

J F Robillard

Konstantinos Termentzidis

Évelyne Martin

Thermal conductivity of deca-nanometric patterned Si membranes by multiscale simulations