François Thiebaut

C Lartigue

E Duc

A certification method for the milling process of free-form surfaces using a test part

It is generally admitted that the manufacturing of free-form surfaces requires the use of a CAD-CAM system. The toolpath accuracy and the dimensional quality of the final shape have to be in accordance with the geometrical specifications. But most of the time, the final parts present deviations from the expected shape. These deviations may be due to either the toolpath calculation (CAM system) or the cutting process itself. In the paper, we propose an analysis of the whole milling process to point out the possible sources of errors. These errors generally lead to geometrical deviations and the final part does not meet the required specifications. As the errors can be linked to geometrical particularities of the shape, we propose a test part associated with check means to bring out problems. The milling of this part using two different techniques of toolpath generation shows that obviously both toolpaths are not errorfree and that errors result from different geometrical particularities of the part surfaces.

Introduction

Nowadays, high speed milling allows the realization of parts with no rectification thanks to the high quality of metal cutting [START_REF] Konig | New technologie in milling, high speed milling, 3 and multi-axis machining working of hardened parts[END_REF], [START_REF] Paro | High-speed milling in tooling production[END_REF]. The envelope of the toolpath essentially defines the final part shape. Thus, the toolpath has to be as accurate as possible since no polishing is made on the part.

The precision of the whole machining process has to lead to a low roughness, to a nice aspect of the part and of course, to the respect of geometrical specifications (form deviation, roughness, .....). The part dimensions are otherwise given by the geometrical model. But as the final quality is linked to the precision of all the steps of the machining process, the control of each of them yields the global control. As a result, the toolpath generation has to be good enough to ensure the competitiveness of the process.

An analysis of the milling process shows the different areas where problems may appear.

The milling process can be divided as follows (figure 1) :

-geometrical description, -surface modeling with a CAD system, -toolpath calculation with a CAM system, -program translation to be understood by the numerical controller via post-processor, -milling of the part on a machine tool.

The model resulting from the surface description by the CAD system is taken as the reference. As we are only interested in the three-axes milling with a hemispherical tool, the program translation via a post-processor is considered as perfect.

Thus, the precision of the milling process results from the accuracy of the toolpath generation and the cutting of the part on a machine tool.

In the paper, we propose a test part that allows to evaluate the accuracy of the calculated toolpath and that of the cutting process separately. figure 1 : the milling process of a free-form surface We essentially present the approach that leads us to elaborate the test part. A list of the different sources of errors or inaccuracies is established, separating errors coming from unexpected behaviors of the toolpath generation with a CAM system from those coming from the machining process, e.g., C 1 discontinuities on the part may provide errors in the toolpath calculation. To each kind of errors corresponds one or more surfaces on the test part. In other words, each surface of the test part is elaborated so as to make either the toolpath calculation or the part cutting potentially fail.

As our objective is to control the milling of various kinds of parts on the basis of the milling of the test part, checking means are set up that allow to evaluate the errors or the inaccuracies. In particular, we have based our approach on the evaluation of the geometrical deviation between the machined surface and the nominal surface at each point. This allows to verify if the machined surface meets the main used geometrical specification for sculptured surfaces : the form deviation.

The sources of error in the global process

The aim of manufacturing with a CAD system is to produce a part that meets the required specifications. The most relevant tolerance specification for sculptured surface is the form deviation (profile tolerance) [START_REF] Austin | Comparison of discretization algorithms for Nurbs surfaces with application to numerically controlled machining[END_REF]. The part is said within tolerance if the machined surface is inside two envelope surfaces given by sweeping a sphere of radius t along the nominal surface (ISO 1101(ISO , 1983) (figure 2). In addition to the form specification, one can specify the roughness of the finished surface too. figure 2 : tolerance of form deviation Commonly, in 3-axis milling with a ball-end cutter tool, the calculated toolpath with a CAM system is a polyhedral line which represents the centers or the extremities of the hemispherical tool. Obviously this sampling of the shape induces dimensional deviations between the final shape and its geometrical representation. But, as we will see next, only a few toolpath generation techniques rely on the geometrical specifications.

CAM system errors

Each calculation technique uses its own approximation criteria to generate the toolpath. This involves different calculated toolpaths for a unique CAD model and hence various machined shapes. Moreover, the discretization techniques may lead to a final shape that does not meet the required specifications. Calculation techniques rely on the evaluation of the successive positions of the tool which are tangent to the original surface. At each point of the surface, the tangent location of the tool verifies (figure 3):

CL   CP   R n  ( 1
)
where CP is the contact point, CL the cutter location and n  the local normal. The whole toolpath is then defined by the set of the CL points for which two successive CL points are linked together by a straight line. To achieve this calculation, authors first evaluate the offset surface of the original surface with the tool radius for offset value. The offset surface may be determined in an explicit manner or an implicit. For that last case, the problem is solved for a discrete number of points.

Classically, toolpath generation techniques can be divided into four main types.

The first one consists in generating an offset of the part whose sampling gives the toolpath [START_REF] Kim | A new machine strategy for sculptured surfaces using offset surface[END_REF], [START_REF] Kim | Verification of NC tool path and manual and automatic editing of NC code[END_REF], (Lai J.Y, Wang D.J, 1994), [START_REF] Tang | Offsetting surface boundaries and 3-axis gouge-free surface machiningcomputer-aided design[END_REF].

If Q  (u , v )
is the surface equation, the offset is given by :

Q off  (u , v )  Q  (u , v )  R n  (u , v ) (2)
where R is the tool radius and n  (u , v ) the normal to the surface at each point.

In that case, the main difficulty is that the offset surface is completely defined only if the normal is known at each point. We can notice different sources of errors presented in figure 4 :

-when the tool radius is locally greater than the curvature radius of the surface, self intersections generally appear, -shapes such as "butterfly" provide an inversion of the normal that makes the calculations fail, -C 1 discontinuities on the surface, for example sharp edges, lead to C0 disontinuities on the offset surface.

figure 4 : geometrical singularities for the normal evaluation Such calculation techniques must be combined with surface analysis to define areas that require corrections but this is not yet completely achieved.

The second consists in determining the tool contact positions and in deducting their center locations [START_REF] Bobrow | NC machine toolpath generation from CSG part representations[END_REF], (Choi B.K, Jun C.S, 1989). In this case, the obtention of the offset surface is implicit. The problem is defined by (1) and has to be solved after the discretization of the original surface in contact points. Of course, a wrong sampling of the tool contact locations may lead to the risk of missing out small dimension surfaces. In addition to the errors linked to a bad definition of the normal at each point, a possible source of errors is the interferences between the tool and the surface (figure 5).

figure 5 : interference tool-surface

The third type consists in machining a polyhedral approximation of the surface [START_REF] Hwang | Interference-free tool-path generation in the NC machining of parametric compound surfacescomputer-aided design[END_REF]. The determination of the cutter locations is then easier for realized from elementary facets. The precision of the toolpath generation is given by the sampling step of the polyhedral surface and if the step is too large small details may be missed out (figure 6).

figure 6 : offset evaluation from polyhedral approximation

The last technique is the determination of the toolpath with the method of the inverse offset [START_REF] Suzuki | Development of the CAD/CAM System based on parallel processing and inverse offset method -Transputing'91[END_REF]. The main interest of such a technique is that the offset calculation is very simple and avoids undercuts since the calculation takes the orientation of the tool axis into account, as well as the whole working part of the tool. The main difficulty is for the milling of vertical walls that are parallel to the tool axis.

Different approaches can be taken to specify the range of the authorized deviations. Some CAM systems allow the machined surface to be included between two offsets of the model shape. More commonly, others allow the specification of criteria to determine the pass computation. These criteria may be the chordal deviation in the driving tool direction and the path interval in the perpendicular direction [START_REF] Kim | A new machine strategy for sculptured surfaces using offset surface[END_REF]. Generally, to compute the distance between two CL data in the driving direction, authors locally approximate the surface by a constant curvature surface, most often a sphere [START_REF] Choi | [END_REF]. This assumption is widely violated when a sharp edge is treated. The path interval is generally used to compute the discretization step between two successive paths. [START_REF] Lin R | Efficient tool-path planning for machining free-form surfaces -Transaction of ASME[END_REF] present the models used by most of the authors. Most of the time, the path interval is calculated from the maximum scallop height allowed. Works rely on the following assumptions : tool locations are tangent to the surface ; the surface is continuous and can be approximated by a circle in a perpendicular plane to the driving direction. These hypotheses are not always verified and, due to the chordal deviation, the actual scallop height is greater than the authorized one [START_REF] Tournier | Evaluation des modes de génération des trajets outil pour l'usinage de formes gauches. Etude de cas[END_REF].

The toolpath will be considered as error-free if the whole set of CL data allows the respect of the geometrical specifications. At the end of this analysis, we can conclude that techniques of toolpath generation may introduce errors due to the different approximations that have to be done; such errors are generally not easy to evaluate or to correct. Moreover, discretization parameters are rarely linked to the geometrical specifications. As a result the final surface may not meet the required specifications. Austin et al propose a discretization technique that takes the geometrical specification into account : the edge length of the triangular facets is calculated so as to comply with the profile tolerance.

Errors of the machining process

The machining process is also a source of the dimensional difference between the final shape and its geometrical definition. The errors may come from the treatment by the numerical controller, the rigidity of the milling machine and the cutting process [START_REF] Makino | Clothoidal Interpolation, a new tool for hish speed continuous path control[END_REF], [START_REF] Weck | Sharp corner tracking using the IKF control strategy[END_REF]. High speed milling requires high speed cutting and high feedrates, thus a high speed treatment of the toolpath.

If the cutting process is not the same on the whole shape, it may induce differences of roughness on the finished surface. It is the case with a three-axes milling because the orientations of the tool axis relative to the shape surface involve a non-constant cutting process.

An abrupt variation of the feed-rate causes marks on the shape due to the inertia of the machine tool. This variation of the feed-rate comes from the numerical controller for different reasons:

-the transmission speed of the data is less than the treatment speed, so the machine slows down while the numerical controller is waiting for information [START_REF] Schultz | High speed machining[END_REF], (Yeung M.K, Walton D.J, 1994),

-the milling time of an elementary toolpath is smaller than that allowed by the numerical controller, so the machine slows down (Yeung M.K, Walton D.J, 1994),

-the rapid variation of the normal orientation on the shape requires changes of the tool movement, hence some accelerations on the different axes.

Unfortunately, these kinds of deviations can only be prevented by a look ahead function of the numerical controller, and we can only notice their effect at the end of the cutting process.

The test part

For an efficient checking of the process, the test part must be representative of the CAM problems and the machining process problems. The main problems are, as seen previously, the treatment of the small curvature radius, the normal singularities on the shape, the length of an elementary toolpath, the rapid evolution of the normal orientation and the changes in the process cutting. We have to consider the details with small dimensions and the undercut shapes too.

These different sources of deviations can be shown using surfaces which present geometrical particularities. They are classified in seven types:

1 -small curvature radius on concave surface, 2 -geometrical discontinuities of the surfaces: gap, hole and overlapping, 3 -singularities in the evaluation of the normal : undefined normal, normal inversion, 4 -undercut shape, 5 -surface with great curvature radius, 6 -small dimensions, 7 -sharp edge between surfaces in convex relation and vertex.

A solution to analyze the machining process is to get a geometrically defined shape which contains all the necessary elements to bring out the behavior of the whole process. Thus, the test part has to be geometrically defined and easy to check. As the existing test parts do not take all of these specifications into account, we propose one that considers these issues [START_REF] Thiebaut | Conception d'une pièce test pour la FAO et l'usinage de formes[END_REF].

The part is mostly composed of canonical but representative surfaces in order to minimize the potential deviations of the CAD system (assumed to be perfect in our work). Of course, the surfaces could be modelled by Nurbs surfaces or other models without limiting the proposed approach. As far as possible, the choice of the surfaces is realized to show only one of the problems for each surface. For most cases, the size of the geometrical characteristics is defined relatively to the used tool radius (R=5 mm). A representation of the resulting part is given with the recapitulation of the different surfaces that compose it figure 7.

figure 7: the test part Type 1 particularities are analysed by the means of three different kinds of grooves in order to essentially determine the behavior of the toolpath generator. The groove located near the details presents a radius that is smaller than the tool radius and allows to know if some interferences occur between the shape and the tool. The width of the second is exactly equal to the diameter of the tool, so the tool is tangent on both sides. The evolutional width of the third one can provide the minimum required curvature radius of a machinable shape (figure 8).

figure 8: a specifical shape, the grooves Type 2 singularities are directly brought out with the crossing faces, the hole, the gap and the overlapping. These discontinuities are generated with planar surfaces in order to minimize the other deviations that may occur during the machining of biparametric surfaces.

Type 3 particularities require surfaces whose the normal is not defined, inversed or not continuous. We take these mathematical degeneration into account on triangular patches where the normal is not defined, on one vertex, by the use of the tetrahedron. The "butterfly"

surface is a patch defined by two parametric segments whose parametric orientations are opposed. The result is a surface presenting an inversion of its normal. The last differentiation concerns the curvature continuity between linking surfaces. The most representative associated surfaces are the sequences from plane with torus and torus with sphere.

Type 4 particularity is analyzed by the means of the undercut shape.

The smooth surface combines all the type 5 particularities. The association of two revolution surfaces with a ruled surface, allows to know the differences of the behavior of both the toolpath generator and the cutting process in function of the curvature and the concavity.

Type 6 particularities are created with cuboids of small dimensions. These cuboids are located either on vertical or horizontal faces, thus it is possible to know if details are either missed by the extremity of the tool or by its side.

The main objective of the type 7 surfaces is to analyze the process cutting and more specifically the behavior of the numerical controller. Indeed, the easiest way to determine the pursuit deviation (characteristic of the numerical controller) is to produce it between planar surfaces: the right cylinder and the parallelepiped. For such surfaces, we expect the CAM system to generate a toolpath without deviation. Thus, if there are deviations on the machined surface, they entirely result from the pursuit deviation of the numerical controller.

The main characteristics of the test part are now defined. The surfaces can also be used to study the influence of the tool axis orientation relative to the surface on the quality of the cutting process, for example: the right cylinder presents some associated surfaces with the same orientation relative to the tool axis but that are milled "upwards" and "downwards".

Moreover, an analysis can be made in five-axes milling to determine the capability of the toolpath generator to propose a toolpath with no useless rotation of the spindle or to know the influence of a low and permanent rotation of the spindle on the final part.

Let us now see the checking means that are associated with this test part and allow to evaluate the different kinds of deviations.

Checking means to assess the manufacturing process

CAM system

The checking of the toolpath is based on the fact that it is known and on the knowledge of its geometrical shape. The existing methods which allow to compare the geometrical part and the toolpath are either approximation methods that use a sampling surface: the Z-buffer method [START_REF] Kim | Verification of NC tool path and manual and automatic editing of NC code[END_REF], the point vector method [START_REF] Chappel | The use of vectors to simulate material removed by numerically controlled milling -computeraided design[END_REF], (Jerard R.B et al, 1989a), (Jerard R.B et al, 1989b) and the visual method [START_REF] Suzuki | Development of the CAD/CAM System based on parallel processing and inverse offset method -Transputing'91[END_REF] or exact methods : CSG representation of the toolpath (Jerard R.B et al, 1989a), [START_REF] Hanada | Development of automated toolpath interference check system by using solids models[END_REF]. This method is interesting but the computations lead to great model sizes, and this approach is only used to analyze local deviations.

In the paper, the assessment of the manufacturing process is realized relative to the form deviation : the machined surface has to lie inside the two envelope surfaces defined by sweeping a sphere of radius t along the nominal surface. The machined surface is only known through a set of points and the assessment is realized from the evaluation of the geometrical deviation at each point, which is the distance between a point on the machined surface and the nominal one.

Relative to the CAM system and the toolpath calculation, the geometrical deviation is calculated in two steps. First the point on the machined surface is determined by the use of the Z-buffer method. From a point on a cartesian grid defined on the nominal surface, the machined point is determined by the intersection between the vertical line and the envelope toolpath (figure 9).

figure 9: calculation of the geometrical deviation For an elementary toolpath (straight line between two CL data), the enveloppe toolpath is materialised by a cylinder limited by two half-spheres. The intersection gives the Z M value of the machined point for which (X M , Y M ) are defined by the grid. M belongs to the vertical line and M belongs to the cylinder, so M verifies :

C m M   C M M  C m C M   R (3)
with M such as :

C m C M   C m M   0 C m C M   C M M   0      (4)
If ( 4) is not verified we have to compute the intersection with the half spheres.

The step of the grid is obviously an essential element since it characterizes the number of assessment points. In the second step, the geometrical deviation is calculated, i.e. the distance of the machined point to the nominal one.

The distance between each point M which belongs to the toolpath and the machining face (which is a re-limitation of the surface) is given by ( 5), where P is the projection of the M point on the face :

e  PM   n  where n  is the surface normal (5)
The toolpath is said error-free if the geometrical deviations are all in the authorized range.

Machining process

At first, a visual analysis allows to detect the potential presence of marks on the final part.

In order to identify the errors due to the machining process. The final part has to be measured by use of a Coordinate Measuring Machine. In order to elaborate a collision-free probe path, the movements of the probe can be built from the toolpath. The position of the measuring points are selected among the CL data and the measuring direction is given by the normal to the surface at each point. At this stage, the comparison with the nominal model provides the errors due to the CAM system and the machining process. The specific CAM system errors may be assessed by the means previously exposed. Thus, the machining process error can be determined by difference. The knowledge of the failing of each chain link during the entire process obviously allows the checking of the whole process.

Validation of the certification method

In order to validate the certification method, we have carried out the entire milling process for the test part. Two different toolpath generation techniques have been tested. The CAM parameters, i.e, the machining direction, the maximum chordal deviation, Mt, and the maximum scallop height, hc, are given in table 1.

Machining direction Maximum chordal deviation

Maximum scallop height

Parallel to the plane xOz 0.01 mm 0.01 mm table 1: CAM parameters These parameters are not directly linked to the geometrical specifications. But a previous analysis [START_REF] Tournier | Evaluation des modes de génération des trajets outil pour l'usinage de formes gauches. Etude de cas[END_REF] has shown that the respect of such parameters implies that the machined surface is included between two offset surfaces of the nominal one, for which the -on the butterfly shape the normal inversion is not correctly treated : the adjoined plane is machined at the same altitude and as a result the geometrical deviations are equal to 0.5mm, -in the case of concave surfaces machined by the tool side we find numerous overcut points essentially in areas where the curvature radius is slightly greater than the tool radius;

the geometrical deviations reach 0.09mm, -on the sides defined by X= Xmin and X = Xmax, the geometrical deviations reach -0.12; this provides undercut points on the part boundaries without apparent reason.

We can notice that for the technique b) the right cylinder presents a CAD modeling error which implies a systematic geometrical deviation of 0.1mm. To conclude, the treatment of the grooves and the undercut shapes is different from the two techniques but both are error-free.

The second step of the validation of the certification method is the milling of the test part and its visual analysis. The realization is carried out from the technique a) toolpath which presents less calculation errors. The visual analysis of the part brings out local marks and systematic marks wherever the curvature radius is less than the tool radius. Such kinds of marks, probably due to the cutting process, should disappear for a high speed cutting. On the other hand, most details are missed out.

Conclusion

The certification method for the milling process of free-form surfaces we propose is

figure 3 :

 3 figure 3 : definition of the cutter location, CL data

  offset values are given by -Mt and Mt + hc. So, the toolpath generation technique is assessed if the geometrical deviation at each machined point belongs to the range [-Mt, Mt + hc].

Figure 10

 10 Figure 10 presents a map of the geometrical deviations that means the points on the part

  figure 1 : the milling process of a free-form surface figure 2 : tolerance of form deviation


M Q Nguyen

O Sename

L Dugard
email: luc.dugard@gipsa-lab.grenoble-inp.fr

An LPV Fault Tolerant control for semi-active suspension -scheduled by fault estimation

Keywords: semi-active suspension, fast fault estimation, LPV/H ∞ control, fault tolerant control

In this paper, a novel fault tolerant control is proposed to accommodate damper faults ( oil leakages) in a semi-active suspension system based on a quarter-car vehicle model. The fault accommodation is based on the Linear Parameter Varying (LPV) control strategy and involved in 2 steps. At first, a fast time-varying fault is estimated by using the fast adaptive fault estimation (FAFE) algorithm and based on an unknown input adaptive observer. Thanks to information about the estimated fault, the dissipativity domain of the semi-active suspension is adapted according to the fault. Then a single LPV fault tolerant controller is developed to manage the system performances. The controller solution, derived in the LPV/H ∞ framework, is based on the LMI solution for polytopic systems. Some simulation results are presented that show the effectiveness of this approach.

INTRODUCTION

Automotive vehicles are extremely complex systems composed of many interrelated subsytems that, in particular enhance the overall driving comfort, stability and safety (see [START_REF] Gillespie | Automotive control systems: for engine, driveline, and vehicle[END_REF], [START_REF] Gillespie | Automotive control systems: for engine, driveline, and vehicle[END_REF]). Among all sub-systems impacting on the vertical vehicle dynamics, the supension systems play a key role since they ensures the link between the wheels and the chassis (see [START_REF] Fischer | Mechatronic semi-active and active vehicle suspensions[END_REF]). A well designed suspension system may considerably improve not only the passenger comfort but also the car road holding. Recently, the semi-active suspensions have been more and more widely used in automobile industry since they achieve the main performance objectives while they are smaller in weight and volume, cheaper in price, more robust and less energy consuming. Several control design problems for semi-active suspension systems have then been tackled with many approaches during the last decades. In the works of [START_REF] Savaresi | Semi-active suspension control design for vehicles[END_REF], [START_REF] Poussot-Vassal | Survey and performance evaluation on some automotive semi-active suspension control methods: A comparative study on a single-corner model[END_REF],the authors presented several control strategies for semi-active suspensions (based on the Skyhook, Groundhook, ADD). Moreover, to cope with the dissipativity constraints of semi-active dampers, some control approaches using the LPV techniques have been presented. In [START_REF] Do | An LPV control approach for semi-active suspension control with actuator constraints[END_REF][START_REF] Do | LPV modelling and control of semi-active dampers in automotive systems[END_REF], the nonlinearities of the semiactive suspension are taken into account and written in LPV form using the saturation of the control input. This method has been extended in [START_REF] Sename | Lpv methods for fault-tolerant vehicle dynamic control[END_REF]) in case of damper faults but using 3 varying parameters which could be conservative. In [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF]), a kind of LPV gain scheduling anti-windup strategy has been proposed to handle such a constraint. This approach will be extended here in the case of some damper malfunctions. About the fault tolerant control (FTC), in general there a 2 main groups of FTC: passive (off-line designed) using the robust control with respect to possible system fauts and active (online control reconfiguration mechanism) using a FDI module and accommodation techniques. In more detail, a FDI strategy can be model-based or data-based and is used to detect, isolate as well as estimate the faults. During the last decades, FDI modules, based on the analytical redundancy for fault estimation (e.g in Zhang and Jiang (2008)), have been received a lot of attention by many reseachers.

Fault estimation is a key step in designing a fault tolerant control. In the litterature, there exist many different approaches to estimate a fault which can be either actuator or sensor malfunction. Let us mention the classical methods, based on the parity space theory (see in [START_REF] Gertler | Fault detection and isolation using parity relations[END_REF]) to generate the residues and approximate the fault or the bank of observers approach (see [START_REF] Varrier | Détection de situations critiques et commande robuste tolérante aux défauts pour l'automobile[END_REF]) as well as by sliding mode observers [START_REF] Edwards | Sliding mode observers for fault detection and isolation[END_REF]). Recently, a new approach (see in [START_REF] Shi | Active fault tolerant control of lpv descriptor systems based on extended state observers[END_REF]) considered the fault element as a state of the augmented system and designed an extended observer to estimation at the same time the state and the fault of system. However, it is limited to constant faults ḟ (t) = 0. Then, Zhang et al. (2008) presented a method allowing to evaluate the time-varying fault by using a fast adaptive fault estimation (FAFE) methodology based on an adaptive observer. But therein, the authors solved the problem with a regular LTI system without considering the disturbances. Next, [START_REF] Rodrigues | Actuator fault estimation based adaptive polytopic observer for a class of lpv descriptor systems[END_REF] proposed an adaptive polytopic observer for timevarying fault estimation in spite of the disturbances for a class of descriptor LPV systems.

The main purpose of this paper is to propose an active FTC for a semi-active suspended quarter-car vehicle which suffers from oil leakage. The paper contributions are twofold:

• First, the fault estimation is based on the FAFE algorithm allowing to determine the magnitude of the fault. This approach is used for the first time to estimate the actuator fault in a semiactive suspension system while taking into consideration the unknown input road disturbance.

• Then an LPV fault scheduled controller is developed to ensure the damper dissipativity property and enhance passengerss comfort and road holding (in normal and faulty conditions). This controller is an extension of the one in [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF] where, using the obtained fault information, the anti-windup kind strategy is adapted according to the estimated fault in order to guarantee the semi-activeness of the damper.

The paper is organized as follows: the next section describes about the problem statement. Section 3 presents a method to estimate the fault on damper. Section 4 is devoted to the design of a fault tolerant LPV semi-active suspension control. Some the results of the proposed method are given in section 5. Finally, some conculsions are drawn in the section 6.

PROBLEM STATEMENT

In this work, a simple quarter vehicle, depicted in Fig. 1 is taken into account: This picture presents a single corner of a vehicle. . Concerning the semi-active suspension control, the main challenge is to take into account the dissipativity of the damper and the saturation in the synthesis step. If this dissipativity constraint is not considered, it is necessary to "saturate" the control input without any performance and stability guarantees, which is referred to as the clipped strategy in [START_REF] Savaresi | Semi-active suspension control design for vehicles[END_REF]). In [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF]), the considered semi-active damper is simply modeled as a static map of the deflection speed/Force, i.e a lower bound and upper bound of the achievable forces as shown on Fig. 2. This static model is thus a saturation function of the deflection speed, denoted as the dissipative domain D(ż de f ). Moreover, a smart parameter is introduced allowing to take the real abilities of the damper into account. This scheduling parameter is indeed defined as a function of the difference "ε" between the computed damper force F d (given by the controller) and the achievable one F ⊥ d , that was used to satisfy the dissipative damper constraints. Therefore, the scheduling parameter depends on D(ż de f ).

Fig. 2. Dissipative domain D graphical illustration

As shown in [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF], this method gives very good results nominal condition i.e without faults. However, if the damper is subject to faults e.g. oil leakages, the method becomes no longer appropriated. Indeed, for a faulty semi-active suspension, the available damping force (low and high level) is lower than a heathy damper and consequently the deflection motion increases. It means that the dissipative domain will change in the presence of fault. Therefore, this domain depends on both the deflection speed and the fault, and so-called D f (ż de f , f ) as in Fig. 3 with f stands for a fault to the damper.

Fig. 3. Dissipative domain D f in presence of fault

As a consequence, if the saturation constraint is not adapted, then the dissipativity condition is not be guaranteed. Indeed, the required force could be outside of the range of the "real" faulty force (even if valid inside of the range of the healthy force). In this case, the control performances are not ensured if some fault information is not included into the control design. To deal with this problem, in this work, a LPV fault tolerant control is proposed to include the dissipativity constraints of the semi-active suspension in case of malfunction in order to achieve good performances.

Problem formulation

Let us consider now a fault on the semi-active damper e.g an oil leakage which induces a lack of force modeled as:

F sa = α(c 0 żde f + u) = c 0 żde f + u + (α -1)(c 0 żde f + u) = c 0 żde f + u + f (3)
where F sa stands for the fault force expressed as a reduction of the nominal semi-active force and α ∈ [0 1] is the oil leakage degree, e.g. α = 0.8 means that the damping force will be of 80% of the nominal damper force F sa due to a lost force f of 20%. f = (α -1)(c 0 żde f + u) stands for the fault.

Then the state space representation of the vertical dynamic using a quarter car model and taking into account a faulty semiactive damper, is given as follows:

ẋ = Ax + B 1 w + B 2 u + E f (4) y = Cx where x = (z s , żs , z us , żus ) T ∈ R n is the state vector, w = z r is the disturbance input, u ∈ R m is the control input, y = [z s -z us , żs - żus ] = [z de f , żde f ] ∈ R p
is the output vector and f ∈ R r stands for actuator faut.

A =                       0 1 0 0 -k s m s -c 0 m s k s m s c 0 m s 0 0 0 1 k s m us c 0 m us - k s + k t m us - c 0 m us                       , B 1 =                  0 0 0 k t m us                  , B 2 = 0 -1 m s 0 1 m us T , E = 0 -1 m s 0 1 m us T , C = 1 0 -1 0 0 1 0 -1 .
This model will be used in the next section in order to estimate the fault actuator f based on an adaptive observer. Then from the estimated fault, a LPV FTC will be developed and allows a fault accommodation.

FAULT ESTIMATION

This part shows the procedure to estimate the fault actuator on a damper using an adaptive observer. This is developed based on the method in Rodrigues et al. ( 2014) (which is used for the descriptor LPV systems).

Consider now the state space repesentation (4), if the following assumptions are satisified:

• Assumption 1: rank (CE)= r.

• Assumption 2: The matrix (A,C) is observable • Assumption 3: The fault f (t) and the derivative of f (t) with respect to time are norm bounded i.e: 0 ≤ f (t) < α 1 and 0 ≤ ḟ (t) < α 2 with 0 ≤ α 1 , α 2 < ∞.

The state space repesentation (4) includes the road disturbance w = z r which is clearly the actual unknown input in the suspension system. Therefore, while estimating the fault actuator f , one needs to take into account this unknown input. To deal with this problem, an unknown input adaptive fault observer is proposed and given as the following structure:

ż = Nz + Gu + Ly + E f x = z + T 2 y (5) ŷ = C
x where z ∈ R n is the state variables of the observer, x ∈ R n the estimated state variables. ŷ ∈ R p is the estimated output vector and f (t) ∈ R r is the estimation of the damper fault f (t). N,G, L, T 2 are the observer matrices to be designed.

Denote : e x (t) = x(t)x(t), e y (t) = y(t)ŷ(t), e f (t) = f (t)f (t) are state, output and fault estimation error respectively. If

T 1 = I -T 2 C, then e x = T 1 x -z and ėx = T 1 ẋ - ż = T 1 (Ax + B 1 w + B 2 u + E f ) -(Nz + Gu + Ly + E f )
After some manipulations, one has:

ėx = Ne x + Ee f + (T 1 E -E) f + + (T 1 A -NT 1 -LC)x + (T 1 B 2 -G)u + T 1 B 1 w
If the following conditions hold:

T 1 A -NT 1 -LC = 0, T 1 B 2 -G = 0, T 1 B 1 = 0, M := T 1 E -E Then,
the estimation error dynamic is given by:

ėx (t) = Ne x (t) + Ee f (t) + M f (t), (6)
e y (t) = Ce x (t) The convergence of the state estimation error (6) can be verified by the following theorem:

Theorem 1: (see [START_REF] Rodrigues | Actuator fault estimation based adaptive polytopic observer for a class of lpv descriptor systems[END_REF]) Under these assumptions 1-3, given scalars σ, µ > 0, if there exist symmetric positive definite matrices Q, P 1 , P 2 , P 3 , and matrices W, U such that the following conditions hold:

         * - 1 σ E T (QT 1 A + WC) -2 1 σ E T QE + 1 σµ P 2 + 1 σµ P 3          < 0 (7) and E T Q = UC (8) where = (T 1 A) T Q + Q(T 1 A) + C T W T + WC + 1 µ P 1 and *
stands for the symmetric elements in a symmetric matrix.

Then the following fault estimation algorithm: ḟ (t) = ΓU(ė y (t) + σe y (t)) (9) can realize e x (t) and e f (t) uniformly bounded, where Γ ∈ R r×r is a symmetric positive definite learning rate matrix.

Then it is easy to show that the estimated fault can be deduced from the expression (9) as follows :

f (t) = ΓU        e y (t) + σ t t f e y (τ)dτ        ( 10
)
where t f is the time since fault occurs. The fault reconstruction in (10) combines a proportional term with an integral one that allows to improve the rapidity of fault estimation. Therefore, to obtain the fault estimation, one needs to solve conditions in Theorem 1 which consists of a Matrix Inequality (7) and an equality matrix constraint (8) which are transformed into following optimization problem:

Min γ subjects to (7) and

γI E T Q -UC (E T Q -UC) T γI > 0 (11)
Solving this optimization problem, we obtain the values of observer matrices N, L,G, T 2 and the estimation of the fault actuator.

FAULT TOLERANT LPV SEMI-ACTIVE SUSPENSION CONTROL

In this section, an LPV fault scheduled state feedback controll is proposed to ensure the damper dissipativity and keep good dynamic performances of the faulty semi-active suspension system. Of course, some gracefull performance degradations are allowed.

Fig. 4. General block diagram

The overall controll structure is presented in Fig ( 4). In this work, the LPV static state feedback controller K receives the state variables x as an input and computes the damping forces u to be added to the nominal damping forces c 0 z de f in order to improve the vehicle performances. This controller is scheduled by the parameter ρ that constraints the control signal or not, in such a way that the required forces F d remain semi-active and adapted to damper ability. The fault estimation algorithm, proposed in the last section, allows to estimate the fault in damper. This estimated fault is used to modify the dissipativity domain D f (ż de f , f ) of the semi-active suspension allowing to schedule parameter ρ. Some weighting functions are then taken into account in the controller synthesis to improve the performances of the vehicle:

• W z r is used to shapes the road disturbances effects z r

• W z s = k z s s 2 +2ξ 11 Ω 1 s+Ω 1 2 s 2 +2ξ 12 Ω 1 s+Ω 1 2 relates to the comfort perfor- mances • W z us = k z us s 2 +2ξ 21 Ω 2 s+Ω 2 2 s 2 +2ξ 22 Ω 2 s+Ω 2
2 relates to the road-holding performances • W u (ρ) = ρ is used to penalize (more or less) the control input signal amplification according to the ρ signal. More specifically, it is used to guarantee the semi-activeness (see next subsection).

Remark:

The weighting functions W z s , W z us are chosen thanks to the genetic algorithm as in [START_REF] Do | Multi-objective optimization by genetic algorithms in h/lpv control of semi-active suspension[END_REF] .

Scheduling parameter

The method proposed in the previous work of [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF] for a quarter car in order to fulfill the dissipativity constraint, aims at increasing or decreasing the gain of the weighting filter W u on the damper control signals, according to a given scheduling strategy. Indeed, if the required force computed by the controller is active, a scheduling parameter allows the controller to enhance or not the performance specifications, so that the required force remains dissipative. This method is extended here in case of a faulty damper. Let us now define the clipping function for the faulty case:

Definition of Fault-scheduled clipping function:

Due to the controlled damper limitations (i.e. the effective force provided by the damper F d should lie in the dissipative domain), the following Fault-scheduled clipping function D f (F d , żde f , f ) is defined (see also illustration in Fig. 3) as:

D f (F d , żde f , f ) → F d = F d if F d ∈ D f F ⊥ d if F d D f (12)
where F d is the required force (given by the controller) and

F ⊥ d is the orthogonal projection of F d on D f .
This definition will inspire the form of the considered scheduling parameter used in the LPV control.

In more detail, the ρ parameter is tuned as following:

• when ρ is low, W u (ρ) is small and it does not penalize the control signal u. • when ρ is high, W u (ρ) is large and it attenuates the control signal u to remain in semi-acitve domain.

For that purpose the following scheduling strategy ρ(ε) is introduced: 14) where ε is the distance between the required force and the force projected on the adapted dissipative domain D f (according to the function of D f (F d , żde f , f ). µ is a design parameter that modifies the dead-zone of the ρ(ε) function (µ is chosen sufficiently low e.g µ = 0.1).

ρ(ε) :=                ρ if ε < µ ρ + ρ -ρ µ (ε -µ) if µ ≤ ε ≤ 2µ ρ if ε > 2µ (13) ε = ||F d -F ⊥ d || 2 (
Remark:

• By this definition, ρ(ε) belongs to [ρ, ρ] which is essential in the LPV framework (ρ = 0.01, ρ = 1). • ε 0(⇔ F d F ⊥ d )
means that the required force is outside the allowed range. Conversely, ε = 0(⇔ F d = F ⊥ d ) means that the force required by the controller is reachable for the considered semi-active actuator.

This varying parameter has been used to schedule the designed static state-feedback vehicle controller.

H ∞ /LPV control design for FTC

It is worth noting that, while the model car is a LTI system, the generalized plant (which consists of the suspension model and weighting functions) is a LPV one because of parameter dependant weighting function W u (ρ). Then the following parameter dependent suspension generalized plant (Σ v (ρ))is expressed in: The generalized plant (15) depends on the varying parameter (ρ), so Σ v (ρ) can be expressed as a polytopic system composed by N = 2 vertices ω i , i=1, 2:

Σ(ρ) :          ξ = A(ρ)ξ + B 1 (ρ) w + B 2 u z = C 1 (ρ)ξ + D 11 (ρ) w + D 12 u y = Cξ
Σ v (ρ) = |ρ -ρ| ρ -ρ Σ v (ρ) + |ρ -ρ| (ρ -ρ) Σ v (ρ) (16)
where Σ i v (ρ) defines the system at i th vertice; Now, let denote K is the state feedback controller of the generalized plant Σ v (ρ) such that the control input: u H ∞ = K.x. Then, thanks to the Bounded Real Lemma (BRL), finding such a controller leads to solve the optimization following (see in [START_REF] Scherer | Multiobjective output-feedback control via lmi optimization[END_REF]):

Min γ such that:

          (A + B 2 K)P + P(A + B 2 K) T B 1 P(C 1 + D 12 K) T * D T 11 * * -γ 2 I           < 0
(17) where P is a symmetric positive definite matrix. The condition (17) is a Bilinear Matrix Inequality (BMI), by denoting Q = KP, one obtains the following LMI:

          (AP + B 2 Q) + (AP + B 2 Q) T B 1 (C 1 P + D 12 Q) T * -I D T 11 * * -γ 2 I           < 0 (18)
Solving this problem, one receives the state feedback controller

K = QP -1 .
It is important to note that the suspension controller is designed using H ∞ /LPV framework for the polytopic system. Then the LPV-FTC controller is a convex combination of the controllers computed at each vertex, so the control input can be expressed as:

u H ∞ =        |ρ -ρ| ρ -ρ K(ρ) + |ρ -ρ| (ρ -ρ) K(ρ)        x (19)
Since the LMI problem is solved at each vertex of the polytope formed by the limit values of the varying parameter, the stability will be guaranteed for all trajectories of the varying parameter.

SIMULATION RESULTS

To The following scenario is used to test the performance of the proposed LPV/H ∞ fault tolerant control:

• The vehicle runs at 30km/h in a straight line on a dry road (µ = 1 stands for the adherence to the road).

• The vehicle has a faulty semi-active damper because of a oil leakage, 50% of reduction of the nominal damping force (α = 0.5), occurs at t=0. • a 3cm bump on the wheel from t = 1 to 1.5s.

First, the actuator failure has been estimated using the procedure presented in section 3. Fig. 5 shows the estimation of oil leakage degree α = F sa -f F sa which represents a lost force of 50% with respect to the nominal damper force. Here the given simulation results aim at proving that the faultscheduled LPV strategy to handle a damper failure, are improving the performances of the previous LPV semi-active control strategy in [START_REF] Poussot-Vassal | A new semi-active suspension control strategy through lpv technique[END_REF]. The result will show a comparaison between a classice semiactive LPV controller without fault tolerance features, denoted as "LPV nominal without FTC" in blue, and the proposed LPV-FTC in red. Figure 6 shows that in case of faulty damper, the LPV-FTC allows to mitigate the sprung mass displacement with respect to the LPV nominal without FTC. Indeed, when the vehicle occurs a bump, the suspension force is required to reduce this motion. Because of a reduction of the nominal damping force, the LPV without FTC cannot compensate the fault. On the other hand, by a FTC strategy, it allows to reconfigure the damper force and attenuates the motion. The comfort of the passengers is then increased in contrast to the not adapted case. The unsprung mass displacement and the suspesion displacement are plotted in Figure 7, 8. The LPV-FTC allows to reduce the unsprung mass motion as well as the relative deflection with

CONCLUSION

The paper has presented a new LPV/H ∞ fault tolerant control for a faulty semi-active suspension system. When there exist a fault in damper such as oil leakages, the fast adaptive fault estimation algorithm is used to estimate the damper fault. Then using the obtained fault information, the LPV FTC is designed to ensure the damper dissipativity constraint. This FTC strategy allows to online reconfigurate provided suspension force according to fault situation, in order to achieve the designed performance objectives in both comfort and road holding. Moreover, the Fault Tolerant Control has been designed in the LPV/H ∞ framework allowing to simplify the implementation procedure. The next step of this work will be the implementation of this strategy on a test benchmark, available at Gipsa-lab Grenoble, developed in collaboration with a high-tech start-up "SOBEN". It consists of a vehicle equipped with four controllable Electro-Rheological dampers, and of 4 DC motors generating separately different road profiles on each wheel. First experimental results on the test-bed are presented in [START_REF] Sename | The inove anr 2010 blan 0308 project: Integrated approach for observation and control of vehicle dynamics[END_REF].

  This work was supported by the French national project INOVE / ANR 2010 BLAN 0308.

Fig. 1 .

 1 Fig. 1. Quarter-car vehicle model In this model, the quarter vehicle body is represented by the sprung mass (m s ), the wheel and tire are represented by the unsprung mass (m us ). They are connected by a spring with the stiffness coefficient k s and a semi-active damper. The tire is modeled by a spring with the constant stiffness coefficient k t . As seen in the figure, z s (respectively z us ) is the vertical displacement around the equilibrium point of m s (respectively m us ) and z r stands for the variation of the road profile. It is assumed that the wheel-road contact is ensured. The dynamical equations of a quarter vehicle model are given by: m s zs = -F spring -F sa m us zus = F spring + F sa -F tire (1) where F spring = k s (z s -z us ) is the dynamical spring force, F tire = k t (z usz r ) is the dynamical tire force. Let us denote z de f = z sz us the damper deflection and żde f = żsżus the damper deflection speed. The semi-active damper force is expressed by a linear model as follows: F sa = c.ż de f = c 0 żde f + u (2) with c min ≤ c ≤ c max and c 0 = (c min +c max) 2

  ξ = [χ quarter χ w f ] T χ quarter , χ w f are the quarter model and weighting function states respectively. z = [z 1 z 2 z 3 ] T are the controlled outputs, w = [z r ]: the disturbance input signal, y = [z de f żde f ]: output signals. u = u H ∞ : the suspension control signal derived from the H ∞ /LPV framework, ρ : the varying parameters, ρ ∈ [0.01 1].

  Fig. 5. Estimation of oil leakage

  Fig. 6. Sprung mass motion

  Fig. 7. Unsprung mass motion

Table 1 .

 1 validate the proposed LPV-FTC, simulations are performed on a non linear Renault Mégane Coupé (RMC) model from the test car available in MIPS Laboratory (Mulhouse, France) (see in Zin (2005)). The model parameters are given in the following table: Parameter values of the Renault Mégane Coupé quarter car model

	Parameter m s [kg]	m us [kg]	k s [N/m]	k t [N/m]	c[Nm/s]
	Value	315	37.5	29500	210000	{700, 5000}


M Q Nguyen

O Sename

L Dugard
email: luc.dugard@gipsa-lab.grenoble-inp.fr

A switched LPV observer for actuator fault estimation

Keywords: Multiplicative fault, fault estimation, switched extended observer, LPV system

In this paper, an actuator fault estimation problem is tackled, based on the LPV approach. First, the actuator fault is modeled in the form of a multiplicative fault by using a constant coefficient representing the loss of efficiency of the actuator's power α. Therefore, the estimation of a time varying faulty actuator can be transformed into the estimation of a constant coefficient α. Then, the faulty system is rewritten in the form of a switched LPV system. The coefficient α and the system states are estimated using an extended switched observer. The stability and performance of the observer is ensured considering a switched time-dependent Lyapunov function. The observer gains are derived based on LMI solutions for polytopic switched systems. Some simulation results are presented that show the effectiveness of this approach.

INTRODUCTION

Fault diagnosis plays an important role in the control and supervision of the industrial systems in order to enhance the safety, reliability and to reduce the loss of productivity. Nowadays, it has been attracted more and more attention from the researchers. In fact, fault estimation is a key step in Fault Detection and Isolation (FDI). More precisely, a FDI strategy can be model-based or data-based and is used to detect, isolate as well as estimate the faults. During the last decades, FDI modules, based on the analytical redundancy for fault estimation (e.g in Zhang and Jiang (2008)), have received a lot of attention by many researchers. In the literature, there exist many different approaches to estimate a fault occurring either on the actuator or on the sensor. Let us mention the classical methods, based on the parity space theory (see in [START_REF] Gertler | Fault detection and isolation using parity relations[END_REF]) to generate the residues and approximate the fault or the bank of observers approach (see [START_REF] Varrier | Détection de situations critiques et commande robuste tolérante aux défauts pour l'automobile[END_REF]) as well as by sliding mode observers (Edwards et al. (2000)). Recently, a new approach (see in [START_REF] Shi | Active fault tolerant control of lpv descriptor systems based on extended state observers[END_REF]) considered the fault element as a state of the augmented system and designed an extended observer to estimate at the same time the state and the fault of the system. However, it is limited to constant faults ḟ (t) = 0. Then, Zhang et al. (2008) presented a method allowing to evaluate the timevarying fault by using a fast adaptive fault estimation (FAFE) methodology based on an adaptive observer. But therein, the authors solved the problem with a regular Linear Time Invariant (LTI) system without considering the disturbances. Next, [START_REF] Rodrigues | Actuator fault estimation based adaptive polytopic observer for a class of lpv descriptor systems[END_REF] proposed an adaptive polytopic observer for time-varying fault estimation in spite of the disturbances for a class of descriptor Linear Parameter Varying (LPV) systems. Most of these works considered the fault as an additive one or transformed a multiplicative fault into an additive one.

The main purpose of this paper is to propose a methodology to estimate multiplicative faults for the actuators. An actuator time-varying fault is considered in the form of actuator power loss. Then, a constant "fault" coefficient is used to model the power loss. In this way, a constant coefficient should be estimated rather than the time-varying fault signal. The paper contribution are twofold:

• First, the faulty system is modeled in the form of an LPV system by considering the control input as a scheduling parameter. Then, the LPV system is rewritten to a switched LPV system. • Second, an LPV extended switched observer is designed to estimate both the "faulty" coefficient α and the system states. Two case studies are considered: the system without and with input disturbance. Some conditions on the decay rate of observer are added in order to increase the convergence's speed of estimation.

The paper is organized as follows: the next section recalls some preliminaries on the stability of the switched systems. Section 3 presents the problem statement. Section 4 gives a full description for an actuator multiplicative fault estimation based on the switched observer. Section 5 shows some simulation results of the proposed method using an academic example. Finally, some conclusions are drawn in the section 6.

PRELIMINARIES

This section is devoted to recall some results on the stability analysis for the continuous time, switched system by using the multiple Lyapunov function. Let us consider the following switched linear system:

ẋ(t) = A σ(t) x(t), x(0) = x 0 , (1)
where σ(t) is the switching signal and A σ(t) ∈ {A 1 , ..., A M } , A i ∈ R n×n , i = 1...M. Obviously, this model imposes discontinuity in A σ(t) since this matrix jumps instantaneously from A i to A j for i j. Now, the stability of the switched system (1) is guaranteed if there exists a family of symmetric and positive Lyapunov matrices {P 1 , ...P M }, each one associated to the correspondent A i such that the Lyapunov function V σ(t) (x(t)) is non increasing for all t ≥ 0 [START_REF] Liberzon | Switching in systems and control[END_REF]). Indeed we need to ensure that the Lyapunov function is non increasing at the switching instants (V(x(t k ) ≤ V(x(t - k )) with t k is the switching instant). This condition is conservative and it is replaced by a weaker condition that the sequence V(x(t k )) is decreasing for t k , k = 0, ..., ∞ is switching instant [START_REF] Geromel | Stability and stabilization of continuous-time switched linear systems[END_REF]) i.e if t k and t k+1 are successive switching times such that σ(t k ) = i and σ(t k+1 ) = j then V j (x(t k+1 )) ≤ V i (x(t k )). To ensure the global stability under slow switching, the concept of "dwell time τ d " is considered imposing that the interval between 2 successive switching instants satisfies t h+1 -t h ≥ τ d , for all switching times. However, in the design step, it is not easy to handle these conditions in terms of convex optimization (and therefore as LMIs), specially for systems with uncertainties or polytopic systems. An interesting solution has been recently proposed by Allerhand et al. (2011), using a time-dependent parametrized Lyapunov function: V σ (t) = x(t) P σ (t)x(t), where P σ(t) (t) is the Lyapunov matrix and is given by:

P σ(t) (t) =              P i,k + (P i,k+1 -P i,k ) t -τ s,k T/K := Pi,k i f t ∈ [τ s,k , τ s,k+1 ) P i,K i f t ∈ [τ s,K , τ s+1,0 ) P i 0 ,K i f t ∈ [0, τ 1 ) (2)
where i = 1, ...N with N is number of subsystems, i 0 = σ(0). τ 1 , τ 2 , ... are the switching instants. T is the dwell time satisfying τ s+1 -τ s ≥ T , and τ s,k = τ s +k(T/K) for k = 0, ...K , τ s = τ s,0 , τ s,K = τ s +T ). P i,k are symmetric matrices of compatible dimensions, where K is an integer that may be chosen a priori. This Lyapunov function allows to solve the conditions at switching instants easier and to deal with uncertainties or LPV systems. Throughout this paper, this Lyapunov function will be used.

PROBLEM FORMULATION

Consider the following single input continuous time linear invariant system:

ẋ(t) = Ax(t) + B 2 u(t) + B 1 w(t) (3)
y(t) = Cx(t) where x(t) ∈ R n , u(t) ∈ R, w(t) ∈ R p and y(t) ∈ R m are the state, the control input, the input disturbance and the measured output vectors, respectively. Matrices A, B 2 , B 1 ,C are real matrices of appropriate dimensions. Assume that the pair (A,C) is observable. Now, assume that the system (3) is in the faulty actuator Fig. 1. System with actuator fault situation, e.g loss of actuator power (Fig. 1). Let us model this failure by a multiplicative fault. In fact, denoting ū is the output of faulty actuator, then:

ū(t) = αu(t)
(4) where α stands for the faulty coefficient i.e if α = 1, the actuator is healthy , and when α = 0.8, the actuator loses 20% of its efficiency. Then, the system (6) becomes:

ẋ(t) = Ax(t) + B 2 αu(t) + B 1 w(t)
(5)

y(t) = Cx(t)
As mentioned previously, the objective of this work is to estimate the actuator fault. To do this, the main problem consists in estimating the coefficient α, which will be tackled, based on an extended switched observer and presented in the sequel.

Remark : It can be seen that if the control input u(t) = 0, then the fault information α in (5) becomes unobservable. The solution in the synthesis step becomes therefore infeasible and the coefficient α cannot be estimated . Thus, in this work, it is assumed that in the synthesis step, the control input u(t) is kept always from zero. To account for the change the sign of u(t), the system will be rewritten as a switched LPV system in the next section.

It is worth noting that the corresponding additive fault f (t) = (1 -α)u(t) is a time variant signal and depends on the value of the control input u(t), which is more complex to handle and estimate.

ACTUATOR FAULT ESTIMATION BASED SWITCHED OBSERVER

The actuator fault estimation in this work is done by estimating a effciency coefficient α of the actuator, and is given step by step in the sequel. Let us consider 2 cases studies: system without input disturbance and system with input disturbance.

System without input disturbance

In this section, it is assumed that the input disturbance w(t) is null or vanishes. Then, from the system (5), one obtains a system with the fault actuator but without input disturbance as follows:

ẋ(t) = Ax(t) + B 2 αu(t) (6)
y(t) = Cx(t)
where α has to be estimated. Since the control input u(t) is known, in this work, one can model the system (6) as an LPV system by choosing u(t) as a time-varying parameter. Moreover, it is assumed that the control input u(t) is bounded, i.e: 0 <| u(t) |≤ u max , u max > 0 (7) Let us rewrite u(t) = |u(t)|sign(u(t)), and denote ρ(t) = |u(t)| as a time varying parameter. Then,

u(t) = ρ(t)sign(u(t)) =      ρ(t) i f u(t)) > 0 -ρ(t) i f u(t) < 0 (8)
Note that from (7), ρ(t) is also bounded: 0 < ρ(t) ≤ u max . Moreover, the bounds of ρ shoud be chosen to be able to apply the polytopic approach in the latter. Here, one chooses ρ(t) ∈ [ρ ρ] where ρ = 0.01, ρ = u max . Then (6) becomes: 9)

ẋ(t) = Ax(t) + B 2 αρ(t)sign(u(t)) = Ax(t) + B(ρ)αsign(u(t))(
y(t) = Cx(t) where B(ρ) = B 2 ρ(t).
As seen before, α is constant, then α = 0. Hence the LPV system (9) can be augmented and rewritten in the following form:

ẋ α = A B(ρ)σ 0 0 A e (ρ)
x α (10)

y = [C 0] C e x α
where σ = sign(u(t))

Furthermore, it can be rewritten in a switched LPV system form [START_REF] Lu | Switching lpv control designs using multiple parameter-dependent lyapunov functions[END_REF]) as follows:

ẋ α = A e,σ (ρ) x α (11) y = [C 0] C e x α
where A e,σ (ρ) is switched between two modes A e,1 (ρ), A e,2 (ρ) according to the switched signal σ(t) = sign(u(t)) as follows:

A e,σ (ρ) =                                A B(ρ) 0 0 A e,1 (ρ) i f u(t) > 0 A -B(ρ) 0 0 A e,2 (ρ) i f u(t) < 0 (12)
It is worth noting that the previous definition of switched system allows to handle control input signals that crosses the 0 value, which cannot be done when choosing u(t) as scheduling parameter. Now, the following LPV extended switched observer is used to estimate the extended state of the switched system (11):

ẋ α = A e,σ (ρ) x α + K σ (t)(y -ŷ) (13) ŷ = C e x α
From ( 11) and ( 13), the estimation error e(t) satisfies:

ė = ėx ėα = (A e,σ (ρ) -K σ (t)C e ) e x e α (14)
The first main result about stability of the estimation error is given below.

Theorem 1. Consider the switched LPV system (11) and the extended switched observer (13), if there exists a collection of matrices P i,k > 0, Y i,k , k = 0, ...K, i = 1, ...M (M = 2: number of subsytems), of appropriate dimensions and K is prescribed integer, such that for all i = 1, ..M and j = 1, ...N, (N = 2: number of the vertices of the polytope), the following LMIs hold:

(P i,k+1 -P i,k ) T/K + A ( j) e,σ P i,h -C e Y i,h + P i,h A ( j) e,σ -Y i,h C e < 0 (15) for k = 0, ...K -1, h = k, k + 1, A ( j) e,σ P i,K -C e Y i,K + P i,K A ( j) e,σ + P i,K A ( j) e,σ -Y i,K C e < 0 ( 16
) P i,K -P l,0 ≥ 0 ∀l = 1, ..i -1, i + 1, ...M. ( 17
)
then

K σ(t) (t) = P σ(t) (t) -1 Y σ(t) (t) =            P-1 i,k Ŷi,k i f t ∈ [τ s,k , τ s,k+1 ) P -1 i,K Y i,K i f t ∈ [τ s,K , τ s+1,0 ) P -1 i 0 ,K Y i 0 ,K i f t ∈ [0, τ 1 ) (18)
is the gain of the extended observer (13) and the error estimation asymptotically converges to zero for a given dwell time T . where

Y σ(t) (t) =              Y i,k + (Y i,k+1 -Y i,k ) t -τ s,k T/K := Ŷi,k i f t ∈ [τ s,k , τ s,k+1 ) Y i,K i f t ∈ [τ s,K , τ s+1,0 ) Y i 0 ,K i f t ∈ [0, τ 1 ) (19)
Proof:

Define the Lyapunov function candidate: V σ (t) = e(t) P σ (t)e(t) where P σ (t

) = P σ (t) > 0. Denoting Y σ(t) (t) = P σ(t) (t)K σ(t) (t).
Then the estimation error in ( 14) is asymptotically stable if: Vσ < 0 ⇔ ė P σ e + e P σ ė + e Ṗσ e < 0 ⇔

(A e,σ (ρ) -K σ (t)C e ) P σ (t) + P σ (t)(A e,σ (ρ) -K σ (t)C e ) + Ṗσ (t) < 0 (20) ⇔ A e,σ (ρ) P σ (t) -C e Y σ (t) + P σ (t)A e,σ (ρ) -Y σ (t)C e + Ṗσ (t) < 0 (21)
From the formula of P σ (t) in (2), we ensure (21)⇔

A e,σ (ρ

) P i,h -C e Y i,h + P i,h A e,σ (ρ) -Y i,h C e + (P i,k+1 -P i,k ) T/K < 0 (22) holds for h = k, k + 1, i = 1, 2, k = 0, ...K -1.
and A e,σ (ρ) P i,K -C e Y i,K + P i,K A e,σ (ρ)

+ P i,K A e,σ (ρ) -Y i,K C e < 0 ( 23
)
hold for i = 1, 2.

The equation ( 22) guarantees that the Lyapunov function V σ (t) decreases during the time interval t ∈ [τ s,0 , τ s,K ). The LMIs (23) ensures that V σ (t) decreases after the dwell time and before the next switching time, i.e t ∈ [τ s,K , τ s+1,0 ).

From the definition of P σ(t) (t), consider that at the switching instant τ k , the system switches from the mode i to the mode l, and

V(τ - k ) = x(τ k ) P i,K x(τ k ), V(τ k ) = x(τ k ) P l,0 x(τ k ).
To guarantee the non-increasing of the Lyapunov function at the switching instants, we must ensure that V(τ k ) ≤ V(τ - k ) which corresponds to the inequality (17), i.e: P i,K -P l,0 ≥ 0 ∀l = 1, ..i -1, i + 1, ...M. In our case, we have 2 subsystems, M = 2, i.e P 1,1 -P 2,0 ≥ 0 and P 2,1 -P 1,0 ≥ 0. Now, in order to resolve the LMIs in ( 22), ( 23), that are parameter dependent, the polytopic approach for LPV systems is considred where the 2-vertices polytope is given by Ω ρ = ρ ρ . The LMIs ( 22), ( 23) is therefore solved at the vertices of the polytope (i.e ρ, ρ). Then, one obtains the LMIs in ( 15), ( 16). Finally, K σ (t) = P -1 σ Y σ (t) is the gain of the switched extended observer (13).

The following theorem extends the previous one imposing a pole placement constraint for the observer design, imposing a prefixed decay rate for the exponential stability of the estimation error.

Theorem 2. Consider the switched LPV system (11) and the switched extended observer (13), if there exists a collection of matrices P i,k > 0, Y i,k , k = 0, ...K, i = 1, ...M (M = 2: number of subsytems), of appropriate dimensions, K is prescribed integer, and a positive scalar β such that for all i = 1, ..M and j = 1, ...N, (N = 2: number of the vertices of the polytope), the following LMIs hold:

(P i,k+1 -P i,k ) T/K + A ( j) e,σ P i,h -C e Y i,h + P i,h A ( j) e,σ -Y i,h C e + 2βP i,h < 0 (24) for k = 0, ...K -1, h = k, k + 1, A ( j) e,σ P i,K -C e Y i,K + P i,K A ( j) e,σ + P i,K A ( j) e,σ -Y i,K C e + 2βP i,K < 0 (25) P i,K -P l,0 ≥ 0 ∀l = 1, ..i -1, i + 1, ...M. ( 26
) then K σ(t) (t) = P σ(t) (t) -1 Y σ(t) (t)
is the gain of the extended observer (13) and this observer converges with a decay rate of β, i.e the estimation error e(t) asymptotically exponentially converges to zero for a dwell time of T .

Proof: The idea of the proof is similar to the "Theorem 1" and is omitted here. We need only to note that in order to increase the convergence's speed of the estimation, one introduces another term of 'decay rate', i.e, if there exists a positive scalar β such that Vσ (t) ≤ -2βV σ (t), then the estimation error converges to zero with a decay rate of β.

System with input disturbance

In this section, the fault estimation in the presence of input disturbance is considered, i.e, the system ( 5) is taken into account:

ẋ(t) = Ax(t) + B 2 αu(t) + B 1 w(t) (27) y(t) = Cx(t)
Similarly, one has u(t) = |u(t)|sign(u(t)), and denote ρ(t) = |u(t)| as a time varying parameter and ρ(t) is also bounded: ρ ≤ ρ(t) ≤ ρ, where ρ = 0.01, ρ = u max . Then, (27) becomes:

ẋ(t) = Ax(t) + B(ρ)αsign(u(t)) + B 1 w(t) (28)
y(t) = Cx(t)
The system ( 28) is augemented as follows:

ẋ α = A e,σ (ρ) x α + B 1 0 w (29) y = [C 0] C e x α
where A e,σ (ρ) is switched between two modes A e,1 (ρ), A e,2 (ρ) according to the signal σ(t) = sign(u(t)), as in ( 12). Then, the following switched extended observer is proposed to estimate the system's state and the coefficient α :

ẋ α = A e,σ (ρ) x α + K σ (t)(y -ŷ) (30) ŷ = C e x α
From ( 29) and ( 30), the estimation error e(t) satisfies:

ė = ėx ėα = (A e,σ (ρ) -K σ (t)C e ) e x e α + B 1e w (31)
where B 1e = B 1 0 and K σ (t) is the observer gain which has to be determined. The extended observer (30) and the estimation error (31) are affected by the disturbance effect w(t). The observer design problem in this case involves the calculation of the observer gain K σ (t) such that:

• when w(t) ≡ 0, the observer is asymptotically stable • when w(t) 0, the following L 2 -induced gain performance criterion is satisfied :

min γ s.t e 2 < γ w 2 (32)
where . 2 stands for L 2 norm. The following theorem solves the above problem.

Theorem 3. Consider the switched system (29) and the switched extended observer (30), if there exists a collection of matrices P i,k > 0, Y i,k , k = 0, ...K, i = 1, ...M (M = 2: number of subsytems), of appropriate dimensions and K is prescribed integer, such that for all i = 1, ..M and j = 1, ...N, (N = 2: number of the vertices of the polytope), the following LMIs hold:

              (P i,k+1 -P i,k ) T/K + A ( j) e,σ P i,h -C e Y i,h + P i,h A ( j) e,σ -Y i,h C e P i,h B 1e I * -γ 2 I 0 * * -I               < 0 (33) for k = 0, ...K -1, h = k, k + 1,            A ( j) e,σ P i,K -C e Y i,K + P i,K A ( j) e,σ + P i,K A ( j) e,σ -Y i,K C e P i,K B 1e I * -γ 2 I 0 * * -I            < 0 (34) P i,K -P l,0 ≥ 0 ∀l = 1, ..i -1, i + 1, ...M. (35) then K σ(t) (t) = P σ(t) (t) -1 Y σ(t) (t) =            P-1 i,k Ŷi,k i f t ∈ [τ s,k , τ s,k+1 ) P -1 i,K Y i,K i f t ∈ [τ s,K , τ s+1,0 ) P -1 i 0 ,K Y i 0 ,K i f t ∈ [0, τ 1 ) (36)
is the gain of the extended observer (30) and the error estimation asymptotically converges to zero for a dwell time of T . where

Y σ(t) (t) =              Y i,k + (Y i,k+1 -Y i,k ) t -τ s,k T/K := Ŷi,k i f t ∈ [τ s,k , τ s,k+1 ) Y i,K i f t ∈ [τ s,K , τ s+1,0 ) Y i 0 ,K i f t ∈ [0, τ 1 ) (37)
Proof: From the Bounded Real Lemma, the condition ( 32) is satisfied if the following condition holds:

V(t) + e e -γ 2 w w < 0 (38) ⇔           (A e,σ (ρ) -K σ (t)C e ) P σ (t) + P σ (t)(A e,σ (ρ) -K σ (t)C e ) + Ṗσ (t) P σ (t)B 1e I * -γ 2 I 0 * * -I           < 0 (39)
From the formula of P σ (t) in ( 2), ( 39) is satisfied if

              A e,σ (ρ) P i,h -C e Y i,h + P i,h A e,σ (ρ) -Y i,h C e + (P i,k+1 -P i,k ) T/K P i,h B 1e I * -γ 2 I 0 * * -I               < 0 (40) holds for h = k, k + 1, i = 1, 2, k = 0, ...K -1. and           A e,σ (ρ) P i,K -C e Y i,K + P i,K A e,σ (ρ) + P i,K A e,σ (ρ) -Y i,K C e P i,K B 1e I * -γ 2 I 0 * * -I           < 0 (41) hold for i = 1, 2.
The equation ( 40) guarantees that the Lyapunov function V σ (t) decreases and (38) holds during the time intervals i.e t ∈ [τ s,0 , τ s,K ). The LMIs (41) ensure that V σ (t) decreases and that (38) holds after the dwell time and before the next switching instant, i.e t ∈ [τ s,K , τ s+1,0 ). From the definition of P σ (t), consider that at instant τ k , the system switches from the mode i to the mode l, to guarantee the non-increasing of the Lyapunov function at the switching instants, we must ensure:

P i,K -P l,0 ≥ 0 ∀l = 1, ..i -1, i + 1, ...M. (42
) In our case, we have 2 subsystems, M = 2, i.e P 1,1 -P 2,0 ≥ 0 and P 2,1 -P 1,0 ≥ 0. Now, in order to resolve the LMIs in (40), (41), we apply the polytopic solution for LPV system where the polytope is given by Ω ρ = ρ ρ and obtain the LMIs in (33), (34). The next result extends the previosu one imposing a prefixed decay rate on the convergence of the estimation error.

Theorem 4. Consider the switched system (29) and the switched observer (30), if there exists a collection of matrices P i,k > 0, Y i,k , k = 0, ...K, i = 1, ...M (M = 2: number of subsytems), of appropriate dimensions, K is prescribed integer, and a positive scalar β such that for all i = 1, ..M and j = 1, ...N, (N = 2: number of the vertices of the polytope), the following LMIs hold:

              (P i,k+1 -P i,k ) T/K + A ( j) e,σ P i,h -C e Y i,h + P i,h A ( j) e,σ -Y i,h C e + 2βP i,h P i,h B 1e I * -γ 2 I 0 * * -I               < 0 (43) for k = 0, ...K -1, h = k, k + 1,            A ( j) e,σ P i,K -C e Y i,K + P i,K A ( j) e,σ + P i,K A ( j) e,σ -Y i,K C e + 2βP i,K P i,K B 1e I * -γ 2 I 0 * * -I            < 0 (44) P i,K -P l,0 ≥ 0 ∀l = 1, ..i -1, i + 1, ...M. ( 45
) then K σ(t) (t) = P σ(t) (t) -1 Y σ(t) (t)
is the gain of the extended observer (30) and the error estimation asymptotically converges to zero for a dwell time of T .

Proof: The proof is similar to the last cases and is omitted here for the simplification. Finally, in order to design the switched observer (30), one has to solve the following optimization problem:

min P i,k ,Y i,k γ 2
subject to ( 43), ( 44), ( 45) and P i,k > 0 (46) By solving this optimization problem, one can derive P σ (t), Y σ (t) and the observer gain is calculated by K σ (t) = P σ (t) -1 Y σ (t).

NUMERICAL EXEMPLE

The following numerical example illustrates the effectiveness of the proposed observer. Let consider a LTI system of the form (3) with the following matrices:

A = -3 -1 1 0 , B 2 = 1 0 C = 1 -1 0 1 , B 1 = 0.1 1
The system is subject to the control input u(t) = 1 + 2 sin(t) (see Fig. 2) and the initial value x(0) = 0 0 . Thus, one has :

-1 ≤ u(t) ≤ 3 and 0 < |u(t)| ≤ 3.

Then, the time varying parameter ρ(t) = |u(t)| is assumed to be bounded by: 0.01 ≤ ρ ≤ 3 and A e,σ (ρ) is switched between two modes A e,1 (ρ), A e,2 (ρ) according to the switched signal σ(t) = sign(u(t)), and is given by ( 12), and at each mode, one has 2 matrices corresponding to 2 vertices of the polytope Ω ρ = ρ ρ as follows:

                             A e,1 (ρ) =                  -3 -2 0.01 1 0 0 0 0 0         ,         -3 -2 3 1 0 0 0 0 0                  A e,2 (ρ) =                  -3 -2 -0.01 1 0 0 0 0 0         ,         -3 -2 -3 1 0 0 0 0 0                 
The figure Fig. 2 shows also the scheduling parmeter and the switched signal. Finamly, as seen in figure 3 the efficiency coefficient of the actuator is such that, for t ∈ [0 -5]s, α = 1, the actuator is healthy, then, for t > 5s, α = 0.5, the actuator loses 50 % of its efficiency. Now, the different results are given for the two case studies (without and with input disturbance) with the observer design of section 3. In this section, the fault estimation for the system without input disturbance is considered. By solving the LMIs (24,25, 26) (where the decay rate is taken into account) for the 2 modes A e,1 (ρ), A e,2 (ρ) to obtain the matrices P i,k , Y i,k , then the observer gain is calculated by K σ (t) = P σ (t) -1 Y σ (t). Fig. 3 and Fig. 4 show the estimation of the faulty coefficient α and the system states, with and without taking into account the decay rate in the observer design step. It shows that the switched observer allows to estimate the faulty coefficient as well as the system state. Moreover adding a decay rate constraint allows to get a faster response of the observer. The disturbance effects on the estimation are considered in this section. See the disturbance appearing at t = 10s in the Fig. 5. By solving the LMIs (43,44, 45) (where the decay rate is taken into account) for the 2 modes A e,1 (ρ), A e,2 (ρ) to obtain the matrices P i,k , Y i,k , then the observer gain is calculated by K σ (t) = P σ (t) -1 Y σ (t). Fig. 6, Fig. 7, show the estimation of the effciency coefficient α and the system states without and with taking into account the decay rate. Obviously, despite of the input disturbance w(t), the switched observer allows to have a good estimation of the coefficient α and the system states, in particular for the observer with decay rate which gives a much better result.

CONCLUSION

In this paper, an actuator fault estimation was proposed within the LPV approach. The actuator fault is modeled in a multiplicative way by using a constant coefficient (α ∈ [0 1]) that represents the fault information. This allows to facilitate the fault estimation. The estimation is based on an extended switched observer and the simulation results show the effectiveness of the proposed approach. In the future work, this approach can be applied on a real system such as the suspension system. Therein, we want to estimate the fault on the damper malfunction (e.g. oil leakage of the damper,...) where some preliminary results are presented in [START_REF] Nguyen | An LPV fault tolerant control for semi-active suspension -scheduled by fault estimation[END_REF]. Moreover, some fault tolerant control strategies can be developped, based on the proposed fault diagnosis. Allerhand, L., Shaked, U., et al. (2011). Robust stability and stabilization of linear switched systems with dwell time.

Automatic Control, IEEE Transactions on, 56(2), 381-386. Edwards, C., Spurgeon, S.K., and Patton, R.J. (2000). Sliding mode observers for fault detection and isolation. Automatica, 36(4), 541-553.

Fig. 2 .

 2 Fig. 2. The control input, varying parameter and switched signal 5.1 System without input disturbance (case 1)

Fig. 3 .

 3 Fig. 3. Estimation of α without decay rate β = 0 (left) and with β = 1 (right)-(case 1)

Fig. 4 .Fig. 6 .

 46 Fig. 4. State estimation without decay rate β = 0 (left) and with β = 1 (right)-(case 1)

Fig. 7 .

 7 Fig. 7. State estimation without decay rate β = 0 (left) and with β = 1 (right)-(case 2)


Cyril Jean

Laurent Belliard

Loïc Becerra

Bernard Perrin

Backward Propagating Acoustic Waves in Single Gold Nanobeams

Femtosecond pump-probe spectroscopy has been carried out on suspended gold nanostructures with a rectangular cross section lithographed on a silicon substrate. With a thickness fixed to 110 nm and a width ranging from 200 nm to 800 nm, size dependent measurements are used to distinguish which confined acoustic modes are detected. Furthermore, in order to avoid any ambiguity due to the measurement uncertainties on both the frequency and size, pump and probe beams are also spatially shifted to detect guided acoustic phonons. This leads us to the observation of backward propagating acoustic phonons in the gigahertz range (∼3 GHz) in such nanostructures. While backward wave propagation in elastic waveguides has been predicted and already observed at the macroscale, very few studies have been done at the nanoscale. Here, we show that these backward waves can be used as the unique signature of the width dilatational acoustic mode.

Probing the elasticity at the nanoscale is a challenge that led a wide community to study confined acoustic modes of nano-objects in the 10 GHz-1 THz range using time-resolved pump-probe experiments [1]. In such an approach, acoustic modes are excited by the absorption of a femtosecond laser pulse and detected in transmission [2,3] or reflectivity geometry in near [4,5] or far [6,7] field by the induced change in the material optical properties. As ensemble studies result in inhomogeneous broadening of the acoustic features [8,9], single particle spectroscopy has been considered to circumvent this drawback and to clarify the acoustic response [3,10,11]. However, the strong damping at these extremely high frequencies has led several groups to isolate the nanoresonators from their substrate to avoid energy leaking through the nanoobject-substrate contact [12][13][14]. These breakthroughs have made possible the observation of an other source of acoustic energy leaks. Indeed, it has been proven that acoustic phonons are also guided along the nanoobjects [15,16]. Other major developments illustrate the possibility to use nano-objects as nanoscale acoustic transducers for hypersonic wave imaging [17]. Non destructive imaging with nanometric resolution in both depth and lateral direction is now one step ahead.

Furthermore, there is a recent and growing interest in both the electromagnetic and acoustic wave communities for materials and waveguides that exhibit backward propagating waves [18]. In backward waves, the phase velocity describing the propagation of individual wave fronts in a wave packet and the energy flux of the wave, characterized by the Poynting vector are anti-parallel. This anti-parallel propagation constitutes the definition of a negative index material and opens the way to a large variety of intriguing physical phenomena [19]. There now exist multiple experimental evidences of this effect, and applications of negative index materials for electromagnetic waves have already been conceived [20]. The pos- * Laurent.Belliard@upmc.fr sibility of such backward wave motion in elastic waveguides has been revealed by the early work of Lamb [21]. Many different elastic waveguides such as cylinders [22] or surface acoustic waves [23] can also exhibit backward wave motion. Recently, negative refraction and focusing of elastic Lamb waves have been investigated but are still confined to the low-MHz range [24,25]. In this letter we provide evidence of backward wave propagation in one-dimensional elastic waveguide in the gigahertz range. We use time-resolved pump-probe spectroscopy with spatially shifted pump and probe beams to excite and detect the backward wave motion. As a first application, unambiguous mode identification using backward wave propagation is investigated.

In order to study guided modes in one dimensional elastic waveguides, gold nanostructures with a rectangular cross section are conceived (Fig. 1). Their thickness, h, is fixed to h = 110 nm and their width, W , ranges from 200 to 800 nm. As shown in Fig. 1(a), the gold stripes are connected to large mesas at both ends. To reduce the influence of the silicon substrate on the elastic confinement, the gold nanostructures are suspended several micrometers above the substrate as illustrated in Fig. 1(b). This geometry results in what we call "gold nanobeams". All boundaries are thus mechanically free. As their length is greater than 15 µm, they are considered as infinitely long and the fixed boundary conditions at the extremities can be neglected. The samples are fabricated by electronic lithography coupled with anisotropic wet etching on silicon substrates. To prepare these structures, single face polished Si(001) substrates are used. First, an electronic resist (PMMA-950K-A6) is spin coated on the Si wafer surface. The exposition is made at 20kV using a Zeiss Supra 40 scanning electron microscope equipped with a Raith Elphy Quantum module. After being developed, rinsed with water and blown dry, the sample is coated with 10 nm of chromium to improve the adhesion, followed by a 110 nm thick gold layer. Then, a lift-off process in aceton with careful ultrasonic agitation is performed to delimit the future gold nanobeams. Finally, the sample is dipped in a hot 40% KOH solution to partially etch the Si substrate and free the desired gold nanobeams. This quick anisotropic Si etching process leads to the formation of suspended gold nanobeams between two big gold pads. The key point of this proces is that the main axis of the beams has been rotated 45 • with respect to the < 110 > direction of Si substrates to allow this underetching.

Our experimental setup working in reflection geometry was described in detail elsewhere [7]. Ultrafast pumpprobe spectroscopy experiments are performed using a mode-locked Ti:sapphire (MAI-TAI Spectra) laser source operating at 800 nm with a pulse duration below 100 fs at the laser output and a repetition rate of 78.8 MHz. Synchronous detection on the sample reflectivity is performed by modulating the pump beam at 1.8 MHz. Both pump and probe beams are focused by means of a microscope objective with a N A = 0.9 and are normally incident on the sample. The laser spots can be focused around 1 µm diameter at 1/e 2 . A telescope is fixed on a XY piezoelectric stage such that the probe beam is laterally positioned with respect to the fixed pump beam. A two-color experiment is performed by doubling the pump frequency (λ = 400 nm) with a nonlinear cystal (BBO) to avoid scattered light coming from the pump beam. A dichroic filter located in front of the diode system suppresses the light of the pump beam, its power is reduced around 500 µW and the power of the probe beam does not exceed 15 µW.

Such experimental conditions place us in the thermoelastic regime and the acoustic signal and optical re- flectivity remain stable during all the average processing. The reflectivity from the sample is measured by an avalanche photodiode and analyzed with a lock-in amplifier. A maximum pump-probe time delay equal to 12 ns is achieved using a mobile reflector system mounted on a translation stage. First, let us consider the case where the pump and probe beams are spatially superimposed. Fig. 2 shows the time domain reflectivity change measurements ∆R/R obtained on a 400 nm ± 20 nm wide (determined by SEM measurement) and 110 nm thick gold nanobeam in a usual reflectometry set-up. After a sharp and sudden rise, a slowly decreasing non-oscillatory background in ∆R/R is observed. This is the signature of the rapid heating of the electron gas induced by the pump pulse absorption, followed by the slow cooling down process. This thermal stress launches the acoustic vibrations of the nanobeam [26]. The time-resolved signature also shows a superposition of high (inset Fig. 2(a)) and low frequency oscillations. By performing a numerical fast Fourier transform (Fig. 2(b)) or by fitting the time oscillation with a damped sinus function, one can extract a 2.8 GHz ± 0.5 GHz and a 16 GHz ± 0.5 GHz acoustic signatures. In the following we consider that the Young modulus, the Poisson's ratio and the density of gold are respectively 79 GPa, 0.44 and 19 300 kg /m 3 [START_REF] Meyers | Mechanical Behavior of Materials[END_REF]. The longitudinal v L and transverse v T sound velocities are then 3.6 µm ns -1 and 1.2 µm ns -1 respectively. Thus, the higher 16 GHz frequency appears to be the thickness vibration of the beam owing to the fact that v L /(2h) = 16.5 GHz. The displacement field of this d 2 eigenmode is plotted in Fig. 3 that this thickness resonance signature is also observed when the pump-probe experiment is undertaken on the gold mesa using an interferometric scheme in this case as FIG. 4. Using a finite elements method, the acoustic mode frequency at wave number k = 0 is plotted versus the nanobeam width. Flexural, shear and dilatational modes are respectively the orange dotted lines, the green dashed lines and the blue solid lines. Experimental measurements are the full black circle with bars. The vertical black dashed line at 370 nm is the specific width where dispersion relations have been calculated (Fig. 3) and propagation has been observed (Fig. 5).

(b). At this point, we verify

photoelastic signals are negligible in gold layers [START_REF] Perrin | [END_REF]. Consequently, the reflectometry measurement undertaken on the gold nanobeam is probably partly an interferometric like detection with the light reflected by the silicon substrate underneath acting as a reference mirror. The case of the lower frequency ∼3 GHz is dealt with by performing an eigenmode analysis by finite elements method. The dispersion relation of a 370 nm wide nanobeam is plotted in Fig. 3. Several acoustic modes with non-zero cut-off frequencies are revealed. Three distinct families emerge : the flexural modes f 1 , f 2 and f 3 , the dilatational modes d 1 and d 2 and the shear modes s 1 , s 2 and s 3 (see Fig. 3(b)). The frequency evolution of s 1 , s 2 , f 1 , f 2 , f 3 and d 1 with the nanobeam width is plotted in Fig. 4. The measured frequencies on different gold nanobeams are also reported. Owing to the measurement uncertainties, the experimental frequencies ranging from 1.7 GHz to 5.5 GHz may be consistent with the frequencies of the modes f 2 , s 2 and d 1 which exhibit close cut-off frequencies. Consequently, additional evidence on the detected mode is needed. Given that the excitation process imposes the relative amplitudes of the excited modes, one can use the initial displacement field projection on the orthogonal basis formed by the nano-object eigenmodes to identify which mode will be efficiently excited [29,30].

Here we propose to discriminate unambiguously between the different acoustic modes by probing their propaga- tion. Very different behaviors are expected owing to the very different slopes of the dispersion relationship.

In order to be sensitive to the propagation phenomenon, the pump and probe beams have to be separated [8]. The transient reflectivity measured on a 370 nm width gold nanobeam, with pump-probe separation ranging from 0 µm to 2.1 µm, is presented in Fig. 5(a). The vertical scale corresponds to the pump-probe spatial separation. The exponential thermal background has been removed. As physical insight on the attenuation process is beyond the scope of this study, each transient reflectivity signal is normalized by its maximum : the eighth experimental signal with 2.1 µm spatial separa-tion is one order of magnitude less intense than the first one when pump and probe beams are spatially superimposed. The signal is mainly composed by a 3.2 GHz oscillation. As the pump-probe separation increases, this wave-packet shifts to a longer time delays as expected for a propagation guided along the nanobeam. This is analogous to what has already been observed in copper nanowires [15]. However, let us have a closer look at the individual wave fronts (Fig. 5(e)). It appears that as the wave-packet propagates (from darker to lighter curve), the individual wave fronts move to smaller time-delay. To investigate this propagation quantitatively, we use the previously developed expression of the reflectivity for a propagating gaussian wave-packet with the parabolic dispersion relation ω(k) = αk 2 + βk + γ [15]

∆r(t) ∝ Re exp(-iγt) √ σ 2 + 32iαt exp - 8(z 0 -βt) 2 σ 2 + 32iαt (1)
where z 0 is the spatial pump-probe separation and σ is the root of the sum of squares of the pump and probe beam diameters at 1/e 2 . A polynomial fit at small wavenumber for d 1 , f 2 and s 2 (see inset Fig. 3(a)) is used. We then evaluate the theoretical reflectivity signal for these three modes d 1 (Fig. 5(b)), f 2 (Fig. 5(c)) and s 2 (Fig. 5(d)) at different pump and probe separations. The progressive shift of the wave-packet is very slow for the f 2 mode (Fig. 5(c)), due to the very flat dispersion relation. Both the f 2 and s 2 modes show parallel propagation of the individual wave fronts and the wave packet (Fig. 5(g),(h)), whereas the d 1 mode exhibits the expected backward propagation behavior (Fig. 5(f)). According to our simplified analytical expression, the parallel (forward) or anti-parallel (backward) behavior is strongly dependent on the sign of the curvature parameter α. The first dilatational mode d 1 is the only one to exhibit a negative curvature α at low wave numbers k which are the only wave numbers excited owing to the laser spot diameters of around 1 µm. To confront the experiment, the envelope of the theoretical signal d 1 is also plotted on the experimental oscillation and shows a good agreement (Fig. 5(a)). Finally, it could be argued that the experimental signal is a superposition of two or more acoustic modes. However, if a contribution of beyond 10% of f 2 or s 2 is added in the simulated signal d 1 , the resulting simulated propagative wave no longer exhibits the backward wave behavior. The backward wave observation in our 1D elastic waveguide is thus the unique fingerprint of the first width dilatational mode d 1 .

In conclusion, two dilatational acoustic modes are investigated on e-beam lithographed gold nanobeams by pump-probe time-resolved spectroscopy leading to the observation of backward wave propagation in the gigahertz range. As a first application, we show unambiguous acoustic mode discrimination using this unique property in the dispersion relations. Furthermore, as already observed in 2D Lamb waveguides, one can imagine to tune the width and thickness of our nanobeam, thus opening the way to negative refraction physics in acoustics at the nanoscale.

FIG. 1 .

 1 FIG. 1. (a) Scanning electron microscope image (top view) of several gold nanobeams suspended a few micrometers above the silicon substrate. (b) Artist view of nanobeams with a rectangular cross section suspended above the silicon substrate. (c) Scanning Electron microscope image (top view) of a 400 nm wide gold nanobeam. The edge roughness is lower than 20 nm.

FIG. 2 .

 2 FIG. 2. (a) Typical time-resolved modulation of the reflected probe beam ∆R/R of a 400 nm ± 20 nm wide and 110 nm thick gold nanobeam that exhibits a low frequency oscillation at 2.8 GHz superimposed to a higher modulation around 16 GHz (the inset is a zoom where the lower modulation is removed, dotted line is a 16 GHz sinus). (b) Fourier Power spectra of the time transient signal (a) showing both 2.8 GHz and 16 GHz acoustic signatures. d1 and d2 refer to the first dilatational modes.

2 FIG. 3 .

 23 FIG. 3. (a) Acoustic modes dispersion relations in a 370 nm wide and 110 nm thick gold nanobeam (with a 10 nm thick adhesion layer of chromium). fi, si and di refer respectively to flexural, shear and dilatational modes. The inset shows the parabolic fit of three specific modes at low wave number used to evaluate the corresponding propagating wave packet. The parabolic fit is respectively f d 1 (k) = -3.4 × 10 -2 k 2 -9.2 × 10 -5 k + 3.2 GHz, fs 2 (k) = 5.0 × 10 -2 k 2 + 5.0 × 10 -4 k + 3.6 GHz and f f 2 (k) = 5.7 × 10 -3 k 2 + 2.8 × 10 -4 k + 3.4 GHz for the modes d1, s2 and f2 with the wave number k in µm -1 . (b) Cross-sectional representation of the displacement field of three flexural (f1,f2 and f3), two dilatational (d1 and d2) and three shear (s1,s2 and s3) acoustic modes. The black arrows show the displacement direction.

2 FIG. 5 .

 25 FIG. 5. (a) Experimental transient relative reflectivity signal at different spatial pump-probe separations (the decreasing exponential background has been removed). The offset corresponds to the spatial pump-probe shift and the reflectivity modulation is normalised by the maximum. The envelope of the eight experimental traces in (a) are the simulated envelope of the first dilatational mode d1 after propagation . (b), (c), (d) are the simulated transient relative reflectivity signal of the modes d1, f2 and s2 according to the parabolic fit of the dispersion relation (Fig.3(a)). (e), (f), (g), (h) are respectively the same signals as in (a), (b), (c), (d) zoomed in a shorter time window. Only half of the eight curves are plotted for readability. The offset is also suppressed in order to observe the propagation direction of the phase of the wave. The black arrow indicates this propagation direction.

ACKNOWLEDGMENTS

The authors would like to thank Michael Rosticher at the ENS for the EBL steps during the maintenance in our facilities, Mélanie Escudier for the gold deposition process and Frances Edwards for her careful proofreading of the manuscript's English.