A Cochard

K Zhu

S Joulié

J Douin

J Huez

L Robbiola

P Sciau

M Brunet

Natural aging on Al-Cu-Mg structural hardening alloys -Investigation of two historical duralumins for aeronautics

Keywords: Aluminum alloys, Aging, Precipitation, Electron microscopy, Mechanical characterization

This article presents the metallurgic study of sixty years old structural hardening Al-Cu-Mg alloys recovered from a double-deck aircraft Breguet 765. This study aims to characterize and understand the evolution of alloys over very long periods thanks to an approach coupling materials characterization and archives researches. The focus is made on two ancient age hardening Al-Cu alloys: A-U4G and A-U4G1, whose compositions are close to the international designations 2017A and 2024 respectively. Multi-scale structural investigation was carry out and mechanical properties were extracted. It is shown that despite natural aging during long periods, these alloys exhibit mechanical properties (yield strength, ultimate tensile strength and rupture strain) close to the expected original specifications.

Introduction

Characterization of the long-term behavior of materials exposed to complex environmental conditions is of great challenge, especially in aeronautics industry where failure can cause irremediable damages. Aluminum alloys are particularly concerned as they constitute the majority of aircraft structures. Developed since the beginning of the 20th century up to now, Al-Cu-Mg alloys or duralumins can be considered as relatively recent metallic materials.

Their good mechanical properties are obtained thanks to age or precipitation hardening. The metallurgical process behind hardening consists in heat treating near solvus temperature the Al-Cu-Mg alloy, quenching it in water then performing an aging treatment. The nanoscale precipitation responsible for the hardening of the material occurs at this step. Aging was done originally by maintaining the alloy at room temperature for several days (natural aging or equivalent T3 and T4 treatments). Other treatments can now be found consisting in heating the products between 150 and 200°C for several hours up to tens of hours (artificial aging or T6 and T8 treatments) to reach the maximum hardness and tensile strength, i.e. peak aging [START_REF] Starke | Heat-treatable aluminum alloys[END_REF]. The precipitation sequences in various alloys and the influence of the process conditions on the microstructure and the mechanical properties are well documented [START_REF]Handbook of Aluminum, Physical Metallurgy and Processes[END_REF].

However, the long-term behavior of these alloys when exposed to real environmental conditions during the lifespan of an aircraft remains to be fully understood. Accelerated tests such as creep and high-temperature fatigue or exposure to corrosive environment have been developed to identify degradation mechanisms and consequently to assess long-term aging [START_REF]Accelerated Aging of Materials and Structures: The Effects of Long-Term Elevated-Temperature Exposure[END_REF]. Elevated temperature exposure under applied stress is known to introduce changes in the microstructure such as growth of new phases, changes in dislocation density, distribution or grain-boundary precipitation, which influence largely the mechanical properties [START_REF]Accelerated Aging of Materials and Structures: The Effects of Long-Term Elevated-Temperature Exposure[END_REF][START_REF] Wiersbinska | The influence of long-lasting annealing on microstructure of AlCu4Ni2Mg2 alloy[END_REF][START_REF] Moy | Influence of heat treatment on the microstructure, texture and formability of 2024 aluminium alloy[END_REF][START_REF] Larignon | Mécanismes d'endommagement par corrosion et vieillissement microstructural d'éléments de structure d'aéronef en alliage d'aluminium[END_REF]. In order to be closer to reality, application of low temperature ranges and long times is a first approach to model long-term aging. For instance, studies were published where heat was applied on age-hardening alloys for up to two months at temperatures between 100 and 150°C [START_REF] Larignon | Mécanismes d'endommagement par corrosion et vieillissement microstructural d'éléments de structure d'aéronef en alliage d'aluminium[END_REF] or up to 8 months (5800h) at 85°C [START_REF] Braun | Mechanical properties and corrosion behaviour of 2024-T351 sheet thermally exposed at slightly elevated temperatures[END_REF]. In these two cases, the effect on tensile properties was proved not significant. Accelerated corrosion exposure tests have also been applied for these alloys [START_REF] Larignon | Mécanismes d'endommagement par corrosion et vieillissement microstructural d'éléments de structure d'aéronef en alliage d'aluminium[END_REF][START_REF] Braun | Mechanical properties and corrosion behaviour of 2024-T351 sheet thermally exposed at slightly elevated temperatures[END_REF][START_REF] Alexopoulos | Experimental and theoretical studies of corrosion-induced mechanical properties degradation of aircraft 2024 aluminum alloy[END_REF][START_REF] Guillaumin | Localized corrosion of 2024 T351 aluminium alloy in chloride media[END_REF] to study the impact on the microstructure and the residual mechanical properties.

Although fully justified and helpful for assessing the materials behavior requested by the aircraft industry, accelerated tests, which simulate long-term service usage, do not represent the exact real aging [START_REF]Accelerated Aging of Materials and Structures: The Effects of Long-Term Elevated-Temperature Exposure[END_REF][START_REF] Dedamborenea | Comparison of accelerated and atmospheric exposure tests for corrosion of aluminium alloys[END_REF]. They also take insufficiently into account the synergistic effects corresponding to complex conditions encountered during aircraft service, which include indeed many factors such as temperature, pressure, loads, pollution or relative humidity. All the combined factors lead to fatigue damage accumulation, which is hard to predict.

This work aims at characterizing aeronautic materials naturally aged on an aircraft. More precisely, two Al-Cu-Mg industrial alloys were sampled on an airplane Breguet dating back from the 1950's. These materials, corresponding to the French AFNOR designation A-U4G and A-U4G1 used at the time, are typical structural hardening alloys used in aeronautical industry at this period in Western Europe. A-U4G would be close to the 2017 alloy according to the international designation system in use nowadays [START_REF]The Aluminum Association, International alloy designations and chemical composition limits for wrought aluminum and wrought aluminum alloys[END_REF]. A-U4G1 would be close to the current 2024 alloy and still employed for some parts of aircrafts fuselage skin. They were submitted to the lifespan of the aircraft and outdoor exposure during several decades.

Thanks to the restoration of the aircraft, the opportunity to collect specific pieces was offered.

These materials represent models for studying long-term aging in real conditions.

Information on mechanical properties durability and microstructural evolution that includes very slow physical mechanisms may be extracted. This investigation of industrial heritage samples is a direct continuation of a previous study focused on aluminum rivets [START_REF] Cochard | Benefits of the complementary use of archaeometry investigations and historical research in the study of ancient airplanes: the Breguet Sahara's rivets[END_REF]. It requires a detailed methodology coupling materials science and archive historical researches. Indeed, original state of the aged materials cannot be measured directly and the only way to study it is, indirectly, through information collected in technical documents contemporary of the aircraft construction. Comparison with modern alloys can be made but with care as composition may have evolved over the years (impurity levels for instance), which can affect their microstructure and their properties.

After the presentation of the method integrating the sampling strategy, the characterization of materials will be reported. Multi-scale structural investigation, as well as mechanical properties based on tension and microhardness testing have been carried out. The results are compared to data obtained from ancient scientific literature and from technical (R&D) archives of the airplane constructor. The discussion section proposes some hypothesis on the aging state of the sampled alloys.

Materials and methods

Materials

The studied materials were sampled on a 60 years old Breguet Sahara 765 "Deux-Ponts" (Double-Decker) under renovation at the association "Ailes Anciennes Toulouse" in Blagnac, France. The plane flew from 1959 until 1969 [START_REF] Cuny | Les avions Breguet[END_REF], after which it was grounded, first on the military base of Evreux, France, then in Blagnac. The aircraft thus experienced two periods of lifetime: 10 years of air service and about 47 years at rest in outdoor condition. Such a long period in outdoor conditions gave rise to severe damages including corrosion on the skin and structural elements. During one of the restoration operations involving metal replacements, two original metallic plates were selected and collected after their removal from this ancient aircraft (Fig. 1). Both plates have a thickness close to 2 mm. They are rolled products, commonly used for fuselage skin and stringers [START_REF] Starke | Application of modern aluminum alloys to aircraft[END_REF]. The first one, called L8 is a large flame-shield plate located inside the engine nacelle. The second one, called L2, is a stringer, located underneath the floor, and acts as a reinforcement at the back of the aircraft.

Methods

Investigation of these industrial heritage materials requires an adapted methodology, which involves at the same time archives research to extract information on original pieces (nature, role, expected properties) and classical metallurgical study on sampled materials hereafter described.

Chemical composition

The two plates L2 and L8 were first analyzed to determine the nominal composition of the alloys. Chemical analysis was performed by Inductively Coupled Plasma -Optical Emission Spectrometry (ICP-OES) by Evans Analytical Group SAS. Table 1 gives the measured elemental composition of each plate and for comparison, the expected compositions for the same alloys from the French standard AIR 3350/C of 1957 [START_REF]Tableaux des alliages d'aluminium et des alliages de magnésium, norme AIR 3350/C, édition n°4 du 15 juin 1957[END_REF]. This standard was applied by the Breguet Company at this period, as it is confirmed in the standards for this aircraft [START_REF]feuille 1, Usine d'Anglet[END_REF]. were indeed planned for the fuselage skin and structure [START_REF]Service Technique de l'Aéronautique -section Avions, Clauses Techniques Breguet 765[END_REF]. It must be recalled that these two alloys, with their high specific strength properties, belong to the Duralumin family. They are age hardened aeronautic alloys, commercialized since 1909 for A-U4G [START_REF] Duparc | Alfred Wilm et les débuts de Duralumin[END_REF] and at the beginning of the 1930s for A-U4G1 [19]. The latter contains a higher magnesium content to improve the tensile strength. A-U4G and A-U4G1 have a composition close to the currently available commercial alloys designated respectively as 2017A and 2024, but with a difference in impurities concentration (Si and Fe). Additional analyses performed by Glow Discharge Mass Spectrometry (GDMS) revealed that the L2 alloy also contains markedly less impurities elements (Zn, Cr, Ti, Ni, Pb, Ga) than the L8 plate.

Structural characterization

In order to observe the metallurgical structure including the different phases distribution, the plates were cut along three perpendicular directions. The directions of cut were: (D1) expected rolling direction, (D2) transversal to rolling surface direction and (D3) thickness direction corresponding to usual Normal Direction (ND). Rolling direction (RD) was hypothesized after preliminary optical observation of cross-sections.

The microstructure was determined from optical and scanning electron microscopies. The plate samples were prepared as usually for metallographic examination. Each specimen was embedded in an epoxy resin and mechanically polished on water-lubricated abrasive papers (silicon carbide) of decreasing granulometry and then on polishing cloths with diamond pastes down to 1 µm particle size. Optical observation was conducted on an inversed microscope Nikon MA200. Etching with Keller's solution (HF 1 ml, HNO3 2.5 ml, HCl 1.5 ml in demineralized H2O 95 ml) was applied revealing the microstructure of the alloy and allowing grain shape and size to be measured with linear intercept method.

Examination with scanning electron microscope (SEM) on a JEOL JSM 6460 LV operating at 20 kV and equipped with an Oxford INCA PentaFETx3 energy dispersive X-ray spectrometer (EDS) was performed. Large precipitates (> 2 µm²) were identified and quantified using pure K-line Al, Si and Fe standards. Results were normalized to 100 wt%.

In order to observe the nanostructure and precipitation at this scale, the plates were also cut according to the three directions (D1, D2, D3). After mechanical polishing up to 2400 SiC paper grade, in order to obtain a thickness of 25 µm, the specimens were then electrochemically thinned using a Tenupol-5 Struers apparatus operating at 60V in a solution of methanol and nitric acid (3:1) at -15°C. The nanostructure was observed using a Transmission Electron Microscopy JEOL 2010 operating at 200kV. Chemicals elements presents in some precipitates were identified with a CM20 FEG TEM/STEM microscope, operating at 200 kV and equipped with an energy dispersive X-ray spectroscopy (Quantax EDS system with silicon drift detector, from Bruker).

Mechanical testing

Tensile specimens were cut in each plate: 6 in D1 direction and 6 others in D2 direction (perpendicular to D1) allowing to get a representative set of data in order to evaluate a possible dispersion of the results for these two surface directions. As depicted in Fig. 2, the gauge length/reduced section and width of the tensile specimen are 20 mm and 3 mm respectively. The gauges were polished mechanically up to SiC 1000 grade in order to eliminate the painting and in some cases some irregular surface flatnesses. The thickness (direction D3) of the tensile test specimen is in the range of 1.2-1.5 mm. The tensile test was performed at a strain rate of 0.02 mm.s -1 , on Synergie apparatus using a MTS® extensometer.

For each plate, 8 gauges were measured.

Fig. 2. Shape and dimension (mm) of the tensile specimen.

Hardness

Microhardness measurements were carried out on each side of each plate using a Buehler Omnimet serie 2100, applying a load of 300g during 10s on polished surfaces. Results are given in Hardness Vickers (HV) numbers.

3 Results It is revealed that the elongation ratio defined as length/width is higher than 2 for (D1, D3) and (D2, D3) planes in both plates, but more pronounced for L2 than for L8.

Metallurgical characterization

For L8, this may result from a cross-rolling of the plate, conducting to similar size grains in the (D1, D2) plane. This process involves changing the rolling direction (usually 90°) after each pass [START_REF] Rout | Cross Rolling: A Metal Forming Process[END_REF]. It reduces grain directional dependency of mechanical properties, residual stress and improve the formability of the metal. However, the use of such a process in 1959 by Breguet Compagny remains to be proven. For L2, the rolling direction corresponds most probably to the D2 direction for which a marked elongation of the grains is evidenced.

Furthermore, as given in Table 2, the average size values of the L8 grains are smaller by a factor 2-3 than those in L2 plate. The maximum elongation is about 45 µm for L8, while it is 114 µm for L2. From the hardness measurements also reported in Table 2, the values of L8 and L2 plates are about 132 HV and 150 HV respectively. There is not an important difference between the different planes excepted that the external surface corresponding to rolling plane (D1, D2) exhibits lower hardness values.

Thus when comparing L8 with L2 plates, the L8 cross-rolled plate exhibits a lower grain size and a lower hardness. The difference in microstructure is linked to the mechanical deformation and to the thermal treatment applied. But the difference of mechanical properties (hardness) finds mainly its origin in the difference of composition in addition elements such as Mg in L2 plate (Table 1), as it will be discussed in analysis and precipitations section.

SEM observation and EDS analysis -Constituent particles and equilibrium phases

Fig. 4 shows typical coarse intermetallic particles detected from backscattered electron image with SEM and quantified with EDS. The two alloys reveal a similar density of micro-particles (around 15,000 particles/mm²) with a particle size ranging from 1 µm up to 18 µm in length. Among the coarse intermetallic particles, the distinction between the insoluble constituent precipitates and the equilibrium precipitates is possible.

In general, constituent particles are the coarsest and first particles formed during solidification of the aluminum alloy ingot [START_REF] Wang | Precipitates and intermetallic phases in precipitation hardening Al-Cu-Mg-(Li) based alloys[END_REF]. They contain iron and/or silicon impurities.

Their size depends on the solidification rate, the content in iron and silicon and the manufacturing process. Mostly insoluble during further heat treatments, they do not participate in the hardening process of the alloy. For the L8 plate (Fig. 4), the fractioned and irregularly shaped grey particles correspond to these constituent particles. Al, Cu, Fe, Mn and Si with various proportions were detected by EDS, which could indicate that they correspond to α-AlMnSi or α-Al(Fe,Mn)Si phases.

On L2 plate, the same types of particles (α-Al(Fe,Mn)Si) were found. Another type of particle, rich in Fe and Cu, was also detected that reveals the presence of Al7Cu2Fe phase. All these phases are expected insoluble constituent particles for 2024 or A-U4G1 alloys [START_REF] Davis | Aluminum and aluminum alloys[END_REF].

In Fig. 4, round and bright particles observed in the L8 plate can be readily identified as θ-Al2Cu equilibrium precipitates. On the L2 plate, similar looking particles appear to be also θ-Al2Cu and S-Al2CuMg equilibrium precipitates. These two types of precipitates form during aging after a solution heat treatment and quenching. The possible precipitation sequences for Al-Cu-Mg alloys follow the following steps [START_REF] Starke | Heat-treatable aluminum alloys[END_REF][START_REF] Dubost | Durcissement par précipitation des alliages d'aluminium[END_REF]:

Supersaturated solid solution (SSS) → GP zone → θ'' → θ' → θ-Al2Cu, → GPB zone → S'' → S'→ S-Al2CuMg
with GP and GPB being Guinier-Preston and Guinier-Preston-Bagaryatsky respectively. θ' and S' are metastable phases forming upon artificial aging and whose presence leads to the maximum hardening state. Equilibrium precipitates (θ-Al2Cu and S-Al2CuMg) are larger and more widely spaced and form by transformation or coalescence of metastable precipitates.

The latter precipitation sequence leading to S-Al2CuMg was of controversial matter until recently [START_REF] Moy | Influence of heat treatment on the microstructure, texture and formability of 2024 aluminium alloy[END_REF][START_REF] Wang | Two types of S phase precipitates in Al-Cu-Mg alloys[END_REF][START_REF] Styles | The coexistence of two S (Al2CuMg) phases in Al-Cu-Mg alloys[END_REF], particularly about the crystallographic structure of the lath shape S'' and S' phases. It is now believed that these two phases have a crystalline structure identical to the equilibrium S phase but with varying lattice dimensions [START_REF] Sha | Nanostructure of aluminium alloy 2024: Segregation, clustering and precipitation processes[END_REF]. The lath-shape metastable precipitates will be named S' here for distinction with coarser particles detected by SEM.

In general, for an industrial Al-Cu-Mg alloy, such as ancient A-U4G and A-U4G1 alloys, the copper is in excess compared to the magnesium. Hence, the sequence of precipitation depends on the Cu/Mg ratio [START_REF] Mondolfo | Aluminum Alloys, Structure and Properties[END_REF]. For Cu/Mg ratio between 4 and 8, the precipitation sequence leads to the formation of both θ-Al2Cu and S-Al2CuMg equilibrium phases.

Precipitation sequence leading to S-Al2CuMg is known to provide a higher hardening at room temperature than the one leading to θ-Al2Cu. However, the sequence of precipitation depends also on the Mg/Si ratio [START_REF] Wang | Precipitates and intermetallic phases in precipitation hardening Al-Cu-Mg-(Li) based alloys[END_REF][START_REF] Mondolfo | Aluminum Alloys, Structure and Properties[END_REF]. When there is an addition of silicon in alloys, hardening precipitates such as Mg2Si (β) or Al5Cu2Mg8Si6 (Q phase) can also appear during the hardening sequence [START_REF] Mondolfo | Aluminum Alloys, Structure and Properties[END_REF].

Thus the difference observed between the L2 and L8 plates has to be directly linked to their different compositions (Table 1). For the L2 plate, the Cu/Mg ratio is equal to 2.9, which means that Al2CuMg precipitates should be the principal hardening agent. The presence of θ-Al2Cu precipitates and black precipitates rich in Al, Mg and/or Si indicates that the hardening process followed in this case several precipitation sequences.

The L8 plate contains less than 1 wt% of Mg and 0.6 wt% of silicon (against 0.167 wt.% in L2), which corresponds to a Cu/Mg ratio equal to 4.7 and a Mg/Si ratio equal to 1.4. The precipitation in L8 is thus different than in L2. Because of the low Mg/Si ratio (<1.7) in L8, S-Al2CuMg precipitation sequence is not favorable: only ϴ-Al2Cu and Al5Cu2Mg8Si6 should be formed. In the L8 plate, the black particles observed by SEM, with a composition of Al(Cu)MgSi could correspond to this phase. As regards to the submicronic scale, different rod-like precipitates are evidenced for L2 and L8 plates. For the L2 plate, having a composition close to 2024 alloy, the precipitates have a regular rod-like shape with a length in the range of 30-290 nm, mainly composed of Al, Cu and Mn and homogeneously distributed throughout the matrix. These particles, also called dispersoids, are formed during a long-term heat treatment, such as the ingot homogenization and have thereafter a low solubility in solid aluminum through some of their constituent elements: Cr, Mn and/or Zr. Their role is to control grains size and resistance to recrystallization avoiding large grains formation and favoring grain boundaries to pin dislocations [START_REF] Wang | Precipitates and intermetallic phases in precipitation hardening Al-Cu-Mg-(Li) based alloys[END_REF]. For an alloy of 2024 type with a typical weight Cu/Mg ratio equal to 2.8, the main dispersoids expected are Al20Cu2Mn3, so-called phase T [START_REF] Wang | Precipitates and intermetallic phases in precipitation hardening Al-Cu-Mg-(Li) based alloys[END_REF][START_REF] Shen | The structure determination of Al20Cu2Mn3 by near atomic resolution chemical mapping[END_REF].

In the L8 plate, ovoid and rod-like shape particles with length and width varying from 40 to 200 nm are observed, thus smaller than the dispersoids detected in L2, and sometimes aligned along perpendicular axes. With TEM-EDS, it is noted that the composition of particles is quite heterogeneous: AlMnFeSi, AlMnSi, AlCuFe and AlCuMgSi can be found. In this case, the dispersoids are supposed to be the ones containing Mn: of the type AlMnSi or AlMnFeSi. As shown on Fig. 5(a), precipitation at grain boundaries (inset) is observed. SiO2 faceted particles are also quite often found due to high silicon rate in this A-U4G alloy.

As regards to the nanometric level, in L8 and L2 plates, smaller particles with round shape in the range of 5 and 20 nm and heterogeneous distribution were observed (Fig. 6).

Fig. 6. Small-size particles and dislocations observed by TEM on L2 with EDS mapping showing the copper rich nature of particles.

Composed of mainly Cu, these small size intermediate particles participate in hardening.

Dislocations are indeed pinned by these small size particles, which proves that they contribute to the hardening of the material by hindering the movement of dislocations. These type of small particles observed in both plates can also be found in equivalent modern 2024 alloys with a T3 temper [START_REF] Moy | Influence of heat treatment on the microstructure, texture and formability of 2024 aluminium alloy[END_REF][START_REF] Lin | Precipitation hardening of 2024-T3 aluminum alloy during creep aging[END_REF].

In L2 piece (A-U4G1), whose precipitation is supposed to be conducted mainly via S-type phase, the hardening S'-phase in the form of laths was not found. As explained in [START_REF] Starke | Heat-treatable aluminum alloys[END_REF][START_REF] Lin | Precipitation hardening of 2024-T3 aluminum alloy during creep aging[END_REF], for wrought alloys tempered at room temperature, corresponding to T4 in modern designation of treatments [START_REF] Starke | Application of modern aluminum alloys to aircraft[END_REF], precipitation is exclusively GPB, which is not observable with conventional TEM. Semi-coherent and coherent precipitates, metastable S'-phase, with a lath shape appear generally upon artificial aging (corresponding to T6 or T8 treatment) as shown in [START_REF] Starke | Application of modern aluminum alloys to aircraft[END_REF][START_REF] Sha | Nanostructure of aluminium alloy 2024: Segregation, clustering and precipitation processes[END_REF][START_REF] Lin | Precipitation hardening of 2024-T3 aluminum alloy during creep aging[END_REF][START_REF] Alani | In situ transmission electron microscopy studies and real-time digital imaging[END_REF].

Table 3 summarizes the various particles found in the two plates at different observation scales, including a comparison with known precipitates from literature.

Dislocations

An important feature concerns dislocations microstructure and localization as they could give indication on the mechanical and heat history experienced by each piece. A lot of dislocations were observed by TEM in the two plates, as shown in Fig. 7. In the L2 plate, the dislocations surround the dispersoids. In the L8 plate, however, the dislocations distribution and type are not exactly the same: helical type dislocations loosely distributed in the matrix are observed. This generally happens when heat is applied and dislocation climb is activated [START_REF] Friedel | Dislocations[END_REF], and could be a sign of aging linked to the flight period of the aircraft.

Mechanical characterization

Mechanical tests performed on L8 and L2 plates revealed a very good reproducibility of the tensile strength curves. This is observed both on the D1 direction as well as on the transversal one D2.

A typical curve for each rolling directions is given in Fig. 8. The measured mechanical properties are given in Table 4, together with data from archives discussed hereafter. For the two plates, the Young modulus and the rupture strain are equivalent: around 67 GPa and 16% respectively. However, an important difference is found concerning the strength properties. It is higher for L2 than for L8. This agrees well with the microstructure patterns and hardness values previously evidenced. As a consequence, the alloy of the L2 plate exhibits higher mechanical properties than that of L8. Although the grain size is more important for L2 than for L8, the reason for higher mechanical properties comes from the higher Mg rate, which is known to increase the tensile strength and hardness of Al-Cu alloys [START_REF] Dubost | Durcissement par précipitation des alliages d'aluminium[END_REF].

Discussion

Comparison of mechanical properties: after natural aging vs original informative data

The mechanical properties extracted from measurements of L8 and L2 plates are summarized in Table 4.

They were compared to historical data extracted from archives documents. Table 4 also presents the standard values recommended by Breguet's aircraft company in 1957, together with a range of values for the same alloys reported in 1954 by a generalist French scientific and technical journal: la Revue de l'Aluminium [START_REF] Chevigny | L'aluminium et ses alliages[END_REF]. Current standards data for commercial 2017A-T4 and 2024-T4 alloys are also given for reference [START_REF] Kaufman | Properties of Aluminum Alloys, Fatigue Data and the Effects of Temperature, Product Form and Processing[END_REF]. The first observation is that the measured yield strength, ultimate strength and rupture strain for L8 (A-U4G) and L2 (A-U4G1) are consistent with the company's standards [START_REF]feuille 1, Usine d'Anglet[END_REF] and fall in the expected range for applications broader than aeronautics. The materials recovered from the Breguet 765 are close to current standards for 2017A-T4 and 2024-T4.

Although the original state of these particular parts remains unknown, one can say that after 10 years of flight service and 47 years of outdoor conditions, uncorroded zones of the alloys in use in the Breguet Sahara kept the required mechanical properties.

These results complement studies reported in the literature on naturally aged alloys directly sampled from aircrafts at the end of their usage life. Salimon et al. [START_REF] Salimon | The evolution of electrochemical, microstructural, and mechanical properties of aluminium alloy 2024-T4 (D16AT) during fatigue cycling[END_REF] proved indeed that bulk material 2024-T4 remained unaffected mechanically by the high cycle fatigue experienced during 11-20 years of service. Changes in chemical composition and in the structure were confined to a thin surface layer making the material more sensitive to corrosion.

However, in our study a difference is noticed concerning the hardness. The values of 132±7

HV for L8 (A-U4G) and 150±7 HV for L2 (A-U4G1) are significantly higher than the recommended and standards values, respectively 117 HV and 129 HV. This marked difference could be related to a work hardening process occurring during the lifetime of the aircraft.

4.2. Effect of long-term aging -structural correlation.

For both alloys, A-U4G1 (L2) and A-U4G (L8), most microscopic constituent particles, mesoscopic dispersoids and nanometric hardening precipitates were identified.

In L2 sheet, they correspond, in nature, shape and distribution to what is expected for 2024-T4 [START_REF] Starke | Heat-treatable aluminum alloys[END_REF][START_REF] Moy | Influence of heat treatment on the microstructure, texture and formability of 2024 aluminium alloy[END_REF][START_REF] Lin | Precipitation hardening of 2024-T3 aluminum alloy during creep aging[END_REF]. Treatment T4 is identical to the treatment described in Breguet's standards [START_REF]feuille 1, Usine d'Anglet[END_REF]:

tempered for 4 days at room temperature. The similarity in matrix precipitates in L2 with literature data is consistent with the mechanical properties here measured.

For A-U4G (L8 plate), precipitates nature and distribution were also identified. We noticed a very heterogeneous distribution of compositions for dispersoid-type particles, which is probably linked to the higher impurities level (Si and Fe mainly) in L8 compared to L2 plate.

It is however difficult to make rational comparisons with literature. The modern equivalent alloys (2017A) are not as thoroughly observed with TEM as 2024 alloys, and as regards to scientific literature from the 50s or 60s, there was no similar study to be found as electron microscopy was not at that time being extensively applied for industrial alloys. Saying that, some peculiar features were observed on L8 plate, as detailed in Fig. 9: precipitation of Mgrich particles on dispersoids (AlMnSi type) which act as nucleation sites. This can appear upon heating [START_REF] Wiersbinska | The influence of long-lasting annealing on microstructure of AlCu4Ni2Mg2 alloy[END_REF][START_REF] Alani | In situ transmission electron microscopy studies and real-time digital imaging[END_REF], and it is worth recalling that L8 was a flame-shield plate located inside the engine nacelle, thus subjected to higher temperatures than L2. To validate this point, additional investigation was performed on another A-U4G piece also sampled on the same airplane, but which didn't experience any heating from the engine. This Al-sample is a flat collected on the inner desk part of the Breguet 765. In this sample, the re-precipitation phenomenon was not observed as shown in Fig. 10. Furthermore, for the L8 plate, several other observations at the nanoscale indicate that the plate experienced some heating. In particular, helical-type dislocations have a loose distribution (Fig. 7), compared to the case of L2 where dislocations are pinned by dispersoids. Another peculiar feature is the numerous precipitates with elongated shape at grain boundaries as shown in Fig. 5(a). The temperature and the exposure time experienced by the L8 plate are unknown but some assumptions can be made. The temperature probably didn't exceed 150°C since no needle shape (θ'-Al2Cu) or lath shape (S'-phase) precipitation occurred. On the other hand, it is highly possible that the plate experienced temperatures higher than 80°C as microstructural changes were detected. In [START_REF] Larignon | Mécanismes d'endommagement par corrosion et vieillissement microstructural d'éléments de structure d'aéronef en alliage d'aluminium[END_REF] indeed, no microstructural changes were observed in 2024 alloys for exposure times equivalent to about 100,000 hours at 80°C, which is a much longer time than the lifetime of the Breguet Sahara. It can be estimated indeed that the Breguet 765 Sahara n°504 flew a maximum of 36,500 hours based on the hypothesis that the flight-time was 10 hours per day for 10 years of service. In the present case, observed microstructural changes occurred due to exposure to temperature higher than 80°C.

However, despite the nanoscale changes in microstructure here evidenced, the mechanical properties for L8 alloy remain adequate to the initial company's specifications. They might have evolved with time. For proving this assumption, one would need to measure the mechanical properties of the original material.

To explain the higher hardness in both alloys, further investigation should be conducted, in particular quantitative studies on precipitates size and on dislocation density.

Conclusion

Even if aluminum alloys and related fabrication processes were developed relatively recently in industry, parts recovered from old aircrafts can provide a lot of information useful for the field of material science. Two vintage aluminum alloys, A-U4G and A-U4G1 in the form of plates were recovered from a Breguet 765. Multi-scale observations of their microstructure and mechanical tests allowed the understanding of their nature and fabrication process. It was indeed possible to confirm:

-The nature of the alloy thanks to chemical analysis but also thanks to TEM observation pointing out specific precipitation phases.

But also several fabrication steps:

-The mechanical treatment, i.e. rolling, thanks to the observation of grain size and orientation.

-The age-hardening treatment thanks to the observation of the presence or absence of hardening precipitates. S'-type hardening precipitates in the form of laths were not observed which confirms that the alloys were matured (T3 or T4 temper) and not artificially aged (T6 or T8 temper). Nano-size copper-rich precipitates pinning the dislocations are similar to ones recently observed in 2024-T3 alloys.

Most interesting, the alloys were naturally aged for about 60 years: 10 years experiencing flights cycles of the aircraft and more than 47 years outdoors without flying. It is shown that both alloys kept good mechanical properties, conformed to standards and most probably close to their original state. Hardness was however significantly increased and this remains to be explained. Some structural changes were observed at the micro-and nano-structure levels for one of the plate (L8 -A-U4G alloy). They are related to the conditions experienced by the plate during the lifetime of the aircraft. This plate located close to the engine most probably experienced some increased temperature. But despite this service aging (temperature cycles, mechanical strain for 10 years) and although signs were detected at the nanoscale, there is no evident change in the mechanical properties at macroscale. It would be interesting to measure the mechanical properties and especially the hardness of the same alloys in a state close to the original state (such as spare-parts kept for this aircraft) to determine a more precise aging of the material.

The presented results show that studies of naturally-aged alloys could thus be considered to complete artificial aging tests for evaluating the long-term aging of structures. Future experiments will consist in comparing the microstructure of these ancient alloys to artificially aged modern alloys (such as recent 2017A-T4 and 2024-T4).

Fig. 1 .

 1 Fig.1. (a) Presentation of the Breguet 765 aircraft, source: ©Olivier Beernaert at http://www.traditions-air.fr; and the two samples location (b) L8, a plate located in the left engine nacelle; (c) L2, a plate located underneath the first stage, into the plane.

3. 1 . 1

 11 Metallurgical structureObservation performed after chemical etching, revealed a typical structure linked to mechanical deformation. As shown in Fig.3, the grains exhibit elongation both in the (D1, D3) and (D2, D3) planes. The values of grain sizes are given in Table2. They exhibit a large variation around the average values.

Fig. 3 .

 3 Fig. 3. Metallographic observation in the (D1, D2, D3) space -structure of (a) L8 and (b) L2 plates after etching (Keller solution).

Fig. 4 .

 4 Fig. 4. Constituent particles and equilibrium phases identified by SEM with backscattering detector and EDS, (a) on L8 plate: θ-Al2Cu, Al(Cu,Fe,Mn,Si), Al(Mg,Si,Cu); (b) on L2 plate: θ-Al2Cu, S-Al2CuMg, Al(Cu,Fe,Mn,Si), Al(Mg,Si,Cu)

3. 1 . 3

 13 Fig. 5 details typical submicronic particles identified with TEM and EDS. At this investigation scale, the two studied alloys reveal different particle shapes and sizes.

Fig. 5 .

 5 Fig. 5. Dispersoids observed by TEM for samples in the rolling plane (D1,D2), (a) on L8 plate, with grain boundary precipitation (inset) and SiO2 faceted particle, (b) on L2 plate

Fig. 7 .

 7 Fig. 7. Dislocations observed by TEM on (a) L8 and (b) L2 plates.

Fig. 8 .

 8 Fig. 8. Tensile strength curves for L8 and L2 plates. Full lines are for specimen parallel to the rolling direction D1 and dashed line are for the other rolling direction D2.

Fig. 9 .

 9 Fig. 9. TEM-EDS mapping on L8 plate showing the nature of precipitates nucleating on intermetallic type particles.

Fig. 10 .

 10 Fig. 10. TEM images of dispersoids (a) on L8 and (b) on unheated A-U4G piece from the Breguet 765 Sahara n°504 showing that re-precipitation observed is specific to L8.

Table 1

 1

	Elemental composition in wt% of sampled aeronautic plates by ICP-OES compared to 1954 literature
	data, standard Air 3350/C [15] and current international standard by The Aluminum Association [11]
	Wt%	Al	Cu	Mg	Mn	Si	Fe	Ti
	L8	95.57	3.88	0.83	0.522	0.60	0.53	0.063
		±1.14	±0.05	±0.01	±0.007	±0.01	±0.04	±0.001
	A-U4G [15]	base	4	0.7	0.5	0.5	0.5 max	0.2 max
	2017A [11]	base	3.5-4.5	0.4-0.8	0.4-1.0	0.2-0.8	0.7	0.15
	L2	94.30	4.36	1.53	0.501	0.167	0.25	0.0097
		±0.41	±0.03	±0.02	±0.006	±0.001	±0.01	±0.0002
	A-U4G1 [15]	base	4.25	1.5	0.7	0.4 max	0.5 max	0.2 max
	2024 [11]	base	3.7-4.5	1.2-1.5	0.15-0.8	0.15	0.20	0.15
	L8 plate contains less magnesium (0.83 wt%) than L2 (1.53 wt%), which shows that L8 is
	made of A-U4G alloy whereas L2 is made of A-U4G1 alloy. For the series of Breguet 765
	Sahara aircrafts, technical clauses and recovered plans indicate that A-U4G and A-U4G1

Table 2

 2

	Grain sizes and hardness of L8 and L2 plates according to the (D1, D2, D3) space.
	Plate	Microstructure	Grain size	Ratio	Number of	Hardness	Number of
	number	-Global shape	(average	L/W	measurements	(HV)	measurements
	/		± std dev.)
	section
			Length	Width
			(μm)	(μm)
	Plate L8	(average)					132 ± 7	75
	(D1, D2)	Equiaxial	35±9	30±9	1.2	50	124±3	20
	external
	surface
	(D1, D3)	Elongated	46±15 15±5	3.1	50	141±8	15
	(D2, D3)	Elongated	44 ±24 19±8	2.3	57	132±5	40
	Plate L2	(average)					150 ± 7	75
	(D1, D2)	Elongated	98±44 64±28 1.5	40	147±4	33
	(D1, D3)	Elongated	69±35 23±8	3	40	161±7	15
	(D2, D3)	Elongated	114±53 29±12 3.9	40	149±5	27

Table 3

 3

	Precipitates identified by SEM-EDX and TEM-EDX within sampled aeronautic L8 and L2 plates.
		L8	L2	Literature known
		A-U4G	A-U4G1	precipitates
	Hardening phase	AlCu	AlCu	θ'-Al2Cu [22]
	Equilibrium phase	AlCu	AlCu	θ-Al2Cu [22]
			AlCuMg	S-Al2CuMg [22]
	Constituent phases	Al(Cu)MgSi	Al(Cu)MgSi	Al5Cu2Mg8Si6 [26]
			AlCuFe	Al7Cu2Fe [21]
		Al(Cu, Fe, Mn, Si)	Al(Cu, Fe, Mn, Si)	Al12(Fe, Mn)3Si [21]
	Dispersoids	AlMn(Fe)Si	AlCuMn	Al20Cu2Mn3 [27]

Table 4

 4 Mechanical properties determined from tensile tests of L8 and L2 plates and corresponding Al-Cu-Mg alloys properties from Breguet Standards (1954)[START_REF]feuille 1, Usine d'Anglet[END_REF], Revue de l'Aluminium (1954)[START_REF] Alani | In situ transmission electron microscopy studies and real-time digital imaging[END_REF] and recent handbook[START_REF] Friedel | Dislocations[END_REF].

			Young	Yield strength	Ultimate	Rupture	Hardness
			modulus	Kg/mm² (MPa)	strength	strain	HB*/HV
			Kg/mm²		kg/mm²
			(GPa)		(MPa)
	Measurements (this work)	L8 plate	6830 (67)	27 (264)	43 (421)	17%	113/132
	Breguet Standards (1957) [16]	A-U4G days (sheet) 495°C Matured for 4	7000 (69)	27 (265)	40 (392)	15%	100/117
	Revue de l'Aluminium (1954) [30]	A-U4G	--	26-45 (255-441)	43-54 (422-529)	11-25%	--
	Kaufman (2008) [31]	2017A-T4	7342 (72)	28 (275)	43 (425)	22%	105/123
	Measurements (this work)	L2 plate	6730 (66)	35 (348)	49 (477)	15%	127/150
	Breguet Standards (1957) [16]	A-U4G1 days (sheet) 500°C Matured for 4	7000 (69)	32 (314)	45 (441)	15%	110/129
	Revue de l'Aluminium (1954) [30]	A-U4G1	--	31-46 (304-451)	44-58 (431-569)	10-21%	--
	Kaufman (2008) [31]	2024-T4	7444 (73)	33 (325)	48 (470)	19-20%	120/141
	* Calculated according to standard ASTM designation E 140: HB = 3.76211 + 0.825368 (HV) for wrought
	aluminum alloys

Acknowledgements

This work was supported by RTRA (Réseau Thématique de Recherche Avancée) -Science and Technologies for Aeronautics and Space, through the VIMA project (Aging of Advanced Materials), by doctoral fellowship PRES (Ref. 5913)/ Région (Ref. 12053330) through the MIPROPAL project and by Programme Investissements d'Avenir under the program ANR-11-IDEX-0002-02, reference ANR-10-LABX-0037-NEXT.

The authors would like to thank the association Ailes Anciennes Toulouse, especially the team "Breguet-2 ponts" for providing the samples and precious information on the aircraft. Dominique Lamirault of the preparation service at CEMES for the thin slices used in TEM is also warmly acknowledged.


Santiago Arroyave-Tobón

Denis Teissandier

Vincent Delos

Applying screw theory for summing sets of constraints in geometric tolerancing

Keywords: Manufacturing defects, Tolerance analysis, Screw theory, Degrees of freedom, Polyhedra, Minkowski sum

In tolerance analysis, approaches based on sets of constraints (also called convex hull techniques) are able to study simultaneously all the possible extreme configurations of a mechanism when simulating manufacturing defects in its components. The accumulation of these defects can be calculated by summing and intersecting 6-dimensional sets of constraints, i.e. polyhedra. These approaches tend to be time-consuming, however, because of the complexity resulting from manipulating sets in R 6 . In this paper, polyhedra are decomposed into a bounded set (a polytope) and an unbounded set (a set of straight lines). The unbounded part of the polyhedra is characterized by the degrees of freedom of the toleranced feature or the joint. Therefore, the decomposition can be performed based on a kinematic analysis of the studied assembly using screw systems. The proposed decomposition is presented for the most common features used in geometric tolerancing. The idea behind this strategy is, instead of summing polyhedra in R 6 , to sum only their underlying polytopes by isolating the unbounded part of the operands. A slider-crank mechanism is used to show the gain in computational time of the proposed method in comparison with the strategy based on complete 6-dimensional sets of constraints.

Introduction

Tolerance analysis consists in studying the impact of manufacturing defects on product behavior. The objective is to determine if the accumulation of defects in the different parts allows the ideal functioning of the product. One way to perform tolerance analysis of three-dimensional tolerance chains is by means of sets of constraints. Using these approaches, it is possible to simulate simultaneously all possible displacements (3 rotations and 3 translations) of a toleranced feature inside its tolerance zone or the relative movements of two features of two parts mated by a contact restriction [START_REF] Fleming | Geometric relationships between toleranced features[END_REF]. They allow, additionally, the treatment of overconstrained (hyperstatic) mechanisms. Some models consider quadratic constraints (such as Domains [START_REF] Giordano | Clearance space in volumic dimensioning[END_REF] and T-Maps [START_REF] Davidson | A new mathematical model for geometric tolerances as applied to round faces[END_REF]), while others are based on linear constraints (such as Polytopes [START_REF] Teissandier | Operations on polytopes: application to tolerance analysis[END_REF]) by considering manufacturing defects as small displacements [START_REF] Bourdet | The concept of the small displacement torsor in metrology[END_REF]. A comparison of these methods is presented in [START_REF] Shah | Navigating the tolerance analysis maze[END_REF]. As all possible configurations of each toleranced feature must be considered, the cumulative defect limits on parts in a serial configuration can be calculated by means of the Minkowski sum of sets of geometric constraints [START_REF] Fleming | Geometric relationships between toleranced features[END_REF][START_REF] Srinivasan | The role of sweeps in tolerancing semantics[END_REF]. In the case of parts with parallel contacts, all the tolerance zones must be satisfied and then the intersection of the respective sets of constraints is required. By means of these two operations, the stack-up of deviations in any tolerance chain can be computed [START_REF] Arroyave-Tobón | Tolerance analysis with polytopes in HV-description[END_REF].

The degrees of freedom of contact features or the degrees of invariance of geometric features imply unbounded displacements [START_REF] Clément | A Declarative Information Model for Functional Requirements[END_REF][START_REF] Giordano | Clearance space and deviation space[END_REF]. In these cases, the respective set of constraints defines a polyhedron in R 6 . When the toleranced feature has no degrees of freedom or invariance (as in the case of a complex surface), the set of constraints represents a polytope (a bounded polyhedron).

On the one hand, Minkowski sums of polytopes have been a research subject explored with applications in biology [START_REF] Bernholt | Constrained minkowski sums: A geometric framework for solving interval problems incomputational biology efficiently[END_REF], robotics [START_REF] Fogel | Exact and efficient construction of minkowski sums of convex polyhedra with applications[END_REF][START_REF] Rakovic | Computation of invariant sets for piecewise affine discrete time systems subject to bounded disturbances[END_REF], tolerance analysis [START_REF] Wu | Improvements to algorithms for computing the minkowski sum of 3-polytopes[END_REF][START_REF] Teissandier | Algorithm to calculate the Minkowski sums of 3-polytopes based on normal fans[END_REF], among others. Because of their complexity, however, the computation of sums in spaces of dimension greater than 3 is still a challenging matter. Mansuy et al. [START_REF] Mansuy | A new calculation method for the worst case tolerance analysis and synthesis in stack-type assemblies[END_REF] propose a method for calculating separately only the most disadvantageous vertices with respect to a functional polytope. Even if this method avoids the computation of Minkowski sums, the set of computed vertices is only representative of a given functional condition. In addition, the authors just consider the case of tolerance chains made up of features of the same invariance class and in a particular relative position (i.e. a set of parallel planes). In these cases, the displacements linked to degrees of freedom or invariance can be trivially isolated from the computations.

On the other hand, for the Minkowski sums of polyhedra, just a few studies have been proposed. Homri et al. [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF] suggested turning polyhedra into polytopes by introducing virtual boundaries, called cap half-spaces. This strategy has to cope with the multiplication of cap half-spaces during the computation of Minkowski sums. As a consequence, the time taken for computing cap facets (facets associated with cap half-spaces) is in general far greater than that needed for computing facets representing real limits of bounded displacements. Fukuda [START_REF] Fukuda | From the zonotope construction to the minkowski addition of convex polytopes[END_REF] presented an algorithm to compute Minkowski sums of polytopes, mentioning the possibility of applying the same procedure for the case of polyhedra with at least one vertex (pointed polyhedra) by treating infinite rays as points at 'infinity'. However, due to the degrees of freedom (or invariance), the polyhedra manipulated in tolerancing usually do not have vertices. Each degree of freedom (or invariance) implies a sweeping operation of a polytope along a straight line, placing the vertices at infinity.

The set of displacements required to completely define the relative position between two features in a tolerance chain can be identified by a kinematic analysis. At a later stage, the objective of tolerance analysis is to determine the limits of these displacements. For this reason it is common to find tolerance analysis techniques in the literature involving kinematic analysis [START_REF] Frisoli | A new screw theory method for the estimation of position accuracy in spatial parallel manipulators with revolute joint clearances[END_REF][START_REF] Kumaraswamy | A unified framework for tolerance analysis of planar and spatial mechanisms using screw theory[END_REF][START_REF] Desrochers | Application of a unified jacobian-torsor model for tolerance analysis[END_REF][START_REF] Laperrière | Tolerance analysis and synthesis using jacobian transforms[END_REF]. However, it has been done just for the case of parametric techniques (which are not able to treat assemblies with redundantly suppressed degrees of freedom). The advantages of applying kinematic analysis in the case of techniques based on sets of constraints remain unexplored, and this is actually the objective of this paper.

One possible way to model rigid body motions is by means of the theory of Euclidean group displacements. This theory states that in the affine Euclidean space of dimension 3, the set of rigid transformations has an algebraic structure of a continuous group, more specifically it is a Lie group of dimension 6 [START_REF] Marie | Intrinsic formulation of problems of geometry and kinematics of mechanisms[END_REF]. Therefore, the union of two successive displacements is also a displacement [START_REF] Hervé | The Lie group of rigid body displacements, a fundamental tool for mechanism design[END_REF]. Fanghella [START_REF] Fanghella | Kinematics of spatial linkages by group algebra: A structure-based approach[END_REF] presented an exhaustive list of products of subgroups. The author produced this classification by analyzing systematically all possible couples of subgroups and their geometric relations (coincidence, parallelism, intersection, etc.). Hervé [START_REF] Hervé | The mathematical group structure of the set of displacements[END_REF] presented a more general procedure by using the exponential maps of the Lie groups enumerating 12 possible kinds of displacement subgroups. However, the 12 displacement subgroups cannot explain all rigid motions in space. Under most conditions, a rigid motion is simply a displacement manifold included in the displacement group but not a subgroup itself [START_REF] Li | Displacement manifold method for type synthesis of lower-mobility parallel mechanisms[END_REF].

In a similar way, using the theory of screws [START_REF] Ball | A Treatise on the Theory of Screws[END_REF] it is possible to perform spatial kinematic analysis. By computing systematically the union and intersection of screws the mobility condition of any two surfaces of a mechanical system can be determined [START_REF] Adams | Application of screw theory to constraint analysis of assemblies of rigid parts[END_REF][START_REF] Gerbino | How to investigate constraints and motions in assemblies by screw theory[END_REF][START_REF] Rico | Screw theory and higher order kinematic analysis of open serial and closed chains[END_REF]. Open serial and closed chains can be analyzed even if the former present redundantly suppressed degrees of freedom. For these purposes, the classification of screw systems plays an important role. Rico Martinez and Duffy [START_REF] José | Classification of screw systemsi. one-and two-systems[END_REF][START_REF] José | Classification of screw systemsii. three-systems[END_REF] proposed a classification based on the reciprocal basis of screw systems. A general approach for mechanisms analysis and synthesis based on this reciprocity is proposed by Dai and Jones [START_REF] Dai | Interrelationship between screw systems and corresponding reciprocal systems and applications[END_REF].

This paper proposes the use of the theory of screws for modeling the mobility conditions of an assembly during tolerance simulation with sets of geometric constraints. It enables us to simplify the sets of constraints, decomposing them into the sum of a polytope (a bounded set) and a set of straight lines (an unbounded set). When the operand sets are decomposed, their sum can be calculated by dealing with the underlying polytopes. This implies a reduction in operand complexity and consequently a reduction in the computational time. In this work, the advantages of using the proposed model (decomposed polyhedra) are explored only in the case of sums, i.e. when the tolerance chain is made up of serial contacts.

The next section of the paper introduces the way to sum sets of constraints based on polytopes and cap half-spaces. In section 3 the new strategy based on kinematic analysis is presented. Next, in section 4, a case study is solved using the proposed method and the results are compared with those obtained by the method based on cap half-spaces. Finally, some conclusions and perspectives are discussed in section 5.

Tolerance analysis with polytopes of R 6

In this section, the tolerance analysis method based on the operation with polytopes of R 6 is presented with its mathematical support and illustrated with an example. For further details the reader can refer to [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF].

Before going into the details of polyhedra and polytopes in the field of geometric tolerancing, let's present some definitions required for this paper. They are taken from [START_REF] Günter | Lectures on polytopes[END_REF]. Definition 2.1. An H-polyhedron is an intersection of many finitely closed half-spaces in some R n . Definition 2.2. A polytope is a polyhedron that is bounded in the sense that it does not contain a ray {x + ty : t ≥ 0, x, y ∈ R n } for any y 0. A V-polytope is the convex hull of a finite set of points in some R n .

According to the Minkowski-Weyl theorem, the H-description (where H refers to half-space) and the V-description (where V refers to vertex) of a polytope are equivalent. It makes the link between the two dual intuitive definitions of a polytope: one based on facets and the other based on vertices (see Figure 1).

Sets of constraints and polyhedra

Let S 1 be a substitute surface related to a nominal feature S 0 (see Figure 2). The position and orientation of S 0 can be restricted by a specification, which in turn defines a tolerance zone T Z, offsetting S 0 from d in f to d sup . This tolerance zone implies a restriction on the position of the points N i ∈ S 0 :

S 1 ⊆ ZT ⇔ ∀N i ∈ S 0 : d in f ≤ t N i • n i ≤ d sup (1)
where t N i is the translation displacement of S 1 in relation to S 0 at point N i .

Figure 2: Toleranced feature and tolerance zone [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF] Expressing the constraints in (1) at any point M of the Euclidean space and linearizing the rotations under the assumption of small displacements we have:

d in f ≤ (t M + N i M × r) • n i ≤ d sup ( 2
)
where r is the rotation vector of S 1 in relation to S 0 . Considering a vectorial base [u, v, w], we define the vectors r

= [r u , r v , r w ] T , t M = [t u , t v , t w ] T , N i M = [d iu , d iv , d iw ] T , n i = [n iu , n iv , n iw ] T .
When expanding the vectorial and scalar products in inequality (2) and setting r u = x 1 , r v = x 2 , r w = x 3 , t u = x 4 , t v = x 5 and t w = x 6 , we obtain:

d in f ≤(n iv d iw -n iw d iv )x 1 + (n iw d iu -n iu d iw )x 2 + (n iu d iv -n iv d iu )x 3 + n iu x 4 + n iv x 5 + n iw x 6 ≤ d sup (3)
These constraints can be modeled with half-spaces. For each point N i ∈ S 0 , two parallel half-spaces are obtained. Then, if a finite set of m points N i is considered, a set of 2m constraints is obtained [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF]. The intersection of these constraints defines a convex H-polyhedron in R 6 [START_REF] Teissandier | Operations on polytopes: application to tolerance analysis[END_REF]:

Γ = 2m k=1 H- k where H- k = x ∈ R 6 : b k + a k1 x 1 + ... + a k6 x 6 ≥ 0 for all k ∈ {1, ..., 2m} (4)
In a similar way, a polyhedron representing the allowable displacements of a couple of features potentially in contact can be defined. In this case, the toleranced feature is defined in the nominal case of permanent contact between the features and the tolerance zone is obtained by offsetting the toleranced feature according to the clearance value.

Finally, a functional polytope can be created with the same procedure that defines a geometric specification between two surfaces of the same part.

Except for the case of complex surfaces, the set of constraints derived from a toleranced feature (Eq. ( 4)) defines an unbounded set in the 6-dimensional space of deviations, i.e. a polyhedron. This is a consequence of the degrees of invariance of the toleranced features or the degrees of freedom of the joints which define unbounded displacements. These polyhedra are generally not pointed (according to definition 2.3). As a consequence, they do not have vertices and thus, the computation of double descriptions and Minkowski sums becomes more complicated.

Polytopes and cap half-spaces

Due to the aforementioned difficulties treating polyhedra, Homri et al. [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF] suggested turning them into polytopes introducing virtual boundaries called cap facets, Hc:

Γ =         2m k=1 H- k         ∩         2d j=1 Hc - j         = Γ ∩         2d j=1 Hc - j         (5)
where d is the number of degrees of invariance (or freedom) of the toleranced (or kinematic) feature. Notice that the strict minimum of cap facets is introduced to each polyhedron. An example of the introduction of cap facets is presented in Figure 3. In this figure a 6D polyhedron derived from the geometric constraints of a plane surface is depicted in two partial 3D representations. The bottom left polyhedron is a partial visualization of the rotations and the right is one of the translations. As the tolerance zone imposes limits on r x , r y and t z , cap facets are required to virtually limit the unbounded displacements Figure 3: Cap half-spaces to virtually limit r z , t x and t y depicted in red, blue and green respectively. r z (in the polyhedron at the right), t x and t y (in the polyhedron at the left). The decision to add cap facets to obtain bounded sets is only for to algorithmic reasons.

The only issue arising from the method based on cap facets is that the topology of the operand polytopes is affected, increasing its complexity. Although such an issue does not influence the results from the tolerancing point of view, it does involve spending time calculating meaningless information. During a simulation, this problem increases after each sum due to the accumulation of the degrees of freedom along the toleranced chains. Therefore, after the first sum, the definition of the set of cap facets of the operand polytopes is no longer minimal. As a consequence, the time for computing cap facets is in general far greater than that for computing significant facets.

Example of tolerance analysis with polytopes of R 6

A tolerance analysis of a sub-assembly of a slider-crank mechanism is presented below in order to illustrate the method based on polytopes and cap half-spaces. Later, in section 4, the whole assembly is addressed comparing the actual approach with the new enhanced method. As depicted in Figure 4, the functional condition FC1 implies the control of the relative position of the inferior cylindrical surface of the connecting rod and the upper plane surface of the piston. These surfaces are designated 1,2 and 3,2 in Figure 5. Contact and geometrical defects in the involved parts are modelled with polytopes for calculating the stack-up of deviations.

The topological model of the mechanism is depicted by the contact graph in Figure 5b. The nodes designated as α, β represent the nominal model of a part when β = 0, and the substitute surfaces when β 0. Each edge of the graph represents some deviations. These may be geometric deviations, in the case of inner edges, or deviations due to contacts, in the case of edges connecting two nodes from different parts. These deviations can be represented by geometric and contact polyhedra respectively [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF].

By analyzing the contact graph, the set of operations required to simulate the relative position of the surfaces involved by the functional condition can be determined. The following relation can then be deduced:

Γ 1,2/3,2 = Γ 1,2/1,0 ⊕ Γ 1,0/1,1 ⊕ Γ 1,1/2,1 ⊕ Γ 2,1/3,1 ⊕ Γ 3,1/3,0 ⊕ Γ 3,0/3,2 (6)
where Γ f /g is the polyhedron describing the relative position of the surface f with respect to the surface g. In Eq. ( 6) some polyhedra are defined over the same surface and consequently they are homothetic. This is the case for polyhedra Γ 1,0/1,1 and Γ 1,1/2,1 and polytopes Γ 2,1/3,1 and Γ 3,1/3,0 . The sum of homothetic polyhedra can be performed directly by homothetic transformations and no numerical computation is needed. Then, Eq. ( 6) can be simplified as: In order to calculate Γ 1,2/3,2 , the strategy proposed by Homri et al. [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF] suggests the introduction of cap facets to the initial operand polyhedra in order to bound them and obtain polytopes. Polytopes Γ 1,0/2,1 , Γ 2,1/3,0 and Γ 3,0/3,2 (from CAD cylindrical surfaces oriented along the x-axis) required two pairs of cap facets to limit rotation and translation along the xaxis (hereafter r x and t x respectively). Three pairs of cap facets were introduced to Γ 1,0/1,2 (defined from a plane surface with normal along y-axis) to impose virtual limits along t x , t z and r y . The truncation of the operands was required to obtain their list of vertices (see [START_REF] Arroyave-Tobón | Tolerance analysis with polytopes in HV-description[END_REF] for further details). After obtaining all the operands in V-description and belonging to the same space, the Minkowski sums of Eq. ( 7) were computed. The computations were performed with the open source library politopix: http://i2m.u-bordeaux.fr/politopix.html.

Γ 1,2/3,2 = Γ 1,2/1,0 ⊕ Γ 1,0/2,1 ⊕ Γ 2,1/3,0 ⊕ Γ 3,0/3,2 (7)
Figure 7 shows a 3D representation of the first sum of the simulation. It is worth mentioning that some of the figures showing 3D polytopes in this paper are just partial representations since their original polytopes belong to spaces with a dimension larger than three. A projection into a three-dimensional subspace was therefore carried out. For this case, the operands were calculated at the centroid of surface 2,1 and were then projected to the subspace spanned by [r x , r z , t y ] for visualization purposes. Notice that from a strict minimum set of cap facets (facets in red) in the operands, originally required just to close the sets, many cap facets and consequently many unnecessary vertices appear in the calculated polytope Γ 1,2/2,1 . The result of this sum, which models the possible deviations of a plane surface oriented along the y-axis (face 1,2) with respect to a cylindrical surface oriented along the x-axis (face 2,1), shows that only two displacements can be controlled: r z and t y (this conclusion is justified further on with a kinematic analysis). With this strategy, however, many vertices are computed along the directions of unbounded displacements: r x , r y , t x and t z in this case. This situation is clearly illustrated in Figure 7, where many unnecessary facets where calculated along r x .

When the result of a sum is used again as an operand for a subsequent sum, this issue becomes worse since the operands are already highly 'contaminated' with vertices from cap facets. Notice that at the lats step of the simulation, for which the results are presented in Table 1, the calculated polytope is made up of 99.9% of cap facets. This can be justified by the accumulation and propagation of the degrees of freedom along the tolerance chain. The cap facets were identified following the strategy proposed in [START_REF] Vincent Delos | How to trace the significant information in tolerance analysis with polytopes[END_REF]. The above results show clearly a situation to be considered to reduce the computational times. The main issue is how to avoid computing worthless information (cap facets) in tolerance analysis with sets of constraints. One possible solution is presented in the next section, where a new method is described for summing polyhedra based on the decomposition of the operand polytopes.

On the sum of decomposed polyhedra

According to the Minkowski-Weyl theorem [START_REF] Günter | Lectures on polytopes[END_REF], a polyhedron can be decomposed into the sum of a bounded part and an unbounded part, namely a polytope P and a polyhedral cone C:

Γ = P ⊕ C (8)
The polyhedral cone of a polyhedron can be identified considering the degrees of invariance of the toleranced feature or the degrees of freedom of the joint. Each of these unbounded displacements characterizes a straight line ∆ j , in R 6 :

C = d j=1 ∆ j with d ≤ 6 (9)
It is along these straight lines that the underlying polytope (derived from geometric or contact constraints) is 'extruded' to generate a polyhedron (see Figure 8).

Identification of unbounded displacements by kinematic analysis

Because of their linear properties, screw systems are a very suitable tool for analyzing the instantaneous kinematics of assemblies for tolerance analysis purposes. According to this theory, an allowable motion can be represented as a twist and a restricted motion as a wrench [START_REF] Su | Type synthesis of freedom and constraint elements for design of flexure mechanisms[END_REF]. Both twists T and wrenches Ŵ are 1x6 vectors, written as [START_REF] Jeffrey | Feature based analysis of selective limited motion in assemblies[END_REF]:

T = [ω | v] = [ω | r × ω] (10) Ŵ = [ f | m] = [ f | r × f ] ( 11
)
where ω is a unit angular velocity vector, v is a unit linear velocity vector, f is a unit force vector, m is unit moment vector about the point r. The advantage of using screws is that the mobility analysis can be performed just by using the geometric parameters of the features. The degrees of freedom of a joint characterize a subspace of R 6 : the twist-space. This subspace is defined by the set of d twists describing the mobility condition of the feature [START_REF] Konkar | Incremental Kinematic Analysis of Mechanisms[END_REF]:

T =                T1 T2 ... Td                (12)
T is usually called twist-matrix and its dual vector space, represented by a wrench-matrix, corresponds to the subspace of bounded displacements: the wrench-space. As the twist-space and the wrench-space are orthogonal subspaces, it is possible to compute a wrench-matrix from a twist-matrix and vice versa. This property, called reciprocity of screws, is one of the most important properties of this theory since changing from one subspace to the other is made easy. This duality is valid if and only if the virtual work of the wrenches on the twists is zero. It can be said if it is assumed that the constrained parts are rigid, the contact surfaces are frictionless and the contact between parts does not break [START_REF] Gerbino | How to investigate constraints and motions in assemblies by screw theory[END_REF].

The mobility conditions of a set of features can describe as well as it is done for a single one. If the features are assembled in a serial configuration, the twist-matrices must be joined together, as in Eq. ( 12). When the features are assembled in a parallel configuration (as in the case of parts mated with multiple contacts), the set of unbounded displacements can be calculated as the intersection of the twists, or more simply, as the union of the corresponding wrenches [START_REF] Konkar | Incremental Kinematic Analysis of Mechanisms[END_REF][START_REF] Jeffrey | Feature based analysis of selective limited motion in assemblies[END_REF]. Then, by means of unions and intersections of screws, the mobility conditions of any couple of faces in a tolerance chain can be characterized.

A screw can be interpreted as a vector in the 6-dimensional space of Euclidean motions [START_REF] Huang | The correspondence between finite screw systems and projective spaces[END_REF]. For a given feature, a twist can be seen as the director vector of a straight line ∆ j in R 6 characterizing a non-restricted displacement, as depicted in Figure 8. By duality, this vector represents a hyperplane. In turn, the intersection of the hyperplanes orthogonal to each of the twists of a feature characterizes the wrench-space. As the geometric constraints obtained in tolerancing are limits in the directions of bounded displacements, they define a polytope in the wrench-space. The sweeping of this polytope along the straight lines of the non-restricted displacement defines an operand polyhedron in R 6 .

Polyhedra decomposition

In accordance with the above, the mobility conditions of a kinematic joint can be identified according to its degrees of freedom. We distinguish three types of kinematic joints according to their type contact feature: surfaces, lines and point. The contact feature is defined in the nominal case of permanent contact between the surfaces of the kinematic joint [START_REF] Homri | Tolerance analysis by polytopes: Taking into account degrees of freedom with cap half-spaces[END_REF].

According to the ISO standards [START_REF]Kinematic diagrams -graphical symbols, part 1[END_REF], we consider in the case of surfacic contact features: spherical, planar, cylindrical, revolute and prismatic pairs. The straight lines characterizing the unbounded displacements for these cases are presented in Table 2. The helical pair is not treated as it is not of interest for geometric tolerancing. In the case of linear contact features, we consider ball-and-cylinder and cylinder-and-plane pairs. Their decomposition is given in Table 3. In the case of punctual contact feature, we consider the ball-and-plane pair (see Table 4). Even if the same reasoning can be applied for the exhaustive list of cases, we only address those of interest in the subject of tolerancing.

As a toleranced feature has no mobility, it is possible to see it as a joint (usually called internal joint) mating the nominal feature and the substitute one. The mobility in this case is characterized by the set of displacements leaving the toleranced feature invariant in the tolerance zone. In other words, the degree of invariance of an internal joint (i.e. a toleranced feature) can be modelled the same way we model the degree of freedom of a kinematic joint. Tolerance chains generally involve internal and kinematic joints alternatively. The general procedure for identifing the unbounded displacements in a tolerance chain is summarized in Figure 9.

The invariance degree of a surface is defined according to 7 invariance classes [START_REF]Geometrical tolerancing datums and datum systems[END_REF][START_REF] Desrochers | A dimensioning and tolerancing assistance model for cad/cam systems[END_REF]. We consider spherical, plane, cylindrical, of revolution, prismatic and complex classes (see Table 2). Like the helical pair, the helical class is not considered because it is not of interest in geometric tolerancing.

Once the straight lines of the polyhedron are identified, the underlying polytope can be obtained. This can be done by intersecting the original set of constraints (which is obtained by applying double inequality (3) to a set of points belonging to a toleranced feature) with a set of hyperplanes H ∆ j such that H ∆ j ⊥ ∆ j ∀ j ∈ {1, ..., d}. Tables 2, 3 and 4 present this decomposition for sets of constraints derived from surfaces, lines and a point, respectively.

Summing the underlying polytopes of two polyhedra

The decomposition of operand polyhedra proposed in the previous section gives way to simplified techniques for calculating the accumulation of manufacturing defects in mechanical systems. Next, instead of operating polyhedra in R 6 or

P = 2m k=1 H- k H ∆ 1 H ∆ 2 H ∆ 3 ∆ 1 = x = λ Tru = λ[u | N i M × u] ∆ 2 = x = λ Trv = λ[v | N i M × v] ∆ 3 = x = λ Trw = λ[w | N i M × w]
Planar pair Plane surface

P = 2m k=1 H- k H ∆ 3 H ∆ 4 H ∆ 5 ∆ 3 = x = λ Trw = λ[w | 0] ∆ 4 = x = λ Ttu = λ[0 | u] ∆ 5 = x = λ Ttv = λ[0 | v]
Cylindrical pair Cylindrical surface

P = 2m k=1 H- k H ∆ 1 H ∆ 4 ∆ 1 = x = λ Tru = λ[u | N i M × u] ∆ 4 = x = λ Ttu = λ[0 | u]
Revolute pair Surface of revolution

P = 2m k=1 H- k H ∆ 1 ∆ 1 = x = λ Tru = λ[u | N i M × u]
Prismatic pair Prismatic surface its corresponding polytopes made up of cap facets, we propose the manipulation of simplified sets of constraints after the exclusion of the straight lines derived from the unbounded displacements (according to Tables 2, 3 and4).

P = 2m k=1 H- k H ∆ 4 ∆ 4 = x = λ Ttu = λ[0 | u] N/A Complex surface P = 2m k=1 H- k N/A
Tolerancing constraints Lines of unbounded displacements Ball-and-cylinder pair

P = 2m k=1 H- k H ∆ 1 H ∆ 2 H ∆ 3 H ∆ 4 ∆ 1 = x = λ Tru = λ[u | N i M × u] ∆ 2 = x = λ Trv = λ[v | N i M × v] ∆ 3 = x = λ Trw = λ[w | N i M × w] ∆ 4 = x = λ Ttu = λ[0 | u],
Cylinder-and-plane pair

P = 2m k=1 H- k H ∆ 1 H ∆ 3 H ∆ 4 H ∆ 5
with m = 2

∆ 1 = x = λ Tru = λ[u, N i M × u] ∆ 3 = x = λ Trw = λ[w | 0] ∆ 4 = x = λ Ttu = λ[0 | u] ∆ 5 = x = λ Ttv = λ[0 | v]

Kinematic joint Tolerancing constraints

Lines of unbounded displacements Ball-and-plane pair

P = 2m k=1 H- k H ∆ 1 H ∆ 2 H ∆ 3 H ∆ 4 H ∆ 5 with m = 1 ∆ 1 = x = λ Tru = λ[u | 0] ∆ 2 = x = λ Trv = λ[v | 0] ∆ 3 = x = λ Trw = λ[w | 0] ∆ 4 = x = λ Ttu = λ[0 | u] ∆ 5 = x = λ Ttv = λ[0 | v]
Table 4: Decomposition of polyhedra derived from punctual contacts, x = {x 1 , ...,

x 6 } ∈ R 6 , λ ∈ R.
Let Γ 1 and Γ 2 be two polyhedra, then their sum can be calculated by summing their underlying polytopes in the smallest possible subspace. Therefore, according to Eq. ( 8) the following can be posed:

Γ 1 ⊕ Γ 2 = P 1 ⊕ P 2 ⊕ C 1 ⊕ C 2 ( 13
)
The idea is then to isolate the two polyhedral cones C 1 and C 2 from the sum to reduce the complexity. The sum P 1 + P 2 is, in general, a sum of two non-full-dimensional polytopes of R n (with 1 ≤ n ≤ 6). From Eq. ( 13) the synthesis of an algorithm to compute the sum of decomposed polyhedra is straightforward:

Algorithm 1 Sum of decomposed polyhedra Require: Γ A , Γ B Ensure: Γ S = Γ A ⊕ Γ B 1: decompose Γ A = P A ⊕ C A 2: decompose Γ B = P B ⊕ C B 3: if dim(C A ⊕ C B ) = 6 then 4: Γ S = R 6 5: else 6: compute V-description V A = conv(a 1 , ..., a r ) 7: compute V-description V B = conv(b 1 , ..., b s ) 8: compute sum V S = V A ⊕ V B 9: Γ S = V S ⊕ C A ⊕ C B 10: end if
Although the "natural way" to perform the sum of two polytopes is to compute the sum of all the vertices of the operands and calculate after its convex hull, this is a very time-consuming strategy. In addition, the points which are not vertices but are located on the boundary of the calculated polytope pose an additional challenge. There are some improved methods which can treat non-dimensional polytopes, such as those presented in [START_REF] Delos | Minkowski sum of polytopes defined by their vertices[END_REF] and in [START_REF] Fukuda | From the zonotope construction to the minkowski addition of convex polytopes[END_REF]. The first algorithm takes advantage of the property of the uniqueness of the Minkowski vertices decomposition for computing just the vertices defining the convex hull of the resulting polytope. The second one exploits the adjacency properties of the operands' vertices for computing just the Minkowski vertices of the resulting polytope.

4. Case study: slider-crank mechanism

FC1 simulation

The same case study already presented in section 2.3 was solved following the strategy based on the sum of decomposed polyhedra. The centroid of surface 2,1 was chosen as the point to express the screws for the mobility analysis. For convenience, it is the same point used to express the polytopes in Section 2.

Applying the decomposition proposed in Tables 2, 3 and 4 to the set of constraints derived from plane surface 1,2, the following can be obtained:

Γ 1,2/1,0 = P 1,2/1,0 ⊕ C 1,2/1,0
where the polyhedral cone C 1,2/1,0 is defined by the space spanned by the rows of the matrix:

T 1,2 =           r y t x t z           =           0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1          
When isolating the 3-dimensional polyhedral cone C 1,2/1,0 from the operand Γ 1,2/1,0 , the 3-dimensional polytope P 1,2/1,0 was obtained. The truncation of P 1,2/1,0 was performed in a 3D subspace in order to obtain its V-description. Figure 10 shows a 3D representation of this operand. Notice in Table 5 that the number of vertices of P 1,2/1,0 decreased 2 3 times with respect to the equivalent operand Γ 1,2/1,0 from the strategy based on polytopes of R 6 (see Table 1). This is due to the sweeping operation that must be performed recursively 3 times to obtain the polytope Γ 1,2/1,0 in R 6 from the polytope P 1,2/1,0 .

Similarly, the polyhedral cone of the operand related to surface 2,1 is characterized by the twist-matrix T 2,1 presented below.

T 2,1 = r x t x = 1 0 0 0 0 0 0 0 0 1 0 0

In Figure 10, a 3D representation of the decomposition of this operand is given. As r x is an unbounded displacement for the cylindrical surface 2,1, this polyhedron can be decomposed putting aside the straight line ∆ r x . Finally, the list of vertices of P 1,0/2,1 was calculated in a 4D subspace. Once we had the V-description of the operands, the Minkowski sum P 1,2/2,1 = P 1,2/1,0 ⊕ P 1,0/2,1 could be performed. A summary of the results is presented in Table 5. As the manipulated operand polytopes were less complex, computing the sum took less time.

For the next sum, the set of straight lines (or the polyhedral cone) of the operand Γ 1,2/2,1 was calculated joining the twistmatrices T 1,2 and T 2,1 . Then, P 1,2/2,1 is the underlying polytope of the polyhedra Γ 1,2/2,1 , of which the polyhedral cone is characterized by T 1,2/2,1 .

T 1,2/2,1 =               r x r y t x t z               =              
1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1

             
The other two Minkowski sums in Eq. ( 7) were performed following the same strategy. A summary of all the simulation is presented in Table 5.

The main contribution of the proposed method can be observed in a reduction of the computation time by 97%. Furthermore, the manipulation of less complex polytopes reduces the probabilities of having numerical problems during the calculations.

It is worth mentioning that the simplification introduced by using the proposed approach keeps the set of constraints representing real limits in the associated tolerancing problem unaltered. The new method merely avoids calculating worthless information (constraints along axes of unbounded displacements). In fact, the equivalence of the results obtained by this method and the results obtained by the method based on polytopes of R 6 was verified. This was done by projecting the results of both methods in the subspace of bounded displacements (the subspace spanned by [r z , t y ]) and checking the equality of the two polytopes, as shown in Figure 11.

FC2 simulation: whole crank slider mechanism

Tolerances along the whole crank slider mechanism were simulated to analyze FC2. This requirement requires control geometrical variations in the combustion chamber between the upper planar surface of the block and the upper planar surface of the piston, surfaces 5,2 and 1,5 respectively in figure 5a. The surfaces belonging to the block, surfaces 5,1, 5,2 and 5,3, are assumed rigidly fixed. As simulation had to be performed for a given position of the mechanism, we have chosen the extreme configuration with the piston in top dead centre.

The fact of having a closed loop implies the computation of an intersection, as it can be noticed in Eq. 14. Assuming some homothetic transformations, the operations required for the simulation of FC2 are:

Γ 5,2/1,2 = Γ 5,2/5,0 ⊕ Γ 5,0/1,0-a ∩ Γ 5,0/1,0-b (14)
where

Γ 5,0/1,0-b = Γ 1,0/2,1 ⊕ Γ 2,1/3,0 ⊕ Γ 3,0/4,0 ⊕ Γ 4,0/5,0 (15)
The computation of Γ 5,0/1,0-b can be performed in a 4D space isolating the straight lines associated to the DOFs of the operands, which in turn are related with the twists: The next operation, the intersection, demands the union of the wrenches associated to the operands for modelling the mobility conditions of the surface 5,0 with respect to 1,0:

T 5,0/1,0-b =                                                    T 1,0/2,1 T 2,1/3,0 T 3,0/4,0 T 4,0/5,0                                                    =                                                    1 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 d 1z -d 1y 1 0 0 0 0 0 0 0 0 0 d 2z -d 2y                                                    (
W 5,0/1,0 =                                    W 5,0/1,0-a W 5,0/1,0-b                                    =                                   
1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0

                                   (17)
As it can be noticed, W 5,0/1,0 can be reduced to the identity matrix, which means that Γ 5,0/1,0 is a 6D operand and can not be decomposed. Therefore, the final twist matrix is composed only of the mobility introduced by Γ 5,2/5,0 :

T 5,2/1,2 =           0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1           (18)
The summary of the simulation for FC2 is presented in table 6. With the proposed strategy the computational time for computing the sums was reduced by 76% with respect to the approach based on 6D operands. The results of FC1 could not be used directly for the simulation of FC2, because the first one is not exactly a subset of the second one. Comparing the time of both simulations, it can be noticed than that for FC1 is greater. This can be justified by the fact that, in this simulation, the operand of the planar surface 1,2 is summed directly with the ones derived from the cylindrical surfaces. For FC2, on the contrary, the first sum is made up of operands derived from cylindrical surfaces whose axes are parallel between them and belong to the same plane. This implies operands polytopes with similar topology, and then, their sum take much less time. Additionally, the intersection performed before the final sum, implies a significant reduction in complexity due to the restriction of all the DOFs in the operand Γ 5,0/1,0 (as shown in Eq. 17).

Conclusions and future work

In this paper, a new method for calculating the cumulative stack up of deviations with sets of constraints was presented. The method is based on the decomposition of geometric and contact polyhedra into the sum of polytopes (the bounded part of the polyhedron) and polyhedral cones (the unbounded part of the polyhedron). When performing a sum of polyhedra, it is possible to exploit this property and take just the bounded parts to calculate the sum. This new way of operating polyhedra was detailed in an algorithm.

The proposed decomposition of polyhedra is based on kinematic analysis of the toleranced features. Using the theory of screws it is possible to model the mobility conditions of a feature with respect to a reference or to another feature. In this way, the set of unbounded displacements can be identified, which in turn characterizes the polyhedral cone of the polyhedron. As the main interest of tolerance analysis lies in calculating the limits of the controllable displacement, which are characterized by the underlying polytopes of the operands, the polyhedral cones of the operands can be isolated from the simulation to reduce the complexity of the manipulated models.

The procedure and the results from applying the method were shown in a case study. A slider-crank mechanism was analyzed. The direct consequence of applying the new method was a 76% reduction in the computational time for the sums. In addition to the reduction in computational time, the method based on the decomposition of polyhedra reduces the complexity of the operands, and hence, the probabilities of having numerical problems during the calculations.

Due to the computational cost of summing of sets of constraints, the advantages of manipulating decomposed polyhedra were explored only for the case of serial assembly. Therefore, further work is required to formalize a strategy when the tolerance chain also contains parallel contacts.

When dealing with unilateral contacts, the polyhedra obtained from the sets of constraints are pointed (according to definition 2.3). As a consequence, the operand polyhedra cannot be completely decomposed as proposed in this paper. In these cases, the strategy based on cap half-spaces can be useful for simulating bounded sets in order to perform the computations.

In product design, the fact of recognizing the most influential features in a tolerance chain is essential. This information enables tolerances to be allocated in an optimal way in order to reduce production costs. In tolerance analysis, sets of constraints ensure that this can be achieved thanks to the traceability of the vertices and facets through the whole computational process. In this respect, the method presented in this paper is a first filter to identify among all the facets of the operands those that really influence the FC. However, further work is still required to develop a strategy with which to qualify the level of influence of the facets of the operands on the facets of the functional polytope (the polytope representing the FC).

Figure 1 :

 1 Figure 1: Two equivalent definitions of a polytope: facets or half-spaces versus vertices

Figure 4 :

 4 Figure 4: Proposed case study: slider-crank mechanism

Figure 5 :Figure 6 :

 56 Figure 5: CAD and topological model of the assembly (CP: cylindrical pair)

Figure 7 :

 7 Figure 7: 3D representation of the first sum in R 6 for FC1. Cap facets are shown in red.

Figure 8 :

 8 Figure 8: Decomposition of a polyhedron into the sum of a polytope and a straight line

Figure 9 :

 9 Figure 9: Flow chart of the decomposition of polyhedra

Figure 10 :

 10 Figure 10: First sum of the computations for FC1 with decomposed polyhedra

Table 5 :Figure 11 :

 511 Figure 11: Comparison of the results obtained by both methods

16 ) where d 1 =

 161 [d 1x d 1y d 1z ] and d 2 = [d 2x d 2y d 2z ] are respectively the position vectors of the local reference system of the surfaces 3,2 and 4,2 regarding the point chosen for the mobility analysis (centroid of surface 2,1).

Figure 12 :

 12 Figure 12: 3D representation of simulation of FC2 in 6D

Table 1 :

 1 Summary of the computations in R 6 for FC1 performed with an Intel Core i7-3740QM.

	Polytope	Total number of facets Percentage of cap facets Vertices Time [s]
	Γ 1,2/1,0	50	12%	192	-
	Γ 1,0/2,1	20	20%	256	-
	Γ 2,1/3,0	20	20%	256	-
	Γ 3,0/3,2	20	20%	256	-
	Γ 1,2/2,1 = Γ 1,2/1,0 ⊕ Γ 1,0/2,1	1144	99.6%	4032	25
	Γ 1,2/3,0 = Γ 1,2/2,1 ⊕ Γ 2,1/3,0	5208	99.8%	16320	3730
	Γ 1,2/3,2 = Γ 1,2/3,0 ⊕ Γ 3,0/3,2	19810	99.9%	55808	58870

Table 2 :

 2

Decomposition of polyhedra derived from surfaces, x = {x 1 , ..., x 6 } ∈ R 6 , λ ∈ R.

Table 3 :

 3 Decomposition of polyhedra derived from linear features, x = {x 1 , ..., x 6 } ∈ R 6 , λ ∈ R.

Table 6 :

 6 Summary of the sums for the simulation of FC2

		Current method	Proposed method
	Polytope	Vertices Time [s] Vertices Time [s]
	Γ 1,0/2,1 ⊕ Γ 2,1/3,0	3304	5,6	512	0,7
	Γ 1,0/3,0 ⊕ Γ 3,0/4,0	10864	69,3	1728	16,3
	Γ 1,0/4,0 ⊕ Γ 4,0/5,0	34480	522,4	4096	160,4
	Γ 4,0/5,0-a ∩ Γ 4,0/5,0-b 4288	-	4288	-
	Γ 4,0/5,0 ⊕ Γ 5,2/5,0	48440	1207	29544	251,3