Leonardo Sanches

Daniel Alazard

Guilhem Michon

Alain Berlioz

Robustness Analysis of Helicopter Ground Resonance with Parametric Uncertainties in Blade Properties

This paper presents a stability robustness analysis of the helicopter ground resonance phenomenon. By using the lifting procedure, the uncertain Linear Time-Periodic (LTP) model of the helicopter is transformed into an augmented uncertain Linear-Time-Invariant (LTI) model that allows the application of µ-analysis tools. The lifting procedure involves a periodic switching LTI piecewise model computed using oversampling of the system period. The representativeness of the lifted model for various oversampling period values and discretization methods is discussed and compared with a Floquet analysis for several parametric congurations. A µ-analysis is then applied to nd the worst case parametric conguration for a given rotor angular rate. The parametric uncertainties taken into account are the dynamic characteristics (stiness and damping) of each blade hinge. A signicant advantage of the proposed approach is that it enables performing ground resonance analysis for a rotor with dissimilar blade properties due to aging eects. Considering uncertainties on the four blade hinge sti-

nesses and damping factors, the µ-analysis performed on the lifted model leads to the conclusion that the worst case for degraded rotor stability corresponds to the symmetric perturbation of all the blades. Helicopter ground resonance is a dynamical phenomenon that has attracted much attention of from researchers over the last ve decades, particularly for rotors with hinged blades and a shift (oset) between the hinge axis and the main rotor axis. The prediction of critical rotor velocities at which the phenomenon occurs was rst studied by Coleman and Feingold [START_REF] Coleman | Theory of Self-Excited Mechanical Oscillations of Helicopter Rotors with Hinged Blades[END_REF] for helicopters with rotors with identical blade properties. The equations of motions were simplied [START_REF] Bir | Multi-Blade Coordinate Transformation and its Application to Wind Turbine Analysis[END_REF] by eliminating their periodical characteristic and Linear Time Invariant (LTI) stability analysis was performed easily.

Major contributions to understanding this phenomenon in hingeless and bearingless rotors have been made since [START_REF] Donham | Ground and Air Resonance Characteristics of a Soft In-Plane Rigid-Rotor System[END_REF][START_REF] Hodges | An Aeromechanical Stability Analysis for Bearingless Rotor Helicopters[END_REF]. Criteria for determining the dimension of viscous dampers have been established and the design of passive control systems has been studied in order to dissipate energies and avoid unstable motions [57]. Semi-active and active control solutions using the pitch angle of each blade have also been proposed to reduce vibrations [START_REF] Arcara | Periodic Control of Helicopter Rotors for Attenuation of Vibrations in Forward Flight[END_REF][START_REF] Bittanti | Periodic Active Control of Vibrations in Helicopter: a Gain-Scheduled Approach[END_REF].

However, the eects of aging on various mechanical elements can induce unbalanced parametric variations from one blade to another and compromise the rotor's nominal behavior, leading to dangerous conditions in extreme cases. In the eld of aeronautics such situations must be mastered to reduce not only human risk but also maintenance costs. Therefore analysis tools for rotors with dissimilar blades are required to assess ground resonance instability.

When considering blades with dierent mechanical properties, the simplications made by Coleman are no longer valid and Floquet theory has been used to study the stability of time-periodic equations of motion [START_REF] Hammond | An Application of Floquet Theory to Prediction of Mechanical Instability[END_REF][START_REF] Skjoldan | On the Similarity of the Coleman and Lyapunov-Floquet Transformations for Modal Analysis of Bladed Rotor Structures[END_REF]. Predicting the ground resonance phenomenon for a wide range of dissimilar blade congurations means analyzing each point individually on a grid of parametric space, generating high computational costs. Furthermore, there is no guarantee that parametric gridding includes the worst-case parametric conguration.

On the other hand, the stability and performance robustness of linear time-invariant (LTI) systems under structural uncertainties have been analyzed by using the standard µ-analysis method and ecient tools are now available [START_REF] Balas | µ-Analysis and Synthesis Toolbox: for Use with MATLAB, User's Guide[END_REF] [1] [START_REF] Bir | Multi-Blade Coordinate Transformation and its Application to Wind Turbine Analysis[END_REF] .

In [START_REF] Kalender | Control Design and Robustness Analysis of Linear Time-Periodic Systems[END_REF], the parametric robustness analysis of LTP systems is considered using a truncated point mapping technique and µ-analysis. the other in such systems was not studied. The robustness of uncertain polytopic discrete-time periodic systems was considered in [START_REF] Farges | Robust H2 Performance Analysis and Synthesis of Linear Polytopic Discrete-Time Periodic Systems via LMIs[END_REF] from the standpoint of periodic state-feedback design. The robustness analysis of ground resonance stability has been addressed more recently in [START_REF] Masarati | Robust Aeroservoelastic Stability of Helicopters: Application to Air/Ground Resonance[END_REF]. The authors considered complex uncertainties embedding uncertainties on the stiness and damping ratio of lead-lag dampers, but under the assumption that these dampers are identical from one blade to the other. However, the direct analysis of ground resonance stability for asymmetric rotor properties was dealt with unsatisfactorily.

Recently, the problem of robustness analysis of linear time periodic (LTP) dynamical systems [START_REF] Verdult | Identication of Linear Parameter-Varying State-Space Models with Application to Helicopter Rotor Dynamics[END_REF] under structured LTI uncertainties was solved [START_REF] Kim | Robustness Analysis of Linear Periodic Time-Varying Systems Subject to Structured Uncertainty[END_REF][START_REF] Rabenasolo | Robust Stability for Linear Periodic Systems under Structured Uncertainties[END_REF], by combining Floquet theory with the lifting technique [2426]. The original uncertain LTP system was cast in the form of a Linear

Fractional Transformation (LFT) using discretization on an oversampling of the system period in order to use µ-analysis methods. The time-lifted LFT model involves an uncertainty structure with highly-repeated parameters which can raise problems for performing µ-analyses. In order to reduce the size of the uncertainty block and the associated computational burden of µ-analyses,

[1] Ferreres, G. and Biannic, J.-M., The Skew Mu Toolbox", http://www.onera.fr/staff-en/jean-marc-biannic/.

[2] Peaucelle, D., RoMulOC: Robust Multi-Objective Control Toolbox",http://spiderman-2.laas.fr/OLOCEP/ romuloc/index.html.

the representativeness of the (lifted) LFT model for various oversampling period values and three discretization methods (zero-order hold (zoh), rst-order hold (foh) and Tustin method [START_REF] Franklin | Digital Control of Dynamic Systems[END_REF][START_REF] Hoagg | First-Order-Hold Sampling of Positive Real Systems and Subspace Identication of Positive Real Models[END_REF]) is discussed in comparison with Floquet analysis for several parametric congurations. A general MATLAB R function was developed to implement this lifting procedure to any LTP model.

The contributions of this paper are:

• to provide practical and low CPU time-consuming tools allowing the application of µ-analysis

to LTP systems and

• to analyze helicopter ground resonance stability under parametric uncertainties using these tools and µ-analysis.

The example studied corresponds to independent uncertainties on (lead-lag) blade hinge stinesses and damping factors. The method proposed can be applied to any kind of parametric uncertainty but the result of analysis on the sensitivity to blade hinge stinesses (damping) is worth mentioning.

It is shown that the worst case parametric conguration corresponds to an identical variation of each stiness (damping), i.e. a rotor with identical blade dynamic properties.

Section 2 describes the lifting procedure for an uncertain LTP system with particular emphasis given to the discretization method. In section 3, the dynamic model used to study the ground resonance phenomenon is derived, the lifting procedure is validated and the results of stability analysis using the method proposed and µ-analysis are presented. Section 4 presents the conclusions.

II. Robustness Analysis of LTP Systems

A.

General background

Consider the uncertain LTP system S(∆) dened by an LFT representation M(s, t) -∆:

M(s, t) :        ẋ(t) = A(t) x(t) + B(t) w(t) z(t) = C(t) x(t) + D(t) w(t) (1)
with: w(t) = ∆ z(t)

where x(t) ∈ R n is the state vector, w(t) and z(t) are input and output vectors of the LFT. The structured uncertainty matrix ∆ is a p×p diagonal matrix of unknown but bounded real parameters:

∆ = diag [δ 1 , δ 2 , . . . , δ p ] (2)
Matrices A(t), B(t), C(t) and D(t) are real, piecewise continuous and periodic with a period T ,

A(t + T ) = A(t), B(t + T ) = B(t), C(t + T ) = C(t), D(t + T ) = D(t) . (3)
The closed-loop LFT representation is

ẋ(t) = A(t) + B(t)∆(I n -D(t)∆) -1 C(t) x(t) = A p (t, ∆) x(t) (4)
where the matrix A p (t, ∆) is also T -periodic.

The nominal system (∆ = 0) is assumed to be stable. The parametric robustness analysis consists in nding the smallest uncertainty ∆ worst matrix which makes the closed-loop system (4)

unstable.

Floquet theory [START_REF] Gökçek | Stability Analysis of Periodically Switched Linear Systems Using Floquet Theory[END_REF][START_REF] Bittanti | Periodic Systems: Filtering and Control[END_REF] can be used to analyze the stability for a particular value of ∆. Considering the transition matrix Φ(t, t 0 , ∆) associated with the closed-loop system (4), the stability analysis is then characterized by the monodromy matrix R(t 0 , ∆) dened as the transition matrix over one period: R(t 0 , ∆) = Φ(t 0 + T, t 0 , ∆) .

(

Without loss of generality, it can be assumed that t 0 = 0.

Then, the system (4) is exponentially stable if and only if R(∆) = R(t 0 = 0, ∆) is Schur, i.e., all the eigenvalues of R(∆), also called characteristic multipliers: λ i (∆) i = 1, 2, . . . , n, have a magnitude less than one.

In most practical cases matrix R(∆) cannot be determined analytically. Nevertheless, R(∆) can be approximated by assuming that the system in Eq.( 4) can be represented in the form of a periodic linear switched system dened by [START_REF] Gökçek | Stability Analysis of Periodically Switched Linear Systems Using Floquet Theory[END_REF]:

ẋ(t) = A p (kh, ∆)x(t) (6)
∀ t ∈ [lT + kh, lT + (k + 1)h[, l = 0, 1, 2, . . . , k = 0, 1, 2, . . . n h -1 where h = T n h
is the oversampling period and n h is a positive integer. That is to say the system is assumed to be LTI during the oversampling period h.

The LFT representation of the switched system ( 6) is:

M(s, kh) :

       ẋ(t) = A(kh) x(t) + B(kh) w(t) z(t) = C(kh) x(t) + D(kh) w(t) (7)
with:

w(t) = ∆ z(t) ∀t ∈ [kh + lT, (k + 1) + lT [, l = 0, 1, 2, . . . , k = 0, 1, . . . , n h -1
and can be represented by the augmented M(s) -∆ interconnection shown in Figure 1. This LFT involves n h blocks ∆ packed in a p n h × p n h augmented uncertainty block ∆.

. . . switch . . . The integration over one period of the n h switched LTI systems allows approximating the

C(0) A(0) A(h) A((n h -1)h) . . . ∆ ∆ . . . ∆ ẋ(t) x(t) . dt B(0) C((n h -1)h) B((n h -1)h) B(h) C(h) . . . t r = remainder(t,T) t t r + + t r ∈ [(n h -1)h, T [ t r ∈ [0, h[ t r ∈ [h, 2h[ . . . + + + + ∆ M(s)
monodromy matrix R(∆) by R a (∆): R(∆) ≈ R a (∆) = n h -1 k=0 e Ap(kh,∆)h = e Ap((n h -1)h,∆)h . . . e Ap(h,∆)h e Ap(0,∆)h . ( 8
)
Thus computation of the monodromy matrix on a p-dimension parametric space gridding would be too CPU time-consuming to characterize the stability in the whole parameter space. Eq.( 8)

will be used to validate the lifting procedure proposed in the next section for several representative parametric congurations and for the worst-case parametric conguration provided by the µ-analysis on the lifted model (section III).

B. The lifting procedure

The lifting method is performed as follows [START_REF] Kim | Robustness Analysis of Linear Periodic Time-Varying Systems Subject to Structured Uncertainty[END_REF]. Initially, the n h continuous-time LTI systems

M(s, kh) (k = 0, 1, • • • , n h -1)
dened by Eq. ( 7) and represented in Figure 1 are discretized with a sampling period h, as discussed in section II C. The resulting discrete-time periodic system is:

M d (z, k) :        x d (k + 1) = A d (k) x d (k) + B d (k) w d (k) z d (k) = C d (k) x d (k) + D d (k) w d (k) (9)
with:

w d (k) = ∆ z d (k) . Matrices A d (k), B d (k), C d (k) and D d (k)
depend on the discretization method and are n hperiodic:

A d (k + n h ) = A d (k), B d (k + n h ) = B d (k), C d (k + n h ) = C d (k), D d (k + n h ) = D d (k) . ( 10
)
In the second step, the system ( 9) is integrated over one period n h . The nal discrete-time-invariant

LFT model M d (z) -∆ (lifted model) is: M d (z) :        x d (k + n h ) = A d x d (k) + B d V T w d (k) z d (k) = VC d x d (k) + VD d V T w d (k) (11)
with:

w d (k) = ∆ z d (k)
where: The third step uses an inverse Tustin transform to convert the system M d (z) (Eq. 11) back to the continuous-time domain in order to apply µ-analysis tools available only in continuous-time.

∆ = diag [δ 1 I n h , δ 2 I n h , . . . δ p I n h ]. Matrices (A d , B d , C d , D d ),
Since the µ-analysis is performed in the frequency domain, the Tustin transform (with a sampling period equal to T ) is selected here for its property of preserving the input-output frequency-domain response [START_REF] Franklin | Digital Control of Dynamic Systems[END_REF]. The nal continuous-time

-invariant LFT model M c (s) -∆ is M c (s) :        ẋc (t) = A c x c (t) + B c w c (t) z c (t) = C c x c (t) + D c w c (t) (12)
with: w c (t) = ∆ z c (t)

where:

• A c = 2 T (I n + A d ) -1 (A d -I n ) • B c = 2 T (I n + A d ) -1 B d V T • C c = 2VC d (I n + A d ) -1 • D c = VD d V T -VC d (I n + A d ) -1 B d V T
The original LTP system S(∆) in Eq.( 1) is now in the standard continuous-time LFT form M c (s) -∆, as given in Eq.( 12). In the next section, the lifting procedure [3] . is applied to the parametric robustness analysis of the ground resonance phenomenon.

C. Discretization methods

From the stability analysis accuracy point of view, errors in the evaluation of the perturbed monodromy matrix R a (∆) and its characteristic multipliers can only be introduced in the rst step (discretization). The approach was to re-use the discretization methods commonly used in the eld of automatic control [START_REF] Franklin | Digital Control of Dynamic Systems[END_REF]: (i) impulse invariance, (ii) zero-order hold, (iii) rst-order hold, (iv) Tustin transformation, (v) matched pole-zero. The impulse invariance and matched pole-zero methods are not considered since the impulse invariance method cannot handle systems with direct feed-through and the matched pole-zero method works only for single-input single-output systems. The three remaining methods are compared here from the angle of monodromy matrix R a (∆) approximation.

For k = 0, • • • , n h the closed-loop LFT representation (9) is:

x d (k + 1) = A d (k) + B d (k)∆(I n -D d (k)∆) -1 C d (k) x d (k) (13) = A d p (k, ∆) x d (k)
and by integration over one period n h , the closed-loop lifted LFT (11) can be expressed as:

x d (k + n h ) = n h -1 k=0 A d p (k, ∆) x d (k) . (14)
[3] The whole procedure for converting an LTP system into a continuous-time lifted system is embedded in a MATLAB R function ltp2lti.m which can be downloaded from http://personnel.isae.fr/daniel-alazard/ matlab-packages/lifting-procedure-for-linear-time.html. The package also contains a tutorial on the Mathieu equation.

The objective is to approximate the uncertain monodromy matrix for any value of ∆:

R a (∆) = n h -1 k=0 e Ap(kh,∆)h ≈ n h -1 k=0 A d p (k, ∆) . (15)
Therefore, for all k, matrices e Ap(kh,∆)h and A d p (k, ∆) must be compared for the various discretization methods. This problem of continuous-time LFT discretization was rst addressed in [START_REF] Imbert | Robustness Analysis of a Launcher Attitude Controller via µ-Analysis[END_REF]. The discretization error according to various methods is also discussed in-depth [START_REF] Tóth | On the Discretization of LPV state-space representations[END_REF][START_REF] Toth | On the Discretization of Linear Fractional Representations of LPV Systems[END_REF] in the more general framework of the Linear Parameter-Varying (LPV) system where the bounds on the discretization error are proposed in terms of approximating the state or output evolution. Here, we focus on the approximation error on the uncertain transition matrix e Ap(kh,∆)h .

A third order Taylor expansion in h of e Ap(kh,∆)h leads to (kh is omitted for brevity)

e Ap(∆)h ≈ I n + (A + ∆ A )h + (A + ∆ A ) 2 h 2 2 + (A + ∆ A ) 3 h 3 6 . ( 16
)
with • Zero order hold (zoh) method [START_REF] Hoagg | First-Order-Hold Sampling of Positive Real Systems and Subspace Identication of Positive Real Models[END_REF]: the input w(t) of system M(s, kh) is assumed to be constant over the oversampling period h:

∆ A = B∆(I n -D∆) -1 C. The expressions of A d (k), B d (k), C d (k), D d (k) from A(kh), B (
w(t) = w d (k), ∀ t ∈ [kh, (k + 1)h[. Then A d =e A h (17a) B d =A -1 e A h -I n B (17b) C d =C (17c) D d =D (17d)
associated with the state x d (k) = x(kh) and

A d p (∆) ≈ I n + (A + ∆ A )h + A(A + ∆ A ) h 2 2 + A 2 (A + ∆ A ) h 3 6 . (18)
• First order hold (foh) method: the input w(t) of system M(s, kh) is assumed to be linear between two consecutive over-samples:

w(t) = w d (k) + t -kh h (w d (k + 1) -w d (k)), ∀ t ∈ [kh, (k + 1)h[. Thus A d =e A h (19a) B d = 1 h A -1 e A h -I n B (19b) C d =CA -1 e A h -I n (19c) D d =D + 1 h CA -1 e A h -I n -A h A -1 B (19d) associated with state x d (k) = (e A h -I n ) -1 (Ax(kh) + Bw(kh)) -1 h A -1 Bw(kh) and A d p (∆) ≈ I n +(A+∆ A )h+(A+∆ A ) 2 h 2 2 +(A+∆ A ) 3 h 3 6 +(A+∆ A )∆ A (A+∆ A ) h 3 12 . (20)
• Tustin method [START_REF] Franklin | Digital Control of Dynamic Systems[END_REF]: the continuous-time integration presented in Figure 1 is approximated by a numerical integration using the trapezoidal rule:

x

((k + 1)h) = x(kh) + h 2 ( ẋ((k + 1)h) + ẋ(kh)) .
Then, with w d (k) = w(kh),

A d =(I n + h 2 A)(I n - h 2 A) -1 (21a) B d =h(I n - h 2 A) -1 B(kh) (21b) C d =C(I n - h 2 A) -1 (21c) D d =D + h 2 C(I n - h 2 A) -1 B (21d) associated with state x d (k) = I n -h 2 A(kh) x(kh) -h 2 B(kh)w(kh) and A d p (∆) ≈ I n + (A + ∆ A )h + (A + ∆ A ) h 2 2 + (A + ∆ A ) 3 h 3 4 . ( 22
)
Thus the approximation of e A(∆)h by A d p (∆) is only a rst order approximation if the zoh method is used, whereas it is a second order approximation with the foh and Tustin methods. Note that for the zoh and foh methods, A d , as dened in Eq.(A3a), is equal to the nominal monodromy matrix (R a (∆ = 0)) obtained from Eq.( 8) while approximations are made when the Tustin method is used. Therefore if ∆ A is assumed to be small, a better approximation with the foh method than the Tustin method can be expected. Indeed, the Taylor expansion ( 20) is equal to ( 16) except for a fourth order term in ∆ A h 3 /12. It is also important to underline the crucial inuence of n h on the accuracy of the result. High values of n h tend to minimize errors with all the methods, but also substantially increase robustness analysis computation time since the number of inputs and outputs of the initial LFT model is multiplied by n h in the new uncertainty block ∆.

This Taylor expansion based analysis allows us to recommend the foh method for the discretization of the LFT representation of uncertain systems instead of the more commonly used zoh method. This will be conrmed through the numerical results on the study of the ground resonance phenomenon in section III B.

III. Ground Resonance Parametric Analysis

A.

Ground resonance modeling

Figure 2 provides a general diagram of the dynamical system. It represents a simplied helicopter model similar to that used in the earliest research of the ground resonance phenomenon [START_REF] Coleman | Theory of Self-Excited Mechanical Oscillations of Helicopter Rotors with Hinged Blades[END_REF].

The fuselage is modeled as a rigid body with mass m f . x f (t) and y f (t) represent the fuselage positions along the longitudinal and lateral directions, respectively. The mechanical impedance between the fuselage and the ground (landing gear) is modeled by two stinesses K f X and K f Y and two damping factors C X and C Y acting in the longitudinal and lateral directions. At equilibrium, the fuselage center of mass (point O) coincides with the origin of the inertial reference frame The origin of the fuselage frame (x, y, z), parallel to the inertial frame, is located at the geometric center of the rotor hub (coincident at point O). The rotor angular velocity is denoted Ω.

(X 0 , Y 0 , Z 0 ).
The fuselage and rotor head are joined by a rigid shaft while aerodynamic forces on the blades are not taken into account. This assumption is quite realistic since the helicopter is on the ground.

Also note that aerodynamic eects can be embedded in the uncertainties on blade hinge damping The position of the k-th blade center of mass, written in the inertial reference frame, is given as:

x b k = a cos (ψ k ) + b cos (ψ k + ϕ k (t)) + x f (t) (23a) y b k = a sin (ψ k ) + bsin (ψ k + ϕ k (t)) + y f (t) (23b)
where a is the hinge oset and ψ k = Ωt + 2π(k-1)

N b , k = 1, • • • , N b .
The expressions of the kinetic energy, the potential energy and the work of dissipative forces are presented separately in Appendix B. By applying Lagrange equations and a rst order expansion of trigonometric terms, the linear dynamic model can be derived:

M q + G q + K q = F ext (24)
where q (t) = [ x f (t) y f (t) ϕ 1 (t) ϕ 2 (t) ϕ 3 (t) ϕ 4 (t) ] T is the generalized coordinates vector. M, G and K correspond to the mass, damping and stiness matrices respectively and are described in Eq.25. These matrices are not symmetric due to the presence of periodic terms. F ext is equal to zero if all the blades have the same inertial and geometrical properties.

M(t) =     1 0 -rm1 sin(ψ1) -rm2 sin(ψ2) -rm3 sin(ψ3) -rm4 sin(ψ4) 0 1 rm1 cos(ψ1) rm2 cos(ψ2) rm3 cos(ψ3) rm4 cos(ψ4) -r b1 sin(ψ1) r b1 cos(ψ1) 1 0 0 0 -r b2 sin(ψ2) r b2 cos(ψ2) 0 1 0 0 -r b3 sin(ψ3) r b3 cos(ψ3) 0 0 1 0 -r b4 sin(ψ4) r b4 cos(ψ4) 0 0 0 1     (25a) G(t) =    r cX 0 -2Ωrm1cos(ψ1) -2Ωrm2 cos(ψ2) -2Ωrm3 cos(ψ3) -2Ωrm4 cos(ψ4) 0 r cY -2Ωrm1sin(ψ1) -2Ωrm2 sin(ψ2) -2Ωrm3 sin(ψ3) -2Ωrm4 sin(ψ4) 0 0 rc1 0 0 0 0 0 0 rc2 0 0 0 0 0 0 rc3 0 0 0 0 0 0 rc4    (25b) K(t) =      ω 2 x 0 Ω 2 rm1 sin(ψ1) Ω 2 rm2 sin(ψ2) Ω 2 rm3 sin(ψ3) Ω 2 rm4 sin(ψ4) 0 ω 2 y -Ω 2 rm1 cos(ψ1) -Ω 2 rm2 cos(ψ2) -Ω 2 rm3 cos(ψ3) -Ω 2 rm4 cos(ψ4) 0 0 ω 2 b 1 +Ω 2 r 2 a 1 0 0 0 0 0 0 ω 2 b 2 +Ω 2 r 2 a 2 0 0 0 0 0 0 ω 2 b 3 +Ω 2 r 2 a 3 0 0 0 0 0 0 ω 2 b 4 +Ω 2 r 2 a 4      (25c) F ext (t) =       N b P k=1 Ω 2 r m k ( a+b a ) cos(ψ k ) N b P k=1 Ω 2 r m k ( a+b a ) sin(ψ k ) 0 0 0 0       (25d)
where:

r m k = b m b k m f + P N b k=1 m b k , r b k = b m b k b 2 m b k +Iz b k , r 2 a k = a r b k , r c X..Y = C X..Y m f + P N b k=1 m b k , ω 2 x = K f X m f + P N b k=1 m b k , ω 2 y = K f Y m f + P N b k=1 m b k , r c k = C b k b 2 m b k +Iz b k , ω 2 b k = K b k b 2 m b k +Iz b k , k = 1, • • • , N b
The aging or failure of mechanical elements comprising the helicopter rotor head, e.g., springs or dampers, has a direct inuence on the dynamical behavior of the whole system. Depending on the degradation of these elements, new critical rotating velocities may be reached at which the ground resonance phenomenon will occur. Thus the robustness analysis of helicopters under structured uncertainties is required to predict the smallest perturbation leading the system to instability. The

ω 2 b k = (1 + δ k ) ω 2 b k (k = 1, 2, 3, 4) (26)
where δ k corresponds to the relative uncertainties related to ω 2 b k (i.e.: square of the nominal blade resonance frequency).

According to Eq.( 1), the LTP model S(∆) of the ground resonance phenomenon taking into account uncertainties on the 4 blade hinge stinesses takes to form:

A(t) = O6×6 I6 -M -1 (t)K(t) -M -1 (t)G(t) , B(t) =      O6×4 -M -1 (t) 2
C(t) = O 4×2 diag[ω 2 b 1 , ω 2 b 2 , ω 2 b 3 , ω 2 b 4 ] O 4×6 , D(t) = O 4×4 (28)
∆ = diag[δ 1 , δ 2 , δ 3 , δ 4 ] (29)
associated with the state vector x = [q T qT ] T .

The model relative to the 4 blade hinge damping factors is described in section III C, equation [START_REF] Imbert | Robustness Analysis of a Launcher Attitude Controller via µ-Analysis[END_REF]. The numerical data are summarized in Table 1.

Table 1: Nominal parameters of a helicopter with identical blades Fuselage m f = 2902.9 • only the 4-th blade hinge stiness is considered to be uncertain (i.e.,∆ = diag[0, 0, 0, δ 4 ]),

• the rotor angular velocity is constant: Ω = 10π rad/s.

The discrete-time lifted model M d (z) in Eq.( 11) is constructed for three dierent values of n h (10, 30 and 100) and the three discretization methods (zoh, foh, Tustin). Then, for each value of the uncertainty δ 4 (from -100% to 100% by steps of 10%), the LFT M d (z) -∆ is resolved and compared with the Floquet monodromy matrix R a (diag[0, 0, 0, δ 4 ]) computed with n h = 100 (Eq. 2. From these results, it can be concluded that:

• for high values of n h (n h = 100), the stability analyses obtained with the three methods converge with the Floquet-based prediction,

• for low values of n h , the stability analysis based on the zoh method is poor,

• the best trade-o between stability-analysis accuracy and the reduction of n h is obtained with the foh method and n h = 30. This value will be adopted in the next section. only for very low values of δ 4 (δ 4 < -90%).

Remark III.1 The CPU-time for computing the lifted model [START_REF] Skjoldan | On the Similarity of the Coleman and Lyapunov-Floquet Transformations for Modal Analysis of Bladed Rotor Structures[END_REF] on a standard desk-top computer is 0.39 s, 1.98 s and 17.7 s, respectively for n h = 10, 30 and 100 (the eect of the discretization method on CPU-time is negligible). The CPU-time for computing the monodromy matrix [START_REF] Flowers | Nonlinear Rotorcraft Analysis Using Symbolic Manipulation[END_REF] with n h = 100 is 1.26 s. for dierent values of n h using the Tustin method in the lifting procedure, and |λ Ra |(δ 4 ).

C.

µ-analysis of ground resonance stability

In this section, all blade hinge stinesses are considered to be uncertain and independent. The µ-analysis toolbox [2] can be directly applied to the continuous-time lifted model M c (s) in Eq.( 12).

The structure of the uncertainty block ∆ 120×120 is therefore 4 real independent parameters repeated 30 times each. At each frequency ω, the µ-analysis computes an upper bound μ(ω) and a lower bound µ(ω) of the structured singular value µ. The µ-upper bound provides a guarantee of robust stability, i.e.

S(∆)

is stable ∀ δ i / |δ i | ≤ 1 max ω μ(ω) , i = 1, 2, 3, 4
while the µ-lower bound provides the worst parametric conguration ∆ worst (ω) [START_REF] Ferreres | A Practical Approach to Robustness Analysis with Aeronautical Applications[END_REF].

For the nominal rotor angular velocity Ω = 10π rad/s, the µ-upper and lower bounds provided by the Skew Mu Toolbox are plotted in Figure 6. It can be concluded that max ω μ(ω) = 12 (i.e.

the parametric robustness margin is 8.3 %) and max ω µ(ω) ≈ max ω μ(ω) (i.e. this margin is not at all conservative). The µ-analysis tools also provide the critical frequency ω c worst = 23.7 rad/s, the frequency of the instability that occurs when ∆ = ∆ worst . The parametric conguration at ω c worst [START_REF] Bir | Multi-Blade Coordinate Transformation and its Application to Wind Turbine Analysis[END_REF] Ferreres, G. and Biannic, J.-M., The Skew Mu Toolbox", http://www.onera.fr/staff-en/jean-marc-biannic/.

is ∆ worst = diag[0.085, 0.085, 0.085, 0.085] .

It is now possible to compute the monodromy matrix R a (∆ worst ) (Eq.15) to validate the µ-analysis result and the representativeness of the lifted model. The 12 eigenvalues (characteristic multipliers) of R a (∆ worst ) associated with the 6 monodromic modes are given in Table 3. Mode # 2 is unstable with a magnitude very close to one. Although rotor stability is quite robust regarding uncertainty on a single blade hinge stiness (see section III B), robustness to uncertainties on all four blades is quite poor. Our analysis showed that the worst case conguration corresponds to a rotor with identical blades. An interesting observation regarding this analysis is that no dissimilar blades congurations are worse than ∆ worst from a stability point of view. This point is conrmed by further analysis considering uncertainties on 2 adjacent blades, 2 opposing blades and 3 blades (see Table 4 for a summary of µ-analysis results).

Once again, considering Figure 4

This analysis conrms that the parametric robustness margin is better for a rotor with dissimilar blade hinges. Thus to improve ground resonance stability it is possible to imagine a mechanism mounted on a single blade hinge and only operated on the ground in order to create asymmetry in the rotor's properties. Of course, this analyze concerns only the ground resonance phenomenon and any conclusion regarding the advantages of a rotor with dissimilar blade hinges cannot be extended to behavior during ight.

Other analyses were performed for various rotor angular rates Ω (from 1 to 10 Hz). In all cases, the worst case conguration corresponds to a symmetric conguration (i.e.: δ 1 = δ 2 = δ 3 = δ 4 = δ).

Assuming that δ is repeated for the 4 blades, the stability analysis can be simplied (1 uncertain parameter repeated 120 times) and it is possible to plot the stability domain according to δ for various rotor angular rates Ω. This stability domain in the Ω -δ plane is shown in Figure 7 and conrms that the stability margin of the nominal conguration (Ω = 10π (rad/s), δ = 0) is very weak. It is also noteworthy that if the stability analysis is restricted to symmetric uncertainties, the discrete-time lifted model M d (z) can be used directly to plot the evolution of characteristic multipliers as a function of δ in the z-plane using a basic root locus [START_REF] Franklin | Digital Control of Dynamic Systems[END_REF]. This plot is presented in A possible physical justication of the worst-case conguration can be formulated: for a rotor with identical blade hinges, it is easy to check in Table 3 that modes # 5 and 6 are identical. This multiple-monodromic mode may be more signicant in the kinematic energy exchange between the rotor and the fuselage than single modes. In other words, the dissymmetry in the rotor breaks this multiple mode and leads to 6 single modes associated with 6 dierent frequencies and lower modal participation factors. This result must be conrmed by further analyses and confronted with an arbitrary number N b of blades (an odd number for instance).

Remark III. [START_REF] Bir | Multi-Blade Coordinate Transformation and its Application to Wind Turbine Analysis[END_REF] The CPU-time to compute the lifted model (with n h = 30) and µ-upper bound is

Once again, the worst case corresponds to a symmetric reduction of all damping factors.

IV. Conclusions

A procedure that combines Floquet theory and the Lifting technique to convert the uncertain linear time periodic (LTP) system into an uncertain linear time-invariant (LTI) lifted system makes it possible to use µ-analysis to determine the smallest perturbation in helicopter blade properties having the worst case eect on rotor stability.

Three discretization methods were compared in order to reduce the complexity (and hence the CPU time) of the lifted model and minimize discretization errors. The conclusion is that the rst order hold method gives signicantly better results than the more commonly used zero-order hold method.

By taking into account uncertainties on the four blade hinge stinesses and damping factors, the µ-analysis shows that the worst case parametric conguration corresponds to a symmetric perturbation for all blade hinges. Thus stability analysis methods restricted to symmetric rotors like the Coleman method are still relevant. Nevertheless, tools are now available to perform further analyses and verify whether this result can be extended to other rotor congurations (for example a rotor system with an odd number of blades). Such tools save signicant CPU-time in comparison with a pure Floquet analysis combined with a parametric space exploration.

V.

NomenclatureA

  (t), B(t), C(t), D(t) = continuous-time state-space matrices A c , B c , C c , D c = state-space matrices of the continuous-time lifted model A d (k), B d (k), C d (k), D d (k) = discrete-time state-space matrices A d , B d , C d , D d = state-space matrices of the discrete-time lifted model a = rotor eccentricity (m) b = blade equivalent length (m) C X , C Y = fuselage damping factors in x and y directions (Ns/m) C b k = k-th blade hinge damping factor (Nms/rad) F ext = normalized external force vector h = oversampling period (s) O n×m = n × m null matrix I n = n × n identity matrixI Zb k = k-th blade lag inertia (kg m 2 ) K b k = k-th blade hinge stiness (Nm rad -1 ) K f X , K f Y = fuselagelongitudinal and lateral stinesses (N m -1 ) M, K, D = mass, stiness and damping matrices m f , m b k = fuselage and k-th blade masses (kg) n = system order (integer) n h = number of over-samples in one period (integer) N b = number of blades (integer) p = number of uncertain parameters q = vector of degrees of freedom r a k = √ a r b k r b k = k-th blade static moment to total inertia ratio (m -1 ) r c k = k-th blade damping factor to total inertia ratio (s -1 ) r c X , r c Y = fuselage damping factor to total mass ratios, (s -1 ) r m k = k-th blade static moment to total mass ratio (m) R = monodromy matrix t = time (s) T = period of the Linear-Time-Periodic (LTP) system (s) V = permutation matrix x,w,z = state, input and output vectors x (t) , y (t) = fuselage longitudinal and lateral position (m) x bk , y bk = k-th blade positions (m) (x, y, z) = frame attached to the rotor hub (X 0 , Y 0 , Z 0 ) = inertial frame ∆ = uncertainty matrix µ = structured singular value ϕ k (t) = k -th blade lead-lag angle (rad) Ω = rotor angular velocity (rad s -1 ) ζ k = k-th blade azimuth angle (rad) ω b k = k-th blade cantilevered frequency (rad s -1 ) ω x , ω y = cantilevered fuselage frequencies (rad s -1 )

Figure 1 :

 1 Figure 1: Sketch of the n h LTI models switched over one period T in positive feedback with the augmented uncertainty block ∆ (direct feed-through matrices are omitted for legibility).

  the re-ordering matrix V, augmented input and output vectors w d (k) and z d (k) are detailed in Appendix A.

  kh), C(kh), D(kh) and h, the mapping between discrete-time state x d (k) and continuous-time state x(kh) and input w(kh), and the third order expansion of A d (k, ∆) are described below for the three discretization methods (k and kh are omitted for brevity):

  The rotor head system is comprised of one rigid rotor hub and an assembly of N b blades. The k-th blade has a mass m b k , a moment of inertia I zb k around the z -axis located at its center of mass and an in-plane lead-lag motion dened by ϕ k (t). The radius of gyration is dened by the length b. Angular spring and viscous damping are considered on each blade hinge (point B). Spring stiness and viscous damping coecients are denoted K b k and C b k , respectively.

Figure 2 :

 2 Figure 2: Diagram of the Mechanical System

  uncertainties introduced on the mechanical model are related to blade hinge stinesses K b k and damping factors C b k . The 4 blade hinge stinesses are normalized with respect to the in-plane lead-lag cantilevered frequency squared:

  [kg] ωx = 6.0 π [rad/s] ωy = 8.0 π [rad/s] CX = 5.71 10 3 [Ns/m] CY = 7.62 10 3 [Ns/m] Rotor a = 0.2 [m] b = 2.5 [m] m b k = 31.9 [kg] ω b k = 3.0 π [rad/s] I zb 1..4 = 259 [kg m 2 ] C b 1..4 = 432 [Nms/rad] B. Validation of the lifting procedure The objectives of this section are: (i) to validate the lifting procedure presented in section II B by comparison with Floquet analyses for various congurations of the uncertain parameters, and (ii) to select the best discretization method while minimizing the over-sample number n h in order to obtain a good trade-o between analysis accuracy and CPU-time reduction. The following assumptions are made:

8 ).

 8 The comparison index is the highest eigenvalue (or characteristic multiplier) magnitude denoted |λ l |(δ 4 ) and |λ Ra |(δ 4 ) for the lifted model and the monodromy matrix, respectively. Note that the nominal angular velocity (Ω = 10π rad/) leads to a weak stability margin (in terms of characteristic multiplier magnitude), which is quite sensitive to parametric uncertainties. Indeed, the maximal characteristic multiplier magnitude of R a (O 4×4 ) is |λ Ra |(0) = 0.982. The results for the three discretization methods are presented in Figures 3 to 5 and summarized in Table

Figure 3 :

 3 Figure 3: Evolution of the magnitude of the highest characteristic multiplier with respect to δ 4 : |λ l |(δ 4 ), for dierent values of n h using zoh method in the lifting procedure, and |λ Ra |(δ 4 ).

Figure 4 :

 4 Figure 4: Evolution of the highest magnitude of the characteristic multiplier with respect to δ 4 : |λ l |(δ 4 ), for dierent values of n h using foh method in the lifting procedure, and |λ Ra |(δ 4 ).

Figure 5 :

 5 Figure 5: Evolution of the highest characteristic multiplier magnitude with respect to δ 4 : |λ l |(δ 4 ),

where ω 2 b 4 ← 1 .5 ω 2 b 4

 4124 , the best stability margin in term of the highest characteristic multiplier magnitude is obtained for δ 4 = 0.5. The stability robustness analysis for an asymmetric rotor leads to the following results: max ω μ(ω) = 8.11, ∆ worst = diag[0.1312, 0.1312, 0.1312, -0.1312] .

Figure 8

 8 Figure 8 in the case Ω = 10π (rad/s). Instability occurs for δ = 0.085 at frequency

Figure 6 :

 6 Figure 6: µ upper bound (line) and lower bound (• marks) -robustness stability analysis with respect to uncertainties on the 4 blade hinge stinesses at Ω = 5 Hz.

  Analyses were also performed to evaluate robustness with respect to the hinge damping factors(parameters r c k , k = 1, • • • , N b in Eq. (25a)) which are certainly also very sensitive to aging eects.

Figure 7 :

 7 Figure 7: Rotor stability domain (dotted area) according to angular rate Ω and an identical uncertainty δ on the 4 blade hinge stinesses.

Figure 8 :

 8 Figure 8: Characteristic multiplier loci in the z-plane according to an identical uncertainty δ > 0 (black) or δ < 0 (grey) on the 4 blade hinge stinesses at Ω = 5 Hz.

  The continuoustime LTP model is transformed into a discrete-time LTI model and the size of the uncertainty block is increased according to the truncation order. Additional uncertainties are taken into account to handle this truncation error but introduce conservatism. In [1418], symbolic methods (e.g., the multiple scale and harmonic balance methods) are used to analyze helicopter ground resonance by considering mechanical nonlinearities. However, the inuence of dissimilarities from one blade to

Table 2 :

 2 Maximal relative error (%) on the highest characteristic multiplier magnitude with respect to the nominal stability margin: max δ4 |λ l |(δ 4 ) -|λ Ra |(δ 4 ) /(1 -|λ Ra |(0)).

	n h	10	30	100
	zoh	323	115	35
	foh	14.5	1.43	0.49
	Tustin	125	14.5	1.49
	It can also be concluded that rotor stability is quite robust with respect to variations of single
	blade stiness. Instability occurs (i.e. highest characteristic multiplier magnitude greater than 1)

Table 3 :

 3 Characteristic multipliers of matrix R a (∆ worst ).

	Mode Characteristic multipliers Magnitude
	# 1	0.1836 ± 0.7602 j	0.7821
	# 2	-0.6981 ± 0.7164 j	1.0003
	# 3	-0.5907 ± 0.4742 j	0.7575
	# 4	-0.7088 ± 0.5466 j	0.8951
	# 5	-0.5962 ± 0.6876 j	0.9101
	# 6	-0.5962 ± 0.6876 j	0.9101

Table 4 :

 4 µ-analysis results for uncertainties on 1 blade, 2 adjacent blades, 2 opposite blades, 3 blades or 4 blades at Ω = 5 Hz.

	Case	maxω μ(ω)	∆worst
	1 blade	1.1	-0.9
	2 adjacent blades	1.41	diag[-0.71, -0.71]
	2 opposite blades	1.31	diag[-0.78, -0.78]
	3 blades	6.2	diag[0.17, 0.17, 0.17]
	4 blades	12	diag[0.085, 0.085, 0.085, 0.085]

Appendix A: LIFTING PROCEDURE

The lifting procedure, as described in [START_REF] Kim | Robustness Analysis of Linear Periodic Time-Varying Systems Subject to Structured Uncertainty[END_REF][START_REF] Colaneri | Constant-Coecient Representations of Periodic-Coecient Discrete Linear Systems[END_REF], is summarized here for the reader's convenience.

The state-space matrices of the discrete-time lifted model are expressed directly from the state-space matrices provided by the selected discretization method (Eqs. [START_REF] Flowers | Chaotic Dynamical Behavior in a Simplied Rotor Blade Lag Model[END_REF] to [START_REF] Verdult | Identication of Linear Parameter-Varying State-Space Models with Application to Helicopter Rotor Dynamics[END_REF]) and do not involve the transition matrix.

From Eq.( 9), let us dene the packed output and input vectors

Then, the integration over one period of system in Eq.( 9) leads to the discrete-time lifted system

where

Φ u l =I n otherwise.

(A3c)

The general expression for the lower triangular terms of matrix D d is

Finally, since each δ i for i = 1, 2, . . . , p appears in each ∆ block, the row re-ordering matrix can be dened as V so that

Let us dene the re-ordered packed input and output vectors

then the LFT representation (A2) is transformed into the LFT representation [START_REF] Skjoldan | On the Similarity of the Coleman and Lyapunov-Floquet Transformations for Modal Analysis of Bladed Rotor Structures[END_REF].

Appendix B: KINETIC AND POTENTIAL ENERGY AND WORK OF DISSIPATIVE FORCE

The kinetic and potential energy expressions and the work expression of dissipative forces of the dynamical system are presented here for the reader's convenience [START_REF] Coleman | Theory of Self-Excited Mechanical Oscillations of Helicopter Rotors with Hinged Blades[END_REF][START_REF] Hammond | An Application of Floquet Theory to Prediction of Mechanical Instability[END_REF]. They are all written in the inertial reference frame.

• The kinetic energy of the whole dynamical system consists of the sum of kinetic energy expression of the fuselage T F us and rotor head T RH systems:

2b(Ω ẏ+ ẏ φk )cos(ψ k +ϕ k )-2b(Ω ẋ+ ẋ φk )sin(ψ k +ϕ k )

• The potential energy of the whole dynamical system consists of the sum of the potential energy of the fuselage, U F us , and rotor head, U RH , systems:

• The work of dissipative forces of the whole dynamical system consists of the sum of the work done by dissipative forces acting on the fuselage δF F us and the rotor head δF RH systems:


Sihem Cherrared

Sofiane Imadali

Eric Fabre

Gregor Gössler

LUMEN: A Global Fault Management Framework For Network Virtualization Environments

Keywords: Fault management, Network Virtualization Environment, Self-diagnosis

The advent of 5G and its ever increasing stringent requirements for bandwidth, latency, and quality of service pushes the boundaries of what is feasible with legacy Mobile Network Operators' technologies. Network Function Virtualization (NFV) is one promising attempt at solving some of those challenges that were widely adopted by the industry and the standardization bodies. At its essence, NFV is about running network functions as software workloads on commodity hardware to optimize deployment costs and simplify the life-cycle management of network functions. However, it introduces new fault management challenges including dynamic topology, multi-tenant fault isolation and data consistency and ambiguity; that we propose to define in this paper. To tackle those challenges, we extend the classical fault management process to the virtualized functions by introducing LUMEN: a Global Fault Management Framework. Our approach aims at providing the availability and reliability of the virtualized 5G end-to-end service chain. LUMEN includes the canonical steps of the fault management process and proposes a monitoring solution for all types of Network virtualization Environments. Our framework is based on open source solutions and could easily be integrated with other existing autonomic management models.

Introduction

The telecommunication industry is facing increased competition as new players are entering with emerging software technologies and open source projects. One such project, Telecom Infra Project [1] was initiated by Facebook and oriented towards an open source and general purpose hardware for a new generation of Mobile Network Operators (MNOs). The current MNOs infrastructures consist of a large range of proprietary hardware appliances [START_REF] Han | Network function virtualization: Challenges and opportunities for innovations[END_REF]: there is thus a need for investment and adoption of new open technologies such as Software Defined Networking (SDN), Network Function Virtualization (NFV), and cloud-native innovations. Therefore, MNOs are investing a lot of resources considering the softwarization of network functions and are joining open source software communities. Most notably, MNOs contribute to OpenStack for infrastructure, OpenDayLight for SDN controllers, and ONAP for network automation [START_REF]Openstack USER SURVEY[END_REF] [START_REF]Open network automation platform[END_REF].

The current impulse to focus MNOs efforts on softwarization comes from the advent of the next mobile network generation (5G). The objective is the projected benefits the MNOs expects from such a paradigm shift when compared to how traditional networks are operated. The benefits include the reduction of costs (management and deployment), optimization to maximum resource Figure 1: The Fault Management Steps utilization, on demand programmable services, infrastructure and services sharing, and allowing a large variety of ecosystems to co-exist including open source systems [START_REF] Manzalini | Towards 5g software-defined ecosystems[END_REF].

Network management is critical for service availability, reliability and resiliency of networks. It includes the X-management axes or FCAPS : Fault, Configuration, Accounting, Performance, and Security. In the current SDN and NFV context, in addition to the possibility of multi-tenancy and slicing in 5G [START_REF] El Hattachi Rachid | [END_REF], new fault management issues and vulnerabilities are expected. The fault management issues in Network Virtualization Environments (NVEs) address scalability of network functions, lack of network visibility, and dynamic network topologies and their impact. The focus of our paper is to tackle those fault management issues in NVEs by introducing a Global Fault Management Framework (LUMEN).

LUMEN includes the canonical fault management steps and provides a new vision for monitoring all types of NVEs (e.g. SDN, NFV, monolithic or distributed functions, single and multitenant environments). One of LUMEN's objectives is to provide availability and reliability for VNFs deployed end-to-end service chains (NAT, firewalling, intrusion detection, DNS, caching and more). Our framework could also be integrated with other fault management Self-X models: Self-modeling, Self-diagnosis and Self-healing. Our contributions in this paper include:

• A definition of the different fault management steps and a classification of main approaches in the state of the art;

• A description of multi-tenant Network Virtualization Environments and the related challenges;

• Our proposal, LUMEN, to tackle the reliability, availibility, and resiliency of VNF chains in multi-tenant NVE;

• A discussion of LUMEN's advantages over related solutions.

The remainder of this paper is as follows. Section 2 discusses related work addressing fault management in traditional and softwarized networks. Section 3 summarizes the fault management issues and challenges for NVEs. Section 4 presents the different steps of LUMEN framework. Section 5 discusses the advantages of using our solution. Section 6 concludes the paper with some perspectives.

Related work

In this Section, we provide a definition of the classical fault management steps and techniques, the recent approaches in NVEs and the efforts done towards the automation of the fault management procedures and their architecture.

Fault Management

Fault management covers the detection and storage of events and alarm notifications, filtering procedures for these alarms, and diagnostic checks for root cause localization, impact analysis and corrective actions. Events can be related to metrics and Key Performance Indicators (KPIs), such as delay, jitter, or response times. Notifications can be SNMP traps [START_REF] Jeffrey D Case | Simple network management protocol[END_REF] or system log-files entries using syslog protocol [START_REF] Gerhards | The Syslog Protocol[END_REF]. The Fault-diagnosis process consists of the detection of the faulty state or the failure, the localization of faulty entities in the network and the identification of the fault types [START_REF] Ma Lgorzata Steinder | A survey of fault localization techniques in computer networks[END_REF][START_REF] Feng | A survey of intelligent network fault diagnosis technology[END_REF][START_REF] Zaytoon | Overview of fault diagnosis methods for discrete event systems[END_REF].

We propose the taxonomy in Figure 1 to help classify some of the main network fault management techniques. In the fault detection step the main approaches use alarms or metrics and once the faulty state is detected comes the fault localization and root cause analysis methods. We classify the state of the art into two main branches: white and black box techniques. As a definition, white box are techniques that give a clear view about the building process such as the dependency graphs methods, while in the black box techniques, it is harder to deduce from the network what exactly has the framework learned. Neural Networks for instance is in the latter case.

The two approaches can be combined. White box models such as Bayesian networks [START_REF] Tembo | A 3-layered self-reconfigurable generic model for self-diagnosis of telecommunication networks[END_REF][START_REF] Hounkonnou | Empowering self-diagnosis with self-modeling[END_REF][START_REF] Bennacer | Optimization of fault diagnosis based on the combination of bayesian networks and case-based reasoning[END_REF][START_REF] Steinder | Probabilistic fault localization in communication systems using belief networks[END_REF][START_REF] Sanchez | Self-modeling based diagnosis of softwaredefined networks[END_REF] use the network dependency graph to pinpoint faulty components and black box Artificial Intelligence (AI) techniques can be applied to detect a faulty state [START_REF] Zhaojun | Statistic and analysis for host-based syslog[END_REF][START_REF] Zheng | System log pre-processing to improve failure prediction[END_REF][START_REF] Kimura | Proactive failure detection learning generation patterns of large-scale network logs[END_REF]. The second phase is the recovery process. In the recovery process, the first step consists of identifying the impact of faults i.e. the affected entities and the more critical faults to be treated first. After that, the framework can proceed to the healing phase. These approaches usually rely on the expertise of network administrators and their knowledge of the underlying infrastructure.

Fault Management in Network Virtualization Environments

Autonomic management is a topic of interest for many organization (ETSI, IETF and SON) and is proposed in the design of the 5G architecture management and orchestration [START_REF] Nfvgs Etsi | Network functions virtualisation (nfv); management and orchestration[END_REF]. Recent efforts [START_REF] Hounkonnou | Empowering self-diagnosis with self-modeling[END_REF][START_REF] Tembo | A 3-layered self-reconfigurable generic model for self-diagnosis of telecommunication networks[END_REF][START_REF] Sanchez | Self-modeling based diagnosis of softwaredefined networks[END_REF][START_REF] Sauvanaud | Anomaly detection and root cause localization in virtual network functions[END_REF][START_REF] Soualah | A link failure recovery algorithm for virtual network function chaining[END_REF], introduced automation in the fault management steps. For instance, Hounkonnou et al [START_REF] Hounkonnou | Empowering self-diagnosis with self-modeling[END_REF] proposed a self-modeling concept to face the problem of scalability of the IP Multimedia Sub-system (IMS). The authors use instances of generic Bayesian networks adapted to the actual topology. Still, the experts knowledge on the IMS architecture helps to extract the generic model, which is not possible in dynamic networks based on NFV where the topology changes by integrating new tenants on the shared infrastructure. The automation of the steps was also introduced in some recent research papers that address SDN and NFV environments with self-diagnosis methods using Bayesian networks for certain efforts [START_REF] Sanchez | Self-modeling based diagnosis of softwaredefined networks[END_REF]. While others [START_REF] Sauvanaud | Anomaly detection and root cause localization in virtual network functions[END_REF][START_REF] Soualah | A link failure recovery algorithm for virtual network function chaining[END_REF], propose Self-healing for NFV chain failures.

Figure 2 summarizes the different fault management steps and presents a new fault management outline with the Self-X concepts.

In the telecommunication industry, the ITU-T's M30xx recommendation series describe a general framework for Telecommunication Network Management [START_REF]Principles for a telecommunications management network[END_REF]. The framework, while applied to certain deployments, has been enriched to accommodate new challenges and overcome the general inflexibility of their hierarchical setups [START_REF] Varga | Integration of service-level monitoring with fault management for end-to-end multi-provider ethernet services[END_REF]. With regards to automation of fault management and the design of fault management architectures, several approaches can be found in the literature. Sánchez et al. [START_REF] Manuel | Softwarized 5G Networks Resiliency with Self-Healing[END_REF], presented a self-healing framework for the resiliency of 5G networks. The proposed Self-healing framework acts in the three planes of SDN (management, control, and data) and in the service plane, by taking observations from the network and launching recovery actions. The Self-healing framework focuses more in the autonomic recovery actions and the SDN environment issues. Abhishek et al. [START_REF] Abhishek | Buddi: Bug detection, debugging, and isolation middlebox for software-defined network controllers[END_REF], presented a bug detection, debugging, and isolation (BuDDI) middlebox architecture for SDN controllers. The proposed solution focuses on software bugs, it includes a common mode failure and dependencies (CMFD) migration module that prevents the system to recur into the same software failures after a fail over.

Authors of [START_REF] Varga | Integration of service-level monitoring with fault management for end-to-end multi-provider ethernet services[END_REF] proposed a fault management framework (ISF) integrated to a service-level monitoring for an end-to-end Ethernet service. The framework handles different management modules. The root cause analysis module is based on Petri nets. While the architecture is general enough to accept heterogeneous use cases, authors of [START_REF] Manuel | Softwarized 5G Networks Resiliency with Self-Healing[END_REF] limited their case studies to network faults and Service Level Agreement (SLA) violation on the network, regardless of the network functions and their operational conditions that might trigger them.

The Open Network Function Virtualization (OPNFV) is an open source project (hosted by the Linux Foundation) for facilitating the development and evolution of NFV components across various open source ecosystems with automated testing to enable accelerated NFV development and reproducible deployments [START_REF]The Linux Foundation Projects[END_REF]. The OPNFV project Doctor proposes Vitrage, a project proposed by NTT DOCOMO and developed by the OPNFV DOCTOR community in collaboration with Nokia [START_REF]Openstack Vitrage[END_REF]. Vitrage is an OpenStack root cause analysis tool that provides rapid alarm notifications using deduced alarms for virtual and physical entities. Vitrage gets its data source from OpenStack modules and external opensource tools such as Nagios and Zabbix. A real time topology mapping is provided. The diagnosis process is based on templates to express the different rules used for the Root Cause Analysis (RCA). The templates represent static detection and alarming scenarios written by the infrastructure administrators and network functions experts.

Previously mentioned tools and research prototypes do not address the fault management of the whole NVEs including NFV deployments (monolithic and distributed) in a multi-tenant environments.

NVEs challenges and issues

NVEs enables distinct network architectures to coexist in a single infrastructure without affecting the network performance [START_REF] Esteves | On the management of virtual networks[END_REF]. NVEs include interactions between NFV and SDN components, physical and virtual functions in multi-tenants and multi-services environment. The network infrastructure connecting Virtual Network Functions (VNFs) into a chain will be handled by the SDN controller while each SDN application (e.g. Virtual Switch) will itself represent a VNF. However, a new architecture implies new challenges and new ways to perform management. In the following, we discuss the main issues of NVEs with a fault management viewpoint. Most of the challenges described in this section are addressed by our framework in Section 4.

• Number and type of managed objects: NVEs features mentioned above enable deployment of more services involving more entities to be managed. The managed entities have different granularities: Logical resources (e.g. container, Virtual Machine (VM) or libOs.), physical and virtual network functions (e.g. Physical Switch and Virtual Switch) and subcomponents (e.g. CPU and network interfaces). Moreover, the number of faults is expected to augment due to the large amount of services and virtual entities. Network administrators will face the problem of huge alarms quantities, alarm loss, delay and consistency. Therefore, managing logs and metrics without automated solutions became impossible.

• Dynamic network topology: 5G will enable deployment of real-time services tailored to customers' requests. The service is represented as a VNF Forwarding Graph (VNF-FG) [START_REF] Gsnfv Etsi | Network functions virtualisation (nfv); use cases[END_REF][START_REF] Haeffner | Service function chaining use cases in mobile networks[END_REF], which is a chain of connected VNFs. Each VNF deploys a specific function on the service. A good example of such chain is the Clearwater project for IMS [START_REF]Project clearwater[END_REF]. Other examples are Firewalls, Traffic Detection Function (TDF), or Traffic Steering Support Function. This real time service deployment makes the evolution of network topology and the entities dependencies in the network unpredictable. The managed system should consider the real time topology changes to pinpoint the faulty component and avoid false or outdated results.

• Lack of network visibility: The distribution of resources and virtual functions in different infrastructures and distinct locations leads to a network visibility problem. Figure 3 depicts how a unique VNF-FG service is distributed in different locations. Moreover, The state of each VNF is related to the state of hardware servers where the VNF is running. Therefore, to keep the service availability, a global network view is necessary.

• Multi-tenancy and fault isolation: MNOs benefits from sharing their infrastructure by maximizing resource usage. Multi-tenancy consists of sharing an infrastructure among multiple tenants or clients, each with a subset view of the owner's infrastructure resources (called a slice). As a consequence, fault isolation problems may occur [START_REF] Del Piccolo | A survey of network isolation solutions for multi-tenant data centers[END_REF]. Moreover, to allow a rapid notification of tenants faults, the infrastructure owners should identify the affected tenant entities to enable tenants to perform recovery actions. Figure 3 illustrates the multi-tenant fault isolation issue. In this example, a hardware server crash affects two VNFs (i.e. VNF1 and 2) of tenant 1 and 3 respectively. In this case, a rapid isolation of faults and notification of the two tenants is crucial to provide necessary recovery actions.

• Ambiguity and consistency of data: In the detection and localization process, events in the form of logs or metrics are collected to identify the failures. These events contain important information about the health and operations of the system. However, the collected data originate from distinct sources with different formats and are most of the time ambiguous and full of insignificant information for the diagnosis of faults. Therefore, an efficient storing and extraction methods are necessary to prepare data before starting the localization process.

LUMEN Framework

In this section, we discuss the scope of our proposed LUMEN framework and describe the functions of each plane.

Overview of the LUMEN Framework

In order to address fault management NVE challenges, we propose LUMEN: a Global Fault Management Framework. LUMEN (Figure 4) is a four step architecture where each plane summarizes the methods that should be deployed to address the different NVEs challenges of Table 1.

The LUMEN framework leverages an open source tool, namely the Elastic Stack [START_REF]Elastic stack[END_REF], to address fault management challenges. The LUMEN framework can also integrate different tools at each layer if necessary. The Elastic Stack ensures a real time data collection, storage, search, analysis, and visualization. The Elastic Stack allows the centralization of logs collected from different locations. This feature is well suited for a multi-tenant environment, where resources of the same tenant can be located on different infrastructures. The Elastic Stack is composed of three tools: Logstash/Beats, Elasticsearch and Kibana. Beats and Logstash are used for collecting data. Beats are lightweight data shippers that can be installed as dedicated agents on managed entities to send specific types of operational data. Logstash allows a more large data collection and enables filtering, enriching, and transforming data from a variety of sources. Elasticsearch is a search engine and analytics NoSQL database designed for storing efficiently the gathered data. Finally, Kibana is used for the data extraction and visualization.

The way LUMEN was designed shows the importance of data in software networks [START_REF] Cui | When big data meets software-defined networking: Sdn for big data and big data for sdn[END_REF]. Information about the network health and the topology changes are gathered in LUMEN. However, this data is frequently mixed with noise from insignificant messages and unrelated events for the fault diagnosis process. For instance, only syslog messages that have severity from 0 to 4 (0: emergency, 1: alert, 2: critical, 3: error and 4: warning) are considered for the RCA process. LUMEN is a four-layer framework that presents efficient techniques to organize and prepare data for the deductions methods.

LUMEN Planes

The functions of the LUMEN framework planes are described in the following.

Source Plane: The first step of our framework consists of gathering all types of data from distinct entities and distributed locations. As an example, the Beats agents of the Elastic Stack can be deployed in every network entity to send real time alarms and metrics to Elasticsearch. The gathered data depend on the deduction method and the process that will be deployed in the decision plane. For instance, one way to answer the problem of dynamic topology for the model-based technique is to model the network entities dependencies using the real time network topology information extracted from data. The network topology information can be found in SDN controllers or Virtual Infrastructure Manager (VIMs) network modules like OpenStack Neutron [START_REF]Openstack neutron[END_REF]. Monitoring data such as logs and alarms are also an important source of inputs. They generally contain relevant hidden information about the network state, faults and root causes.

The collected data can originate from different management levels with distinct types and formats and depend on:

• the management entity: if the collected data concern a specific tenant or the whole infrastructure.

• the granularity: if the collected data concern virtual or physical entities, VMs or containers or sub-components (VM CPU or Network interface).

• the type of data: if the collected data is syslog alarms, metrics or topology information; application or configuration logs.

• the decision method: if the collected data is for RCA or SLA or fault prediction.

Sink Plane: Since the collected data originate from different sources, entities and technologies, the data will have different formats (e.g. syslog , JSON or CSV). Therefore, a unification of the data format and the organization of data is important to facilitate the next steps of fault management. Moreover, the data should be filtered and the insignificant messages should be dropped.

Figure 5-D presents an example of a controller syslog message structured by Logstash in a JSON format. The example illustrates how Logstash can unify messages from disparate sources and normalize the structure of the message before storing it in the Sink. While Figure 5-A and 5-B exposes how Filebeat configuration files enable message filtering and dropping irrelevant information, respectively. Adding supplementary indications before the storing process is also crucial to localize the tenants and the network slices in further fault management steps. Filebeat enables transforming the data with adding new fields (Figure 5-C).

After the data transformation process is completed, the collected data can be stored in a Sink, Elasticsearch in our case. Elasticsearch is a robust search engine to centralize data and enable efficient extraction. Elasticsearch can also be distributed in different infrastructures.

Extraction Plane: This step is highly dependent of the decision process. In fact, the extraction of data depends on the inference engine used in the diagnosis process and the decision methods. For instance, some KPIs (VM cpu or disk load) can be mined to calculate VNF SLAs violation. Application Programming Interfaces (APIs) of the Sink plane can be used for the data extraction process.

Decision Plane: The last step of LUMEN is the decision plane which is the most important step in our fault management framework. In this step we create an educated guess through one or more deduction methods and diagnosis approaches such as Bayesian network and machine learning techniques. This step is enabled by the first three planes that provide the necessary information in a unified and organized way for a real time, rapid and efficient decision process.

Discussion and Analysis

We here discuss the use of the LUMEN, ISF [START_REF] Varga | Integration of service-level monitoring with fault management for end-to-end multi-provider ethernet services[END_REF], and Vitrage OpenStack [START_REF]Openstack Vitrage[END_REF] frameworks applied to three case studies and give a qualitative analysis.

LUMEN enables an autonomic collection and organization of data to prepare the fault management decision process. Plane 1, 2 and 3 represent the detection step in fault management, while Plane 4 represents the localization and recovery step. Further, Self-X frameworks such as the Self-healing framework in [START_REF] Manuel | Softwarized 5G Networks Resiliency with Self-Healing[END_REF] can be integrated in LUMEN. Table 1 presents how LUMEN planes addresses the issues of NVEs discussed in Section 3. Figure 6 illustrates how LUMEN can be customized to respond to fault management needs, in this figure three case studies are presented:

• case study 1: fault prediction using machine learning in logs in a time slot.

• case study 2: detection of SLA violation with collecting real time VNFs metrics.

• case study 3: RCA of dynamic multi-tenant environment using Bayesian networks.

Qualitative Comparison of Fault Management Frameworks

Table 2 presents a qualitative comparison of the three fault management frameworks: LUMEN, ISF [START_REF] Varga | Integration of service-level monitoring with fault management for end-to-end multi-provider ethernet services[END_REF] and Vitrage [START_REF]Openstack Vitrage[END_REF]. The three architectures are compared with different parameters related to the chosen environment and the methods and steps used in fault management. We focused on the LUMEN case study 3 in this analysis.

The conception of LUMEN and the choice of the opensource tools was made to accept distinct physical and virtual resources distributed in different infrastructures. Vitrage is designed for OpenStack environments, while the ISF framework addresses only Ethernet networks.

The three architectures filter data before processing, but only LUMEN proposes to augment the consistency of data with additional indications such as slice ID. Both LUMEN and Vitrage consider multiple tenants, while in ISF there is no concept of vitalization and network sharing in the considered environment.

The three frameworks consider the information about network topology and entities dependencies. A topology database was defined in ISF framework that stores three types of information: the node address, the node functionality and a list of connecting nodes. However, no topology builder was defined in the ISF framework since the end-to-end Ethernet services has a fixed topology, i.e. only new added entities to the network can be detected but not real time changes such as in the case of an on-demand NFV chain service. In Vitrage an entity graph builds the network topology using the information collected locally in the OpenStack modules. While in LUMEN we consider two sources to collect topology information: the VIM modules and SDN controllers.

The decision methods should be adapted to the scalability and dynamicity of NVEs. Templates or rule based techniques are used for static networks with the presence of expert knowledge . This knowledge is then encoded in rules. Acquiring this knowledge in a dynamic network such as NFV is difficult. Petri nets are dynamic systems that can only encode sequences of events. To progress towards a more ambitious diagnosis, one must adopt more appropriate methods such as Bayesian networks. The later model the network using constraints and statistical dependencies. They provide explanations about failure propagation, deal with multiple failures and capture false or lost alarms. Still, the definition of the model in a scalable and dynamic network such as NVEs is a difficult step.

The Vitrage project was developed to detect the root cause of problems in an OpenStack running environment and notify clients. The defined architecture focuses on the detection and localization of faults. The ISF framework addresses the different management steps with an integration of a service-level monitoring module. The module measures the service performance and compare these measurements with the requirements set in the SLA.

In the case of fault management, LUMEN considers all the steps and is open to the different NVEs. However, the decision module needs further investigation which is one of our future goals.

Conclusion and future work

The 5G deployment consists on a variety of NVEs: SDN, NFV, VNF chains and multi-tenant environments. The coexistence of this distinct and complementary NVEs models will open up new fault management challenges and issues where classical fault management methods are limited.

In this paper, we presented significant issues to be considered in the fault management of NVEs. Such as the dynamic network topology and the large number and variety of managed entities. These issues weaken diagnosis solutions that were designed for limited and fixed network topologies. We referenced relevant studies illustrating this fact.

We highlighted the necessity to provide automation to the canonical fault management steps and proposed LUMEN: a Global Fault Management Framework. LUMEN summarizes the fault management steps in four planes (Source, Sink, Extraction and Decision). The first three planes compose the detection step. In these three phases, the data is remodeled to prepare the decision plane where deduction methods can be deployed. In our future work we will focus in the Decision Plane. Some of the planned activities relate to the third example of Figure 6.

To implement LUMEN, we plan on using OpenStack Ocata and the Docker-engine as a Virtual Infrastructure Manager (VIM), The Clearwater virtual IMS with docker and virtual machines [START_REF]Project clearwater[END_REF] for the virtual environment and R or Python language for the RCA inference engine.

Our next steps consist on remodeling the Bayesian network technique to pinpoint faulty network components in NVEs. Still, some open issues need to be investigated for our future contribution. We will focus on the following points: the problem of scalability of Bayesian networks that was addressed by [START_REF] Hounkonnou | Empowering self-diagnosis with self-modeling[END_REF] for IMS networks, introduce VNF chains in the managed entities plus the SDN model that was addressed by [START_REF] Sanchez | Self-modeling based diagnosis of softwaredefined networks[END_REF] with a finer granularity. This case study should not affect the diagnosis reliability, address in more details the problem of tenants fault isolation and propose a real time Self-diagnosis engine for NVEs.

Figure 2 :

 2 Figure 2: Autonomic Fault Management Steps for NVEs

Figure 3 :

 3 Figure 3: Distributed VNFs in a Multi-tenant environment

Figure 4 :

 4 Figure 4: LUMEN Framework

Figure 5 :

 5 Figure 5: Elastic Stack Logs transformation

Figure 6 :

 6 Figure 6: Case studies applied to the LUMEN Framework

Table 1

 1

	: LUMEN proposed solutions for the NVEs issues and challenges
	Issues	LUMEN Solutions
		-Data filtering before storing to reduce
	Number and type	insignificant data.
	of managed objects	-Efficient storing engine and clustering methods
		to face the scalability. (Plane 1, 2)
	Dynamic network	-Network topology changes can be extracted
	topology	from collected data. ( Plane 1, 2, 3)
	Lack of network	-The centralization of logs from different
	visibility	distributed infrastructures. (Plane 2)
	Multi-tenant fault isolation.	-Additional fields (e.g. tenant ID). -Efficient Self-diagnosis to notify clients. (Plane 1 ,2,3 and 4)
	Ambiguity and consistency of data	-Unification of data format. (Plane 1, 2 and 3)

Table 2 :

 2 A qualitative comparison of fault management frameworks

		Environments	Data filtering	Data type	Fault Management steps	Decision	Dynamic topology network	Multi-tenancy
	LUMEN	NVEs	yes	-logs, -metrics, -topology	Fault management	Bayesian networks	yes	Centralization of tenants data
	ISF Framework	End-to-end Ethernet services	yes	-logs, -metrics, -topology	Fault management and Service-Level Monitoring	Petri nets	no	No resources sharing
	Vitrage OpenStack	OpenStack environments	yes	-logs, -metrics, -topology	RCA : detection and localization	Templates	Yes	OpenStack tenants

Acknowledgment

Authors would like to thank their colleague, Ayoub Bousselmi, for his early review of the paper.


Peng Lv

Sebastien Prothin

Fazila Mohd-Zawawi

Emmanuel Benard

Joseph Morlier

Jean-Marc Moschetta

Study of A Flexible Blade for Optimized Proprotor

Keywords: tilt-body MAVs, flexible proprotor, composite laminate, FSI, LDS

In the present study, a passive twist control is considered as a potential way to improve the overall flight efficiency for proprotor of Micro Air Vehicle (MAV). This paper will focus on the aerodynamic performance and deformation behaviour of a flexible laminate blade. Incorporated with a database of airfoil characteristics, Blade Element Momentum Theory (BEMT) is implemented for performance prediction of proprotor at low Reynolds numbers. The preliminary procedure is based on finding optimum twist distributions for hover and forward flight, but keeping a given chord distribution. A numerical model is developed using a combination of aerodynamic model based on BEMT, and structural model based on anisotropic beam finite element, in order to evaluate the coupled structural and the aerodynamic characteristics of the deformable proprotor blade. The numerical model -Fluid Structure Interaction (FSI) was then validated by means of shape reconstruction from LDS (Laser Displacement Sensor) outputs. It can be concluded that the proposed experiment technique is capable of providing a predictive and reliable data in blade geometry and performance for rotor mode. The FSI approach is also valid as a reliable tool for designing and analyzing the MAV proprotor made of composite material.

Introduction

Tilt-rotor aircraft has been developed be multifunctional in order to offer a wide range of services. It can fly in both of hover and forward flight. In early 1950s, tilt-rotor aircraft was started to be developing. XV-3 of Bell company operated first transition from hover to forward flight. In 1970s, XV-15 was developed by Bell company to demonstrate the feasibility of tilt-rotor concept. Then it was moved to NASA and the U.S. Army for further study in aeroelasticity. The successful XV-15 led to the project of V-22. V-22 is the world's first production tilt-rotor aircraft. With respect to a proprotor of tilt-rotor aircraft, in hover, the inflow velocity is small and the proprotor must provide high thrust to support aircraft weight. By contrast, in forward flight, the inflow velocity is relatively large and the low thrust is just to overcome the drag. The difference in the inflow and thrust requirement between the two flight modes suggests different blade twist and chord distributions. In 1983, McVeigh obtained the twist of XV-15 proprotor through linear interpolation of twist between rotor and propeller by a compromise [START_REF] Mcveigh | Aerodynamic design of the XV-15 advanced composite tilt rotor blade 39[END_REF]. Although this trade-off provided acceptable performance on XV-15, the stiff proprotor with certain twist cannot maximize the efficiency for both flights. In 1988, Nixon proposed a passive blade twist control method for the proprotor of XV-15 [START_REF] Nixon | Improvements to tilt-rotor performance through passive balde twist control[END_REF]. The study demonstrated successfully the feasibility of the passive blade control on conventional tilt-rotor aircraft. The tilt concept of typical aircraft attracts interest of researchers who are working on MAVs. In 2008, Shkarayev and Moschetta introduced the efforts on the aerodynamic design of a tilt-body MAV named miniVertigo, which had a tilt-body configuration [START_REF] Shkarayev | Aerodynamic design of micro air vehicles for vertical flight[END_REF]. The wind tunnel measurements were conducted for a motor, a wing, and an arrangement of a wing with a motor. The results were realized in the design of a prototype of tilt-body MAV which was successfully tested in flight. The small proprotors also suffer the problem caused by different twist between hover and forward flight. However, due to the small size of MAV, the complex tailored blade cross section for passive twist control based on conventinal tiltrotor aircraft is not available any more on it. Therefore, composite laminate is explored to be a more practical method for proprotor blade of MAVion, which is a tilt-body MAV developed by ISAE. In this study, a passive twist control is considered as a potential way to improve the overall flight efficiency of MAV proprotor. The proprotor blade made of composite material is preferred to be used. It is due to their potential benefits such as aeroelastic tailoring, ability to manufacture, more refined aerodynamic designs, significant enhancements in fatigue performance and damage tolerance of the blade. The blade is expected to be deformed in torsion under different airloads and structural loads. This paper is aimed at developing an evaluation of design techniques of a composite flexible proprotor.

The key issue to study flexible blade is to observe the deformation accurately. Optical measurement techniques have been developing for a couple of years in applications of aerodynamics, materials and structure, such as Holographic Interferometry (HI), Electronic Speckle Pattern Interferometry (ESPI), Projection Moiré Interferometry (PMI) and Digital Image Correlation (DIC) [START_REF] Williams | Optical methods in engineering metrology[END_REF]. In 1998, Fleming obtained the 3-D deformation of rotor blade using PMI technique [START_REF] Fleming | Measurement of rotorcraft blade deformation using projection moire interferometry[END_REF]. However, it has low sensitivity for in-plane deformation and moderate for out-of-plane deformation. By contrast, DIC has a relatively high sensitivity that can reach 1/30,000 of the test field [START_REF] Schmidt | Full-field dynamic displacement and strain measurement using advanced 3D image correlation photogrammetry: Part 1[END_REF]. In 2011, Lawson demonstrated the deformation of a rotating blade using DIC [START_REF] Lawson | Measurement of deformation of rotating blades using digital image correlation[END_REF]. The technique was found to have many advantages including high resolution results, non-intrusive measurement, and good accuracy over a range of scales. However, DIC needs a preprocessing which is to apply a stochastic speckle pattern to the surface by spraying it with a high-contrast and non-reflective paint. This complex painting will probably affect the stiffness of blade. Hence, in this study, LDS was developed to measure blade deformation and validate the FSI model.

Aerodynamic modelling

The aerodynamic model based on BEMT is used as a tool to compute the aerodynamic loadings. In the classical approach of rotor analysis, lift polar is a linear function [START_REF] Leishman | Principles of helicopter aerodynamics[END_REF]. Hence, in order to consider non-linear airfoil characteristics prevalent in low Reynolds number regime (generally, Re<70,000), a database of airfoil char-acteristics was incorporated into modified BEMT. To solve the BEMT equations, the blade should be numerically discretized into a series of small elements. In the classical equation,

4F λ 2 rdr = σC lα 2 (θr 2 -λr)dr, (1)
where λ is inflow ratio, F is Prandtl's tip-loss factor, r is nondimensional radius, dr is the nondimensional length of each element, C lα is 2-D lift-curve-slope of airfoil, and θ is the pitch angle at the midspan of each element. By contrast, in the modified version, the incremental thrust coefficients are described using lift coefficient C l directly instead of lift-slope C lα ,

Eq. 2 allows for the solution of induced flow ratio along the blade radius for a given rotor with certain chord distribution, collective pitch, twist distribution and airfoil section. The result will be converged through a few iterations. In propeller analysis, the form of BEMT developed by Adkins was employed [START_REF] Adkins | Design of optimum propellers[END_REF]. In this model, the inflow ratio or inflow angle is iteratively computed until the convergence criteria are reached. The database of airfoil characteristics are obtained from Xfoil, which is an airfoil design and analysis code developed by Drela [START_REF] Drela | An analysis and design system for low Reynolds number airfoil, low Reynolds number aerodynamic[END_REF]. BEMT model was validated using existing performance data of rotor and propeller. The detailed geometries of them are described in [START_REF] Bohorquez | Rotor hover performance and system design of an efficient coaxial rotary wing micro air vehicle[END_REF] and [START_REF] Smedresman | Design, fabrication, analysis, and testing of a micro air vehicle propeller[END_REF]. The airfoil of rotor defined by NACA0012 while the propeller airfoil is based on a modified MA409. Xfoil can not always give a converged solution of airfoil performance for given airfoils. Hence, in order to obtain a smooth database for BEMT convergence, it is necessary to conduct interpolation for lift coefficient C l , drag coefficient and C d and moment coefficient C m with variation to angle of attack and Reynolds numbers. The Ncrit to predict the flow transition is 0.1 for rotor and 7.0 for propeller corresponding to different turbulence intensities [START_REF] Smedresman | Design, fabrication, analysis, and testing of a micro air vehicle propeller[END_REF][START_REF] Schafroth | Aerodynamics, modeling and control of an autonomous micro helicopter[END_REF]. Fig. 1 and Fig. 2 show the interpolated results on C l and C d . Fig. 3 and Fig. 4 exhibit the validation of performance predictions of the rotor and propeller.

The BEMT model over predicts the thrust of rotor slightly while the torque matches well with experiment data. This is probably from that, as the wake becomes turbulent for heavy loading, air with momentum is transported from the outer flow region into the wake. The accuracy of BEMT model is decreased by the strong self-induction wake. Also, under heavy disk loading, the assumption of rigid wake sheet for Prandtl tip-loss function becomes not accurate enough. With respect to propeller, it is found that the developed BEMT method incorporated with database agrees better with the experimental data compared to the linear model Qprop [START_REF] Drela | Qprop user guide[END_REF]. However, for high disk loadings, the results predicted by Qprop which is capable of considering self-induction wake, show potentially better trends. Overall, the BEMT is an efficient tool for at least preliminary analysis of the spanwise distribution of airloads, especially due to the negligible computational cost of this method. 3 Twist effect on the efficiency of small proprotor

In order to study the twist effect on the efficiency of small proprotor, five blades with varied built-in twist are defined first, as shown in Tab. 1. For example, twist -10 represents that twist of blade tip is lower than the counterpart of blade root by 10 • with linear distribution. The airfoil of the blade is flat plate with a thickness of 2.5%. Tab. 2 exhibits the operation conditions for both hover and forward flight.

FSI Model

In order to compute the blade deformation under airloads and centrifugal loads, the coupled FSI model of the proprotor combines the aerodynamic model -BEMT and structural analysis model model -FBEAM, which was developed by Mohd-Zawawi [START_REF] Moha-Zawawi | Study of a flexible UAV proprotor 48th International Symposium of Applied Aerodynamics Saint-Louis[END_REF]. To begin with, the basic geometry (zero deflections) is assumed. Then the structural model calculates the blade deformations under airloads and centrifugal loads. The deformed blade shape is used to update the aerodynamic model. The quasi-steady aerodynamic loads are computed for the update geometry, using the 2D aerodynamic theory based on BEMT. The aerodynamic loads are transferred to the beam nodes as concentrated forces. A new structural analysis is performed to calculate the deformed shape of the blade under the influence of aerodynamic and centrifugal loads. The variation of the blade twist angle along the blade is monitored for convergence. The interaction is repeated until equilibrium between deformation and loadings is achieved. After convergence the proprotor performance characteristics are computed. The approach described above was applied to the constant chord untwisted 2-bladed system for small proprotor made of laminate composite. The deflection results in two basic modes of deformation; spanwise bending and torsion. The effects of camber changes are not included in this study due to limitation in the modeling. Firstly, a static analysis on a rectangular planform with dimension (length 190.5mm, width 12.7mm, thickness 3.175mm) was performed. Carbon/Epoxy with material properties E11=129GPa, E22=9.4GPa, E33=9.4GPa, G12=5.16GPa, G13=4.3GPa, G23=2.54GPa, µ=0.3 and ρ=1550kg/m 3 was used.

The ply angle with respect to pitch axis was 30 • . The result obtained in FBEAM was validated by commercial structural analysis program namely MSc. Nastran and the validation result is illustrated in Fig. 6. to reference plane. Meanwhile, square wave in time domain is measured by optical RPM sensor. Then, average RPM and angular speed at each blade section can be extracted from the square wave. Furthermore, with the azimuth of feathering, polar coordinates is possible to be transferred to Cartesian coordinates X (along spanwise of blade) and Y (along chordwise of blade). Combining coordinates X, Y and Z, a polynomial surface fitting is performed to obtain the bending and torsion of rotating blade. As we can see, coordinates Z can be measured in a direct way by laser. However, coordinates X and Y are instead calculated by local radius of each scanning, angular speed of blade and sampling time,

X = rcos(wt) = rcos( RPM 60 2πt), (3)
Y = rsin(wt) = rsin( RPM 60 2πt), (4)
where r, w and t are local radius, angular speed and sampling time, respectively. The uncertainties U X and U Y are defined by bias limits and precision limits using the Root-Sum-Square (RSS) method [START_REF] Taylor | An introduction to error analysis: the study of uncertainties in physical measurements[END_REF],

U X = (B 2 X + P 2 X ) 1 2 , (5)
U Y = (B 2 Y + P 2 Y ) 1 2 , (6)
where B X , P X , B Y and P Y are bias limits of X, precision limits of X, bias limits of Y and precision limits of Y . The are defined by,

( B X X ) 2 = ( 1 X ∂X ∂r B r ) 2 + ( 1 X ∂X ∂RPM B RPM ) 2 + ( 1 X ∂X ∂t B t ) 2 , (7)
( P X X ) 2 = ( 1 X ∂X ∂r P r ) 2 + ( 1 X ∂X ∂RPM P RPM ) 2 + ( 1 Y ∂X ∂t P t ) 2 , (8)
( B Y Y ) 2 = ( 1 Y ∂Y ∂r B r ) 2 + ( 1 Y ∂Y ∂RPM B RPM ) 2 + ( 1 Y ∂Y ∂t B t ) 2 , (9)
( P Y Y ) 2 = ( 1 Y ∂Y ∂r P r ) 2 + ( 1 Y ∂Y ∂RPM P RPM ) 2 + ( 1 Y ∂Y ∂t P t ) 2 , (10)
With respect to the real Z of blade to reference plane, as shown in Fig. 8, it is defined by,

Z = Z ′ cos(θ), (11)
where coordinates Z ′ are the distance data to measuring plane obtained from LDS directly, θ is the measuring angle of lasers. This is to say, coordinates Z are needed to be transformed from measuring system to the reference system for the reconstruction of blade shape. U Z , it is defined by uncertainty of laser measurement directly,

U Z = (B 2 laser + P 2 laser ) 1 2 , ( 12
)
where B laser and P laser are bias limits and precision limits of laser data. All of bias limits are constants and they will not be affected by sampling. Local radius r is totally determined by track system while time t is from quartz crystal reference.

Here, the precision limits of RPM are considered for uncertainties of coordinates X and Y . The precision limit of RPM is defined as,

P RPM = 1.96σ RPM √ N RPM , (13)
where σ RPM is the standard deviation of RPM and N is the RPM sampling number. The distance to a reference was measured by laser sensor first to study precision limits of coordinates Z, the precision of displacement data is defined as,

P laser = 1.96σ laser √ N laser , (14)
where σ laser is the standard deviation of measured distance by laser and N is the laser sampling number.

Results

The laminate blade was tested at RPM 1,500 in hover mode. The comparison results between BEMT model and experimentsin performance are shown in Fig. 9. The errorbars represent the stand deviation of the mean of sampling data. In simulation, the input blade geometry for BEMT was assumed without deformation. As can be seen, the BEMT agrees well with experimental data of torque. However, it overpredicts the thrust. This is probably caused by the obvious bending deformation of rotating blade since it can generate the thrust vector in spanwise of blade. The precision limits of reference measurement was given by -35.8452±0.0008, as shown in Fig. 10(a). Besides, according to user's manual of KEYENCE LK-G502, when the detection range is between -250mm and 500mm, the bias limits is given by ±0.1%. Based on Eq. 11, the uncertainty of coordinates Z is ±0.1%. This is to say, compared with bias limits, precision limits of laser can be negligible. The normal density of RPM samples from optic sensor are shown in Fig. 10(b). Finally, the uncertainty of RPM was given by 1,517.5±11.4. The uncertainty U X and U Y were calculated using Eq. 5 and 6, as shown in Fig. 11 and 11. The blade was reconstructed through 3-D polynomial surface fitting of second order (Fig. 13(a)). Fig. 13(b) shows the bending distribution of rotating blade, and the errorbars are defined by the standard deviation of surface fittings along spanwise of blade. FSI simulation and experiment are capable of obtaining the bending distribution with approximately same order of magnitude. However, the curvatures of bending in FSI model and measurement are different especially around the blade root. The difference could be from composite properties in FSI model. Hence, it is necessary to measure the properties instead of using ideal characteristics. Also, initial scanning position of LDS may cause the curvature difference in bending.

Conclusion

In conclusion, the evaluation of design techniques, both for the aerodynamic performance and for the structural behavior of a composite flexible proprotor has been presented in rotor mode. The numerical model BEMT has capability of prediction in proprotor aerodynamics. The results show high twist of blade is suitable for propulsive efficiency while low twist is beneficial for hovering efficiency. FSI model has also been validated by LDS technique. The difference of bending distribution is expected to be studied further in terms of materials properties and experiment procedure. It can be concluded that the developed numerical model is a reliable tool for designing and analyzing the proprotor made of composite material. Additionally, current laminate was not observed deforming in beneficial torsion. In future, in order to obtain obvious torsion, the study will focus on the configuration design of laminate blade. Likewise, the experiments will also be expanded to propeller mode.

Figure 1 :Angle

 1 Figure 1: Airfoil characteristics of NACA0012

Figure 2 :

 2 Figure 2: Airfoil characteristics of modified MA409

Figure 3 :

 3 Figure 3: BEMT validation for rotor

Figure 4 :

 4 Figure 4: BEMT validation for propeller

Fig. 5 (

 5 a) shows the twist effect on hovering efficiency -Figure of Merit (FM) with the variation to thrust coefficient. Thrust coefficient is adjusted by collective pitch for the five blades. The blade with built-in twist 10 • has the maximum efficiency at C T =0.01. In Fig. 5(b), the propulsive efficiencies of the five blades are all analyzed under C T =0.003. They vary with the advance ratios in forward flight. The blade with twist 30 • exhibits highest efficiency of all at advance ratio J=1.25. Overall, as we can see from Fig. 5, high twist of blade is beneficial for forward flight while low twist can improve hovering efficiency. According to the mission requirement of MAVion, the solution of proprotor for MAVion is given in Tab. 3. Twist effect on propeller efficiency

Figure 5 :

 5 Figure 5: Twist effect to small proprotor Hover Forward flight Built-in twist (degree) 10 30 Collective pitch (degree) 31 54

Figure 6 :

 6 Figure 6: Static validation of FBEAM

Figure 7 :Figure 8 :

 78 Figure 7: Experimental rigs

Figure 9 :

 9 Figure 9: Aerodynamic performance of laminate blade

Figure 10 :

 10 Figure 10: Normal density distribution

Figure 11 :Figure 12 :

 1112 Figure 11: Uncertainties on upper surface

Figure 13 :

 13 Figure 13: LDS results

Table 1 :

 1 The blades with varied built-in twist

		Hover Forward flight
	RPM	1500	1200
	Velocity (m/s)	0	10
	Advance ratio, J	/	1.25
	Thrust (N)	2	0.3
	Thrust coefficient 0.01	0.003

Table 2 :

 2 Operation conditions of MAVion

Table 3 :

 3 Optimized proprotor of MAVion

Acknowledgements

The authors thank Rémy CHANTON, Patrick MOREL, Philippe BARRICAU, Serge GERARD, Patrick CAZENAVE, Philippe MOULIGNE and Xavier FOULQUIER for advice and assistance. The research is supported by China Scholarship Council (CSC) and Malaysian government.