A Batailly

Q Agrapart

A Millecamps

J.-F Brunel

Experimental and numerical simulation of a rotor/stator interaction event localized on a single blade within an industrial high-pressure compressor

Keywords: rotor/stator interaction, rubbing, nonlinear dynamics, unilateral contact dynamics, friction, abradable coating, high-pressure compressor

This contribution addresses a confrontation between the experimental simulation of a rotor/stator interaction case initiated by structural contacts with numerical predictions made with an in-house numerical strategy. Contrary to previous studies carried out within the low-pressure compressor of an aircraft engine, this interaction is found to be non-divergent: high amplitudes of vibration are experimentally observed and numerically predicted over a short period of time. An in-depth analysis of experimental data first allows for a precise characterization of the interaction as a rubbing event involving the first torsional mode of a single blade. Numerical results are in good agreement with experimental observations: the critical angular speed, the wear patterns on the casing as well as the blade dynamics are accurately predicted. Through out the article, the in-house numerical strategy is also confronted to another numerical strategy that may be found in the literature for the simulation of rubbing events: key differences are underlined with respect to the prediction of non-linear interaction phenomena.

1 -Department of Mechanical Engineering, École Polytechnique de Montr éal, P.O. Box 6079, Succ. Centre-Ville, Montr éal, Qu ébec, Canada H3C 3A7 2 -Snecma, site de Villaroche, Moissy-Cramayel, 77550, France 3 -Laboratoire de M écanique de Lille, Universit é de Lille, 59655 Villeneuve d'Ascq, France

Simulation num érique et exp érimentale d'un cas d'interaction rotor/stator sur une aube de compresseur haute-pression industriel

A. Batailly 1 , Q. Agrapart 2 , A. Millecamps 2 , J.-F. Brunel 3

R ésum é

Cet article porte sur la comparaison entre des résultats expérimentaux et numériques relatifs à un cas d'interaction rotor/stator provoqué par des contacts structuraux. Contrairement à ce qui a été observé dans le cadre de travaux précédents pour un compresseur basse-pression, cette interaction n'est pas divergente : de hautes amplitudes de vibration sont prédites numériquement et observées expérimentalement sur un court laps de temps. Dans un premier temps, une analyse détaillée des données expérimentales permet de définir l'interaction observée comme un cas de frottement aube/carter localisé sur une seule aube impliquant son premier mode de torsion. Les résultats numériques sont en accord avec les observations expérimentales : la vitesse d'interaction, les motifs d'usure sur le carter ainsi que la réponse vibratoire de l'aube sont prédits avec précision. Tout au long de cet article, la stratégie numérique utilisée est confrontée à une autre stratégie utilisée dans la littérature pour la simulation de ce type d'interactions. Les différences fondamentales entre les deux approches sont soulignées.

Mots-cl és

interaction rotor/stator; frottement aube/carter; dynamique non-linéaire; dynamique du contact; frottement; revêtement abradable; compresseur haute-pression

Introduction

The need for the improvement of modern aircraft engines efficiency motivates designers to consider reduced clearances between rotating and stationary components. As a counterpart, structural contacts between these components may occur more frequently and must now be accounted for as part of normal engine running conditions in non-accidental configurations. Subsequent vibratory phenomena may involve a single blade-it is then usually referred to as rubbing [1, 2, 3]-, a full bladed disk and the surrounding casing-with possible modal coincidence yielding very high amplitudes of vibration [START_REF] Schmiechen | Travelling wave speed coincidence[END_REF][START_REF] Batailly | Assessment of reduced models for the detection of modal interaction through rotor stator contacts[END_REF]-or precessional shaft motions such as backward or forward whirl motions [START_REF] Muszynska | Whirl and whip-rotor/bearing stability problems[END_REF][START_REF] Sinha | Rotordynamic analysis of asymmetric turbofan rotor due to fan blade-loss event with contact-impact rub loads[END_REF][START_REF] Salvat | Two-dimensional modeling of shaft precessional motions induced by blade/casing unilateral contact in aircraft engines[END_REF][START_REF] Parent | Whole engine interaction in a bladed rotor-to-stator contact[END_REF]. A usual strategy employed by engine manufacturers in order to mitigate such vibratory phenomena involves the deposition of an abradable coating along the casing contact surface [START_REF] Borel | The wear mechanisms occurring in abradable seals of gas turbines[END_REF][START_REF] Mandard | Interacting force estimation during blade/seal rubs[END_REF]. This coating acts as a sacrificial material that is worn out in the case of blade/casing contacts. However, even when this coating is employed, interaction phenomena leading to blade failure have been reported [START_REF] Millecamps | Influence of thermal effects during blade-casing contact experiments[END_REF].

Because they may threaten the engine components structural integrity, these interactions have been the focus of many recent studies, both numerical and experimental.

In this contribution, the focus is made on an experimental test bench set-up in Snecma facilities. It features a full-scale high-pressure compressor stage and aims at simulating contact induced interactions between one of the blades (slightly longer than the other ones) and the surrounding abradable coating that is deposited along the casing circumference. For this experimental set-up, it is found that the witnessed interaction involves a single blade-thus it should be analysed as a sequence of rubbing events-, more specifically its first torsional mode, which is its second free-vibration mode. The numerical simulation of rubbing events has been the centre of attention in a significant number of recent publications [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF][START_REF] Williams | Simulation of blade casing interaction phenomena in gas turbines resulting from heavy tip rubs using an implicit time marching method[END_REF][START_REF] Legrand | Numerical investigation of abradable coating removal through plastic constitutive law in aircraft engine[END_REF][START_REF] Ma | Vibration response analysis of blade-disk dovetail structure under blade tip rubbing condition[END_REF]. Several strategies may be employed to simulate such occurrences: a brute-force approach using a commercial finite element software with full 3D finite element models for the rotating components [START_REF] Arzina | Vibration analysis of compressor blade tip-rubbing[END_REF] is sometimes considered. This approach advantageously allows for a very precise modelling of the problem (possibly featuring nonlinear deformations, centrifugal and gyroscopic effects. . . ) but inevitably yields very long computation times and only a handful of rubbing events-at most-may be numerically simulated, a phenomenological approach that first aims at characterizing the contact forces during the sequence of rubbing events before applying it on the blade tip as a function of time [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF][START_REF] Sinha | Non-linear dynamic response of a rotating radial timoshenko beam with periodic pulse loading at the free end[END_REF]. This approach is also commonly used for the modelling of machining operations such as cutting and milling [START_REF] Campocasso | A generalised geometrical model of turning operations for cutting force modelling using edge discretisation[END_REF][START_REF] Gao | Mechanistic identification of cutting force coefficients in bull-nose milling process[END_REF]. Using this phenomenological approach, a very long sequence of rubbing events may be simulated but the restriction of rubbing induced contact forces to a pulse load [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF][START_REF] Sinha | Non-linear dynamic response of a rotating radial timoshenko beam with periodic pulse loading at the free end[END_REF] fully known a priori filters contact related nonlinearities thus making it impossible to predict some potentially critical interaction phenomena. Additionally, when a pulse load is considered, strong assumptions are often made on the way the load is applied along the blade tip, simplified analytical models are also frequently used for the simulation of rubbing events [START_REF] Yuan | Compact-impact analysis of a rotating geometric nonlinear plate under thermal shock[END_REF][START_REF] Ma | Vibration response analysis of blade-disk dovetail structure under blade tip rubbing condition[END_REF][START_REF] Shen | Numerical analysis of a rub-impact rotor-bearing system with mass unbalance[END_REF]. In this case, contact conditions are typically treated by a penalty law: on one hand residual penetrations are tolerated and wear removal cannot be physically accounted for but on the other hand, sophisticated bifurcation analyses may be carried out under the assumption that interaction motions are periodic, and interesting qualitative results may be obtained in order to characterize the nature of interaction phenomena, a reduced-order model based [START_REF] Legrand | Numerical investigation of abradable coating removal through plastic constitutive law in aircraft engine[END_REF][START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF][START_REF] Almeida | Numerical analysis of bladed disk-casing contact with friction and wear[END_REF] approach, which is used in this article, is the last known method to simulate rubbing events. It combines the use of precise 3D finite element models of the rotor with a full 3D contact treatment procedure with no a priori assumptions on the rubbing forces. This quantitative approach aims at providing accurate results but remains more computationally costly than the use of simplified analytical models for instance. This strategy has already been applied to the simulation of a bending induced interaction [START_REF] Millecamps | Influence of thermal effects during blade-casing contact experiments[END_REF] for the last stage of a low-pressure compressor [START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF]: the critical numerical speed, wear lobes as well as blade dynamics were accurately predicted. It advantageously accounts for abradable coating removal and contact treatment is based on a Lagrange multiplier computation [START_REF] Carpenter | Lagrange constraints for transcient finite element surface contact[END_REF] thus no residual penetrations can occur between the structures. So far, there is no existing comparison between these methods and only a few of them have been calibrated or confronted with experimental observations. In the context of this article and the numerical simulation of several hundreds of rubbing events for a given blade, both the brute-force approach and the use of simplified analytical models seem irrelevant: the former is not admissible due to the cumbersome computation times it would imply while the latter can not capture accurately the blade dynamics. Accordingly, this article gives the opportunity to confront the phenomenological approach with the reduced-order model based approach. Beside of the presentation of experimental and numerical results, different sections of this article focus on the distinction between these two strategies.

The first section of the article describes the test bench and the targeted interaction scenario, it also provides details regarding the blade instrumentation: the type of gauges used during the experiment is briefly exposed. Experimental results and observations are then listed in the second section of the article. In particular, the second section also contains an in-depth analysis of the stress signals retrieved on the longer blade during distinct phases of the interaction. Following sections are related to the numerical investigations carried out in order to numerically predict the experimental results. The third section recalls the key points of the reduced-order model based strategy and the focus is made on the computation of contact forces that are compared with pulse loads used in the phenomenological approach. Finally, in the last section, numerical results are scrutinized both in time and frequency domains and confronted with experimental observations. Note Amplitudes are normalized with respect to reference data for confidentiality purposes.

Experimental set-up

Test bench

The experimental set-up presented in this article is part of the European research program NEWAC1 dedicated to the reduction of CO 2 emissions and fuel consumption. More specifically, the experimental set-up is used to investigate the interaction between a full scale bladed disk of an axial high-pressure compressor and its surrounding casing as pictured in Fig. 1. The bladed disk of interest is machined from a single solid piece of material, also named blisk. It is rotationally driven by an electric motor via a flexible shaft and sealed bearings. The blisk is surrounded with a thick simplified casing on which an abradable coating is thermally sprayed, see Fig. 2. The alignment is a crucial point of the experimental set-up, it has been the subject of an accurate calibration at different angular speeds to eliminate parasitic vibratory phenomena. Adjustment screws and laminated shims allow for compensating angular misalignment and concentricity deviation between rotor and stator.

high-pressure compressor stage of interest

The rotor is cantilevered at the end of the shaft, in such a way that it allows space for an additional infrared thermal imaging camera, and easy access to the casing for clearance adjustment. The test bench is housed in a vacuum chamber in order to limitate aerodynamic effects. All tests are conducted at room temperature.

Blade of interest

For all blades, clearances at rest are identical from the leading edge to the trailing edge before the engine starts. As mentioned above, one of the blades is slightly longer than the other ones (additional length is about 1 % of the blade height) so that only this blade may impact the surrounding abradable coating deposited on the casing. The clearances are adjusted so that contacts only occur between the longer blade and the abradable coating once the targeted interaction speed Ω c has been reached. In order to better apprehend the vibratory behaviour of the blade during the interaction, it is instrumented with two dynamic strain gauges located on its pressure side as depicted in Fig. 3. Each gauge position and orientation is carefully chosen in order to isolate the first bending mode (1B) and the first torsional mode of the blade (1T). Accordingly, these gauges will be referred to as gauge 1B and gauge 1T in the following.

In order to ensure that the measured interaction is located on the longer blade only, the two adjacent blades and the disk are also instrumented with dynamic strain gauges. Time signals acquired by these gauges are not depicted in this article as they only feature negligible amplitudes of vibration. Also, it is assumed that the casing-which is particularly thick-is perfectly rigid thus casing amplitudes of vibration are not monitored.

Interaction scenario

The Campbell diagram pictured in Fig. 4 shows the evolution of the blade first three eigenfrequencies (the first bending mode 1B, the first torsional mode 1T and the second bending mode 2B). The dashed lines refer to the eigenfrequencies of the real blade measured experimentally at rest. The solid lines refer to the eigenfrequencies of the blade numerical model used in the following. The increase of the computed eigenfrequencies (solid lines only) over the angular speed range of interest is a consequence of centrifugal effects that are accounted for in the modal analysis. Small differences visible between the eigenfrequencies of the numerical model and the experimental blade stem from possible experimental measurements uncertainties as well as the simplification of the boundary conditions used for the numerical model. The latter is simply clamped on its root with a perfectly rigid disk.

It can be seen in Fig. 4 that the experimental angular speed Ω c corresponds to the intersection of the 22 nd engine order line with the frequency line of the first torsional mode. It is also located in the vicinity of the intersection of the 8 th engine order line with the frequency line of the first bending mode. The same observations stand for the numerically predicted angular speed Ω * c (see section 5) for the numerical blade model.

The crossing between an engine order line and an eigenfrequency line is typically an area of interest when investigating the vibratory response of a blade subjected to unilateral contacts conditions. As a matter of fact, this crossing implies the following relation between the eigenfrequency of the blade f i and the angular speed Ω:

f i = kΩ, k ∈ N (1)
meaning the blade can oscillate exactly k times over its i th free-vibration mode during one revolution. The blade is thus excited over a sub-harmonic of its i th eigenfrequency which may lead to a resonance in a highly nonlinear context [START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF][START_REF] Batailly | Conjectural bifurcation analysis of the contact-induced vibratory response of an aircraft engine blade[END_REF].

The targeted angular speed Ω c is approached step by step as shown in Fig. 5(a). A first acceleration quickly brings the rotor to 90 % of Ω c over the first 70 seconds of the experiment. Then, the angular speed is progressively increased and Ω c is reached around t = 508 s as pictured in Fig. 5(b). The first contacts are detected around t = 508 s and self-sustaining interaction stops by itself at t = 620 s. Between these two instants, which define the time window of interest in this article, the vibratory behaviour of the blade is decomposed in six distinct phases-phase 1 to phase 6-marked in Fig. 5(b) and detailed in the next section. shape of these signal does not allow for distinguishing the contact phases mentioned above. Indeed, these phases have been identified based on the analysis of the signals acquired by both gauges in the time domain and plotted in polar coordinates as shown below.

Experimental results and observations

Analysis in the time domain: polar coordinates visualization

In order to go further in the understanding of the interaction, the time signals plotted in Fig. 6 and Fig. 7 are represented in polar coordinates: time is replaced by the corresponding angular position of the blade. This way, the signals obtained from a blade revolution to another may be compared, the observation of distinct patterns leads to the definition of the contact phases mentioned above. In the following, a different scale is used for the representation of stresses values from a phase to another with the sole aim to show distinguishible signals during each phase. For all the figures depicted in this section, the zero stress level is indicated by a black circle (σ = 0) in polar coordinates. Phase 6: t ∈ [618 ; 620] s, the stress levels have reached a maximum, it is thus assumed abradable removal is also maximum during this phase. Contrary to phase 5, the time signals obtained for each gauge do not clearly show 22 oscillations anymore, hinting that the perfect synchronization with the angular speed has been lost, see Fig. 14. This will ultimately lead to contact loss. Phase 7 : t ∈ [620 ; 900] s, contact has been lost, amplitudes of vibration are progressively damped until the blade comes back to its equilibrium position, the observations made during this phase are similar to those made during phase 0.

Frequency analysis for gauge 1B

The time signal x(t) acquired by gauge 1B is pictured in Fig. 6, its empirical mode decomposition (EMD) [START_REF] Huang | The empirical mode decomposition and the hilbert spectrum for nonlinear and nonstationary time series analysis[END_REF][START_REF] Ming | Multi-fault diagnosis of rotor system based on differential-based empirical mode decomposition[END_REF] is shown in Fig 15 with :

x(t) = 3 i=1 γ i (t) + r(t) (2)
where γ i (t) is an intrinsic mode function (IMF) and r(t) is a residual monotonic function featuring negligible amplitudes. The EMD of the time signal allows for locally discriminating the frequency content of a signal: at all time, the highest frequency content of the signal is retrieved in the first IMF γ 1 , second highest frequency content appears on IMF γ 2 and so on, until amplitudes become negligible. The IMF of a time signal are particularly well suited for a Hilbert transform [START_REF] Huang | The empirical mode decomposition and the hilbert spectrum for nonlinear and nonstationary time series analysis[END_REF] and the Hilbert spectrum of x(t) is plotted in Fig. 16. This spectrum is a representation of the amplitude of each IMF in a time/frequency plane using a colour code from white (lowest amplitudes) to dark blue (highest amplitudes). Over time, the relative intensity of the colour lines allows for determining which frequency content is dominant.

As one can see in Fig. 16, free vibration modes 1B and 2B are dominant right after the initiation of the interaction: the first contacts yield bending motions. During phase 2, the participation of the first torsional mode becomes non-negligible and mode 1B is clearly dominant. Starting from the beginning of phase 3, the blue line around mode 2B vanishes and the frequency content captured

Frequency analysis for gauge 1T

The signal acquired by gauge 1T is pictured in Fig. 7 and its EMD in Fig. 17. The comparison of the EMD pictured in Figs. [START_REF] Arzina | Vibration analysis of compressor blade tip-rubbing[END_REF] and 17 reveals that non negligible amplitudes of vibration for the signal acquired by gauge 1T are solely located on its first IMF γ 1 while at least three non-negligible IMF are visible in Fig. 15. This is an indication that non-negligible amplitudes of vibration are predicted for only one frequency during each phase of the interaction. The Hilbert spectrum plotted in Fig. 18 confirms this observation : almost all non negligible amplitudes of vibration are detected along the first torsional mode 1T. During the first phase, both modes 1T and 2B feature non-negligible amplitudes. This participation of mode 2B is key in understanding why the blade dynamics cannot be synchronized with its angular speed in Figs. 9(a) and 9(b). Indeed, it is remarkable that the contribution of mode 2B vanishes at the end of phase 1, and during phase 2 the frequency content of signal acquired by gauge 1T is dominated by the first torsional mode 1T. As the blue line around

Abradable wear

The surface condition of the abradable coating is scrutinized after the interaction using a coordinate measuring machine. Wear levels are measured using the following resolution: 180 equidistant points in the tangential direction and 44 equidistant points in the axial direction. Measured data yield the wear map depicted in Fig. 19: a colour code from light blue (no wear) to dark blue (maximum wear) is used. Black areas correspond to locations that could not be measured due to the use of anchored capacitive sensors in the casing, see Fig. 2. The wear map underlines that the coating underwent severe wear during the interaction. 22 lobes are evenly located along the casing circumference in front of the blade leading edge and trailing edge. This observation is in agreement with the time domain analysis carried above. Remarkably, no significant wear level is measured in front of the middle of chord. When having a close look at the angular location of the wear lobes, one may notice Based on this observation and the time domain analysis carried out previously. One can only speculate that contacts must have happened two times (one contact at the leading edge followed by one contact at the trailing edge) per oscillation of the blade along its first torsional mode. In order to corroborate this hypothesis, linear and nonlinear static computations of the blade/casing clearances as the blade vibrates along its first torsional mode are depicted in Fig. ). The opposite behaviour is predicted along the negative alternation: clearance is closed at the trailing edge while it is open at the leading edge. It is noticeable that clearance closure is much higher at the leading edge than it is at the trailing edge which is consistent with the observation of deeper wear lobes at the leading edge. Also, these calculations reveal that clearance variation at the middle of chord is negligible in agreement with the fact that no wear is observed at this location experimentally.

Following the disassembly, the blisk exhibited neither damage nor substantial friction wear. Contrary to what may have been reported in the literature [START_REF] Millecamps | Influence of thermal effects during blade-casing contact experiments[END_REF] within a low-pressure compressor, this interaction phenomenon stopped before initiating any significant structural damage on the contacting blade.

Refined stress analysis

A close look at the stress signals plotted in polar coordinates in Figs. 10(a), 11(a), 12(a) and 13(a) reveals that they are all in phase. However, a significant phase shift occurs between phases 5 and 6, see Figs. 13(a brightness, increasing contrast and using their negative) helps visualizing the actual rotation of the contact areas between phases 5 and 6. The average of the angular shift measured for each contact areas lead to the same value δθ 0.065 rad. Because phases 5 and 6 are of extremely short duration (a few seconds only), it is difficult to accurately identify when the angular shift starts and when it stops. In particular, it is unclear either this shift happens progressively through out phase 6 or happens more suddenly at the transition between phases 5 and 6.

In any case, this angular shift is coincidental with the blade dynamics desynchronization that leads to contact loss. It will thus be one of the point of interest in the following numerical investigation. As a matter of fact, an angular shift during a blade/casing interaction simulation was previously identified as a possible cause for the desynchronization of the blade dynamics with respect to its angular velocity [START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF].

Numerical modelling

As mentioned in the introduction, the numerical strategy employed in this article is the reduced-order model based approach [START_REF] Legrand | Numerical investigation of abradable coating removal through plastic constitutive law in aircraft engine[END_REF][START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF]. This strategy has been extensively detailed in previous contributions and is only briefly recalled here for the sake of brevity. It is based on an explicit time-marching algorithm [START_REF] Carpenter | Lagrange constraints for transcient finite element surface contact[END_REF] with a Lagrange multiplier-based contact treatment procedure ensuring that no penetration can occur between the impacting components. When contacts occur between a blade and the abradable coating, a plastic consitutive law [START_REF] Legrand | Numerical investigation of abradable coating removal through plastic constitutive law in aircraft engine[END_REF] is used for modeling abradable coating removal. While this strategy allows for taking into account a full bladed disk and a flexible casing, it is here limited to a single blade and the casing is assumed perfectly rigid in agreement with experimental observations.

Numerical framework

Same as suggested by experimental measurements, the casing is slightly ovalized. Numerically, this ovalization is used in order to initiate contacts with the blade. Numerical simulations are carried out with a constant angular speed Ω. At t = 0 s, the casing is perfectly circular and there is no penetration between the blade and the casing. Blade tip displacements due to centrifugal effects and casing ovalization are progressively accounted for during the first 20 % of the simulation during which contact is initiated. After this point, the casing ovalization remains constant and the blade is free to interact with the abradable coating and the rigid casing.

The parameters of the simulation are as follows:

• friction is accounted for, µ = 0.15, this coefficient is used for both blade/abradable coating contacts and blade/casing contacts if the abradable is locally fully removed, • the amplitude of the casing ovalization is equal to the blade/casing clearance at rest, • the abradable coating is 5 mm thick, • simulations are carried out over 100 blade revolutions, • the blade modal damping coefficient is ξ = 5 • 10 -4 , • empirical values are considered for the mechanical properties of the abradable coating, these values were obtained from the calibration of the employed numerical strategy with various experimental results [START_REF] Cuny | Experimental measurement of dynamic forces generated during short-duration contacts: application to blade-casing interactions in aircraft engines[END_REF][START_REF] Mandard | Interacting force estimation during blade/seal rubs[END_REF] : -Young's modulus is E = 20 GPa, the plastic modulus is K = 5 GPa, the yield limit is σ y = 1.5 Pa. The numerical convergence of the simulations has been checked with respect to the time step, the number of abradable elements as well as the blade reduction parameter but is not detailed here for brevity.

Targetted angular speed range

In order to numerically predict a rotor/stator interaction, it is first required to estimate the angular speed range over which it may occur. The Campbell diagram plotted in Fig. 4 combines the eigenfrequencies of the physical model and the numerical model. It appears that the targetted interaction is anticipated for Ω * c < Ω c since the eigenfrequency of the first torsional mode is slightly lower than the one of the physical model.

However, the stiffening effect stemming from unilateral contact constraints typically shifts the predicted interaction speed toward higher values. For that reason, numerical simulations are carried out over an angular speed range Ω ∈ [0.97Ω * c ; 1.10Ω * c ] centered on a higher value than Ω * c . This angular speed range also corresponds to Ω ∈ [0.9Ω c ; 1.02Ω c ].

In total, 101 simulations are carried out over this angular speed range, the angular speed increment between two simulations is δΩ = 0.0012Ω c .

Modelling

Blade

The reduced-order model of the blade used for the simulation is computed from a 3D industrial finite element model clamped on its foot. Eight nodes are used along the blade tip in order to manage contact with the abradable coating and the casing. A Craig-Bampton based component mode synthesis method [START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF] is used in order to obtain a numerically efficient model. The reduced-order model contains 115 degrees of freedom, its convergence with respect to the inital finite element model with respect to the reduction parameter has been carefully checked, its first eigenfrequencies are shown in Fig. 4.

Rub induced contact forces

As mentioned in the introduction, a typical way to account for rub induced contact forces is to consider a pulse load [START_REF] Sinha | Non-linear dynamic response of a rotating radial timoshenko beam with periodic pulse loading at the free end[END_REF][START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF] along the blade tip. In this study however, a full 3D nonlinear contact treatment procedure [START_REF] Carpenter | Lagrange constraints for transcient finite element surface contact[END_REF] is employed and this section aims at precising in which context the two approaches may lead to similar results or not. First of all, the pulse load model assumes that the contact force has a specific shape, be it a sine [START_REF] Sinha | Non-linear dynamic response of a rotating radial timoshenko beam with periodic pulse loading at the free end[END_REF] or a quadratic function [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF] of time. Numerically, both definitions are actually very close-as pictured in Fig. 23. They are compared with a typical contact force pulse numerically predicted by the numerical approach considered in this article. Even if significant differences are visible, in particular the predicted contact force is not perfectly ) with the pulse load model as proposed in [START_REF] Sinha | Non-linear dynamic response of a rotating radial timoshenko beam with periodic pulse loading at the free end[END_REF] ( ) and in [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF] ( ) symmetric with respect to the time t Fmax where F = F max , the pulse load approach provides an acceptable approximation of the rub induced contact force. Yet, in order to get this acceptable approximation, one must be able to quantify very accurately: (1) the duration of the rub event ∆t and (2) the maximum related contact force F max . Such quantification is unachievable for the experimental test bench considered in this article. In the scientific literature, typical values of ∆t given in the literature are around 1 ms [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF] but are related to large fan blades. For smaller, and stiffer high-pressure compressor blades, no estimation of the rub event duration could be found. Numerical simulations carried out in the following suggest a significantly lower duration: ∆t 0.15 ms. As mentioned earlier, the casing is slightly ovalized due to assembly conditions and two symmetric privileged contact areas are thus expected. A typical definition of the pulse load would then feature two rub events per revolution as pictured in red in Fig. 24. To the contrary, at the critical angular speed, the employed numerical strategy predicts that a much larger number of rub events occur, see the blue plot in Fig. 24. This observation shows that the interaction experimentally witnessed could not be detected with a simplified pulse load and justifies a posteriori the need for a more versatile approach of the contact forces computation. In addition, the noticeable variation of the amplitude of the contact force numerically predicted from a rub event to another hints that considering a unique F max value over time may lead to inaccurate results. for Ω ∈ [0.900Ω c ; 0.902Ω c ]: at the very beginning of the angular speed range, significant wear levels with 9 lobes are numerically predicted. This interaction is likely to be related to the crossing of the ninth engine order line with the first bending mode frequency line, see Fig. 4 around 0.856Ω c . Indeed, it has been previously observed that bending initiated interactions typically feature high amplitudes of vibration and deeper wear lobes than torsional initiated interactions [START_REF] Batailly | Conjectural bifurcation analysis of the contact-induced vibratory response of an aircraft engine blade[END_REF]; for Ω ∈ [0.907Ω c ; 0.909Ω c ]: a narrow band of wear is noticeable, the relatively small wear levels combined with the localization of this phenomenon indicates that it may be related with a high frequency mode (bending or torsion) or potentially a sub-harmonic of a high frequency mode.

Numerical results

The impossibility to distinguish a precise number of lobes makes it arduous to analyze further; for Ω ∈ [0.941Ω c ; 0.950Ω c ]: a narrow interaction area with a high number of lobes is predicted, zoom over this area are provided for both the leading edge and the trailing edge in Fig. 25(c) and in Fig. 25(d). These figures allow for a precise count of the number of lobes: 22 wear lobes evenly spaced along the casing circumference are visible. The colour code indicates that higher wear levels are observed at the leading edge; for Ω ∈ [0.983Ω c ; 0.993Ω c ]: eight lobes of relatively low depth are predicted; for Ω ∈ [0.993Ω c ; 1.05Ω c ]: eight deep lobes are predicted through out this angular speed range, which coincidates with the crossing of the eighth engine order line with the first bending mode frequency line in Fig. 4. The analysis of such interaction goes beyond the scope of this article, the reader may refer to [START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF] for more details. In the following, all the numerical results are obtained for Ω * c,a = 0.946Ω c .

Wear profiles at the critical speed Ω *

c,a

Numerical/experimental confrontation

A polar representation of the wear profiles at the leading edge, the middle of chord and the trailing edge is pictured in Fig. 26. The wear profiles are pictured with respect to time (the darker profile ( ) is the profile at the end of the simulation) and superimposed with the casing deformation ( ) so that one may see the privileged contact areas. The final wear profiles in Fig. 26 correspond to the ones plotted in Fig. 25 for Ω * c,a : wear lobes are deeper at the leading edge than at the trailing edge. In agreement with experimental observations, a very small wear level is predicted in front of the middle of chord with only two wider lobes.

The wear profiles obtained from the leading edge to the trailing edge for Ω * c,a are combined in Fig. 27. The colour code used in Fig. 27 is identical to the one used in Fig. 19, it is relative: from highest wear levels pictured in dark blue to lowest wear levels in light blue. The comparison of Figs. 19 and 27 underlines that numerical results and experimental observations are consistent in terms of relative wear levels: in both cases the deepest lobes are located at the leading edge (LE) while the lobes at the trainling edge (TE) are thinner and less deep. Also, both experimental observations and numerical results suggest low to no wear facing the middle of chord. Finally, the same angular shift is found between the location of the wear lobes in front of the trailing edge and the location of the wear lobes in front of the leading edge.

All these observations confirm the fact that the linear-in the sense that the small perturbation framework is considered-finite element blade model used for the numerical simulations provides an accurate approximation of blade/casing clearances, as pictured in Fig. 20.

Refined analysis of abradable removal

The wear lobes pictured in Fig. 26(a) and the zoom in Fig. 26(b) are drawn in Fig. 28 and 29 using cartesian coordinates for the sake of readability. The colour code used for the wear profiles is modified in order to highlight distinct phases of the simulated interaction : for t ∈ [0 ; 0.12] s: the amplitude of vibration of the blade is low and only two areas are worn out on the casing surface ( ). Those lobes can only be seen in Fig. 28; for t ∈ [0.12 ; 0.28] s: as the blade amplitudes of vibration grow, 22 lobes are worn out along the casing circumference. It is noticeable that both the angular location and the depth of these lobes grow linearly with time ( ); for t ∈ [0.28 ; 0.36] s: over this short period of time, the 22 wear lobes become twice deeper (from about 0.4 mm deep to 0.8 mm deep in Fig. 29) while they are still rotating along the casing circumference. The angular shift of the lobes during this phase is about δθ 0.065 rad ( ); The observation of an angular shift δθ of the wear lobes soon before the interaction stops reinforces the consistency of numerical results with experimental observations. At this point, one can only be very careful with respect to the angular values obtained experimentally and numerically: the limited experimental data does not allow for drawing conclusions. The great similarity of the values of δθ mentioned in sections 3.3 and 5.2.2 is intriguing but it is too early to conclude that experimental observations and numerical predictions are indeed related to the same phenomenon. From a numerical standpoint, among the possible explanations, one may wonder if this angular shift originates from a bifurcation-similar to the ones witnessed in [START_REF] Batailly | Conjectural bifurcation analysis of the contact-induced vibratory response of an aircraft engine blade[END_REF] with the same blade-of the system due to the progressive wear of the abradable coating. Whereas from an experimental standpoint, one may wonder if a local alteration of the abradable mechanical properties-due to thermal effects-could be yielding a desynchronization of the blade dynamics with its angluar speed resulting in such angular shift, as anticipated in [START_REF] Batailly | Numerical-experimental comparison in the simulation of rotor/stator interaction through blade-tip/abradable coating contact[END_REF].

Frequency analysis

The radial displacement of the blade leading edge is pictured in Fig. 30. Soon after blade/abradable coating contact is initiated-around t = 0.05 s-significant amplitudes of vibration are predicted. These amplitudes remain steady over a few revolutions before they are amplified and reach a maximum for t = 0.55 s. At this precise moment, the blade did about 70 revolutions. Soon after, for t > 0.6 s, the amplitude of vibrations starts decreasing due to structural damping which indicates that contact has been lost. In order to get a better understanding of the blade dynamics during the predicted interaction, the EMD of the above time response is carried out. The first three intrinsic mode functions γ i are depicted in Fig. 31. Other IMF feature negligible amplitudes and are not considered in the remainder. As mentioned in section 3.1.3, these IMF are particularly well-suited for computation of the Hilbert spectrum which is depicted in Fig. 32. Significant amplitudes of vibration are represented in a time/frequency plane using the following colour code: from white (negligible amplitudes) to dark blue (highest amplitudes). Overall, two lines appear in this spectrum : (1) the first one is located slightly above the frequency of the first torsional mode 1T, it features the highest amplitudes, and (2) a second one is located slightly above the frequency of the second torsional mode 2T. This spectrum underlines that the numerically predicted interaction is almost of the highest amplitudes of vibration are slightly over the frequencies of the torsional modes for t < 0.6 s because of contact stiffening. As a matter of fact, for t > 0.6 s, after contact has been lost, these amplitudes are found closer to the frequency lines as contact stiffening disappeared. The comparison of the Hilbert spectrum in Fig. 32 with the two spectra pictured in Fig. 16 and 18 underlines that the numerically predicted interaction seems to be consistent with experimental observations during phases 3 to 7 essentially. From a numerical standpoint, and in agreement with experimental data, it is predicted that the interaction starts along bending modes but related amplitudes of vibration are so small that they are hardly visible in Fig. 32 for t ∈ [0 ; 0.12] s.

However, numerical simulations fail to predict any participation of the second bending mode 2B. The participation of a mode or another in the transient phase following the first contacts is closely related to the contact configuration. It is possible that the idealized configuration used numerically is responsible for this approximation. Regarding the strong participation of the second torsional mode 2T, this numerical prediction cannot be confronted to experimental data since the eigenfrequency of this mode is higher than the acquisition frequency of gauges 1B and 1T. Because experimental and numerically predicted wear profiles are consistent, the participation of this mode in the blade dynamics on the experimental set-up cannot be ruled out.

Contact forces

The contact forces are computed in a local frame specific to each boundary node from the leading edge to the trailing edge as pictured in Fig. 33. A global cartesian coordinate frame (O, x, y, z) is considered for the whole bladed disk which rotates around the Oz direction. From this global frame, a local contact frame is defined by accounting for the casing conicity (angle φ) and the angular location of each contact node (angle ψ). As an example, the local frame used at the trailing edge is represented in Fig. 33.

Contact forces computed on the blade leading edge during the interaction in the three directions of space-radial r c , tangential θ c and axial z c -are pictured in Fig. 34 drive the blade tip off the abradable coating when contacts occur. Thus, based on the convention depicted in Fig. 33, this implies that the radial force is positive at all time. Also, permanent sliding is assumed between the blade and the abradable coating due to the very high bladed disk angular speed. It may be shown that the blade tip displacement is essentially driven by the rotation of the bladed disk: the blade tip local vibrations stemming from contact interactions generate displacement speeds in the tangential direction that are one order of magnitude lower than the displacement speed due to the bladed disk rotation. For that reason, the tangential contact force is negative at all time.

Regarding the axial force however, it may be positive or negative depending solely on the orientation of the axial blade tip displacement as there is no bladed disk motion in this direction.

In the phenomenological approach, an assumption sometimes made in the numerical simulation of rubbing events [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF] is to consider that the amplitude of the axial force is much lower than the one of the tangential force, so that the axial force may be neglected. While this assumption seems valid in our case for the trailing edge, on the leading edge the amplitudes of both the tangential and axial forces are found to be of the same order of magnitude. A finer analysis of the axial force goes beyond the scope of this study but these results do not necessarily go against the assumption made in [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF]. As a matter of fact, it is possible that the non negligible axial force predicted numerically is inherent to the blade design (high compressor blade while a large fan blade is considered in [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF]) or the torsional induced simulated interaction (while bending motions are expected in [START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF]). In any case, the contact forces depicted in Figs. 34 and 35 are consistent with the displacement shown in Fig. 30: their amplitudes grow until t = 0.45 s. Also, the loss of contact soon before t = 0.6 s is confirmed both on the leading edge and the trailing edge as all contact forces become null. It is noticeable that the very last contacts are witnessed on the leading edge.

Contact initiation

At the beginning of the interaction, over the revolutions following the very first impacts in the abradable coating (around t = 0.1 s), the blade dynamics is dominated by bending motions eventhough their relative low amplitudes make them hardly noticeable in Fig. 32. Six snapshots of the Von Mises stress field within the blade at the beginning of the interaction are depicted in Fig. 36. For this particular blade profile, it is found that the dominance of bending motions implies that highest stresses are located along the leading edge close to the blade root (area A). Other significant stresses are predicted within the blade closer to the trailing edge at about a third of the blade height (area B).

Conclusion

This article describes experimental and numerical investigations carried out for the simulation of a rotor/stator interaction case localized on a single blade within the high-pressure compressor of an aircraft engine. This interaction case is initiated by structural contacts between a blade and the surrounding casing on which is deposited an abradable coating. Because only one blade of the bladed disk is involved in the interaction, the interaction may be qualified as a rubbing event. Experimental observations underline the domination of the first torsional mode-which is the second free-vibration mode of the blade-in the blade vibratory response. Also, an extensive analysis, in both time and frequency domains, of time signals acquired during the interaction allows for the definition of distinct phases that reveal a sophisticated blade dynamics through out the interaction. It is remarkable that the simulated interaction does not diverge: it does not lead to blade failure contrary to bending induced interactions previously observed.

A reduced-order model based numerical strategy is employed for the simulation of this rubbing event. The confrontation of numerical results with experimental observations highlights the relevance of this strategy. In agreement with experimental observations, numerical results feature (1) a blade dynamics dominated by bending modes following the first contacts before torsional modes become dominant, (2) an abradable coating significantly worn out with 22 wear lobes in front of the blade trailing edge, 22 deeper lobes in front of the blade leading edge and practically no wear in front of the middle of chord, and (3) a non-divergent type of interaction motion where the blade amplitudes of vibration decrease at some point.

Key differences between the employed numerical strategy and a phenomenological approach relying on the experimental characterization of contact forces are pointed out through out the article. It is underlined that the identification of a pulse load applied linearly on the numerical model of a blade will filter contact induced non-linearities and make it impossible to simulate the rubbing event mentionned in this article. Though, for a single rubbing event, it is pointed out that traditional sine or quadratic approximations of the contact force are in very good agreement with the contact force numerically predicted with the reduced-order model based numerical strategy.

The confrontation of numerical and experimental results presented in this article implicitly raises the question of the influence of thermal effects in such rotor/stator interactions. As a matter of fact, an apparent paradoxical situation is described in this article: on one hand a significant raise of temperature is experimentally observed with the infrared camera but on the other hand a numerical strategy that does not account for thermal effects provide accurate results. Accounting for thermal effects within the employed numerical strategy is under investigation and goes beyond the scope of this article. One may only speculate that there is some consistency between the consequences of thermal effects on the experimental set-up and the interaction scenario numerically considered. In particular, it is observed that thermal effects are significant all along the casing contact surface. For instance, such thermal load on the casing leads to non-negligible deformations that may accentuate its ovalization in agreement with the idealized numerical interaction scenario.

Figure 1 .

 1 Figure 1. schematic cut view of an aircraft engine

  Fig. 3

Figure 2 .

 2 Figure 2. front view and cut view of the high-pressure compressor second stage used in the experimental set-up (most blades are hidden for confidentiality)

Figure 3 .

 3 Figure 3. detailed view of the longer blade, location of gauges

Figure 4 .

 4 Figure 4. Campbell diagram : experimental values at rest ( ) and numerical modal analysis ( )

3. 1 Figure 5 .Figure 6 .

 156 Figure 5. evolution of the bladed disk angular speed during the simulation: (a) whole time history, (b) zoom over the time window of interest

Figure 7 .Figure 8 .Figure 9 .

 789 Figure 7. time signal acquired by gauge 1T

Figure 10 .Figure 11 .Figure 12 .

 101112 Figure 10. superimposition of the stress signals during phase 2: (a) gauge 1T, (b) gauge 1B

Figure 13 .

 13 Figure 13. superimposition of the stress signals during phase 5: (a) gauge 1T, (b) gauge 1B

Figure 14 .

 14 Figure 14. superimposition of the stress signals during phase 6: (a) gauge 1T, (b) gauge 1B

Figure 15 .Figure 16 .

 1516 Figure 15. EMD of the time signal acquired by gauge 1B

Figure 17 .Figure 18 .

 1718 Figure 17. EMD of the time signal acquired by gauge 1T

Figure 19 .

 19 Figure 19. wear profiles measured on the casing

Figure 20 .

 20 Figure 20. clearance variation along the first torsional mode 1T, nonlinear computations: leading edge (LE) ( ), middle of chord (MC) ( ), trailing edge (TE) ( ) and linear computations: leading edge (LE) ( ), middle of chord (MC) ( ), trailing edge (TE) ( )

Figure 21 .Figure 22 .

 2122 Figure 21. average of the stress signals retrieved during phases 5 () and 6 ( ) from gauge 1T

Figure 23 .

 23 Figure 23. zoom over one rub event : superimposition of the predicted contact force ( ) with the pulse load model as proposed in [3] ( ) and in [2] ( )

Figure 24 . 5 . 1 Figure 25 .

 245125 Figure 24. contact force predicted numerically () superimposed with a pulse load ( )[START_REF] Sinha | Non-linear dynamic response of a rotating radial timoshenko beam with periodic pulse loading at the free end[END_REF][START_REF] Turner | Airfoil deflection characteristics during rub events[END_REF]

Figure 26 .Figure 27 .Figure 28 .Figure 29 .

 26272829 Fig. 26(b)

Figure 31 .Figure 32 .

 3132 Figure 30. leading edge radial displacement

Figure 33 .

 33 Figure 33. conventions used in to define local contact frames

Figure 34 .

 34 Figure 34. contact force at the leading edge in the three directions of space r, θ and z

Figure 35 . 5 . 5

 3555 Fig. 24

Figure 36 . 5 . 5 . 2 Figure 37 .

 3655237 Figure 36. Stress field in the blade in the vicinity of t = 0.1 s where bending motions are dominant

http://www.newac.eu

Acknowledgement

Thanks go to Snecma for its technical and financial support.