A Rabhi

H Imine

N M' Sirdi

Y Delanne

Observers With Unknown Inputs to Estimate Contact Forces and Road Profile

Keywords: Nonlinear observers, Sliding Modes, Vehicle-Road Interaction Models, State Estimation, Tire Forces, Road Profile

This paper presents sliding mode observers designed to estimate tire forces and road profile. Tire forces affect the vehicle dynamic performance and behavior properties. The tire forces and road friction are difficult to measure and their modelling is rather complex. In this work we deal with a simple model of vehicle combined with sliding mode approach to develop robust observers.

I. Introduction

Knowledge of tire forces is essential for systems such as antilock braking systems (ABS), traction control systems (TCS) and electronic stability program (ESP). Vehicle dynamics depend largely on the tire forces which are nonlinear functions of wheel slip and slip angles and depend on some factors such as tire wear, pressure, normal load tire road interface properties [START_REF] Ackermann | Robust control prevents car skidding[END_REF][2] [START_REF] Gipser | Dynamical tyre forces response to road unevennesses[END_REF]. Recently, many analytical and experimental studies have been performed on estimation of the frictions and contact forces between tires and road [START_REF] Canudas De Wit | Dynamic Friction Models for Road/Tire Longitudinal Interaction[END_REF] [START_REF] Clover | Longitudinal Tire Dynamics[END_REF] [START_REF] Delanne | Tire/Road Friction Perfor-mance Models from on-site Measurements[END_REF]. The tire forces affect the vehicle dynamic performance and behavior properties. Thus for vehicles and road safety analysis, it is necessary to take into account the contact force characteristics. However, tire forces and road friction are difficult to measure directly and complex to be precisely represented by some deterministic equations. Vehicle dynamics depends largely on the tire forces represented by the nonlinear functions of wheel slip. The tire models encountered are complex and depend on several factors (as load, tire pressure, environmental characteristics, etc.). This makes on line estimation of forces and parameters difficult for vehicle control applications like detection and diagnosis for driving monitoring and surveillance. In this paper, modelling of the contact forces and interactions between a vehicle and road is revisited in the objective of on line force estimation using robust observers coupled with a robust and adaptive estimation of contact forces. We propose a robust observer to estimate the vehicle state and an adaptive estimator for tire forces identification [START_REF] Nacer | Observateurs robustes et estimateurs pour l'estimation de la dynamique des véhicules et du contact pneu -route[END_REF]. The designed observer is based on the sliding mode approach. The main contribution is on-line estimation of inputs (the tire forces and road profile) needed for control. In this work, we deal with a simple vehicle model coupled with an appropriate wheelroad contact model in order to estimate contact forces. Then, we develop a method to observe tire forces and road profile.

This paper is organized as follows. Section 2 deals with the vehicle description and modelling for estimation of contact forces. The design of an observer and an adaptive tire force estimation is presented in section 3. Section 4 is devoted to develop an observer with unknown inputs to estimate the road profile. Some results about the states observations and estimation of the two kinds of unknown inputs are presented in section 5. Finally, some remarks and perspectives are given in a concluding section.

II. VEHICLE MODELLING

In the literature, many studies deal with vehicle modelling [START_REF] Dugoff | An analysis of tire traction properties and their influence on vehicle dynamic performance[END_REF][9] [START_REF] Lee | Adaptative vehicle traction force control for intelligent vehicle highway systems[END_REF]. The objective may be either analysis for better design and features enhancement or increase of safety and maniability. The vehicle is a complex mechanical system that exhibits nonlinear behaviors. Commonly, the proposed and used models are not very complicated and give partial representation of the system dynamics. It would be relatively difficult and intricacy to involve more complete models and to define the size of different parameters. Several models have been considered in literature for analysis of the road -vehicle interaction and its consequence on the behavior (see eg Figures 1 and2). The motions (longitudinal, lateral, and vertical) depend on interaction between the wheels and the road, the disruptions and the gravity. We generally we can distinguish three types of vehicle models:

-Longitudinal (Figures 1) -Lateral -longitudinal + lateral (Figure 2).

Let us consider, for our case, the last type of representation (Figure 2). We can define as dynamic equations of the vehicle:

m   . v x . v y . v z   =   F x F y F z   (1) J    .. θ .. φ .. ψ    =   M x M y M z   (2)
The wheel angular motion is described by:

J f i
.

ω f i = T f i -r f i F xf i J ri . ω ri = T ri -r ri F xri (3)
where

v = [v x , v y , v z ] T : vehicle velocities along x, y, z,
The subscripts f , r stand for front and rear wheels respectively.

ω f i and ω ri : rotation speeds of wheels (front and rear), r f i and r f i :wheel radius at front and rear respectively, J f i and J ri : front and rear wheels inertia respectively, F xf , F yf : longitudinal, lateral forces on front wheels, F xr , F yr :longitudinal, lateral force on the rear wheels, T f i and T ri : torques applied to front and rear wheels, F x , F y , F z , M x , M y and M z : forces and moments with respect to x, y and z. and θ φ ψ.

We propose therefore, with these dynamic equations, to use a longitudinal-lateral (tire-road) contact model, in order to take into account the slip effects and ground forces. This will provide the latter variables which inputs of the depicted equations.

A. Simplified model

Known as the bicycle model ( [10][11]), the structure described by Figure 2, gives a fairly good representation of the vehicle behavior in the (x,y).plane. This representation is obtained by replacing the two front wheels by an equivalent one. A similar approach is applied for the two rear wheels. This representation considers the following assumptions :

1) The epicenter is assumed to be on road level 2) Neglect the roll, pitch, and vertical motion;

3) The road is assumed to be perfectly flat 4) Neglect influence of aerodynamic side forces.

The dynamic model represents the longitudinal, lateral, yaw motions and the rotation of the wheels as shown in figure 2. The resulting equations of the simplified vehicle model are then

m . V x = F xf (cos(δ F ) -F yf sin(δ F ) + F xr m . V y = F xf sin(δ F ) + F yf cos(δ F ) + F yr (4) J z .. ψ = F xf l f sin(δ F ) + F yf l f cos(δ F ) -l r F yr
δ F is the front wheel steering angle (δ r = 0); l r and l f distances between the center of gravity of the vehicle and the rear and front axis respectively.

B. Tire modelling

The forces are produced by contact between the road and tires. They are transmitted trough the dynamics of wheels and vehicle. They are of major importance for the dynamic behavior of a road vehicle. Hence, accurate tire models are necessary components of models aimed at analyzing or simulating vehicle motion in real driving conditions. The driver can then control the vehicle trough this dynamic. A lot of work has been done in the area of tires model fitting and estimation [START_REF] Gim | Analytical model of pneumatic tyres for vehicle dynamic simulations part1: Pure slips[END_REF][13]. Many models have been previously used to describe the tire-forces. Some of them are theoretical in the sense that they aim at modelling the physical processes that generates the forces. Other ones are empirically oriented and their aim is to describe observed phenomena in a simple way.

1) Lateral dynamics and transient phenomenon: The dynamic behavior of the transverse motion was the subject of several works [START_REF] Drakunov | ABS control using optimum search via sliding modes[END_REF]. In [ [START_REF] Pacejka | Magic Formula Tyre Model with Transient Properties[END_REF], Pacejka describes by a first order model, the variations of the lateral strength and the moment of auto-alignment in presence of weak values of the slip angle, while using the notion of the relaxation length.

To illustrate the concept of the relaxation length, let us consider the dynamic variations of strength lateral F y in the case of weak rates of slip. Suppose that the variation of the vertical strength is weak as . F z = 0, then the variation F y is associated to a variation of the lateral speed of point of contact represented mainly by V cy . So to describe the transient, the variation of F y is represented by a differential first order equation:

σ yi . F yi + V xi F yi = C y V y i = f, r (5)
where: C y is the rigidity of the lateral slip, and σ y represents the length of relaxation. We can extend this equation ( 5) to large slips to get:

σ yi . F yi = -V x (F yi -F yi0 ) + C y V y i = f, r (6)
the unknown parameters F yi0 is intersection of the tangent

F xi = -V x (F xi -F xi0 ) + C x (V x -r i ω i ) b) during acceleration , i = f, r σ xi . F xi = -r i ω i (F xi -F xi0 ) + C x (V x -r i ω i ) (7)
Finally, during acceleration, the model can be written:

.

F xf = - V x σ lf F xf + V x σ lf F xf 0 + C x σ lf V x - C x σ lf r f ω f . F xr = - V x σ lr F xr + V x σ lr F xr0 + C x σ lr V x - C x σ lr r r ω r . F yf = C y σ tf V y - V x σ tf F yf + V x σ tf F yf o (8)
.

F yr = C y σ tr V y - V x σ tr F yr + V x σ tr F yro

III. Adaptive Estimation of Tire forces A. Expression of the robust observer

The system equations can be rewritten in the following state space form:

         • x 1 = x 2 • x 2 = Ω(u)x 4 • x 3 = Gu -[H ; 0] x 4 • x 4 = Ψ(X)θ (9)
where x 1 = (x 11 , x 12 , x 13 ) = (x, y, ψ) To estimate both state an unknown parameters we propose the following sliding mode based observer [START_REF] Hadri | Non-linear longitudinal tire force estimation based sliding mode observer[END_REF]:

x 2 = (x 21 , x 22 , x 23 ) = (V x , V y , . ψ) x 3 = (x 31 , x 32 ) = (ω f , ω r ) x 4 = (x 41 , x 42 , x 43 , x 44 ) = (F xf , F xr , F yf , F yr ) x 4x = (x 41 , x 42 ), x 4y = (x 43 , x 44 ) Ω(u) =   1 m cos(δ F ) 1 m 1 m sin(δ F ) 0 1 m sin(δ F ) 0 1 m cos(δ F ) 1 m 1 Jz l f sin(δ F ) 0 1 Jz l f cos(δ F ) -1 Jz l r   G =   0 0 0 0 1 J f 0 0 0 1 Jr   ; H = r f J f 0 0 0 0 r f J f 0 0 Ψ(X) = Ψ 1 (X) 0 0 Ψ 2 (X) Ψ 1 (X) = C x x 11 -
           . x2 = Ω(u) x 4 + Λ 2 sign( x 2 ) . x3 = Gu -H x 4x + Λ 3 sign( x 3 ) . x4 = Ψ( X) θ + Λ 4 sign( x 2 ) + Λ 5 sign( x 3 ) . θ = η (10)
where x i represents the observed state vector, x 2 = x 2 -x 2 , x 3 = x 3 -x 3 are the state estimation errors.

The observer gains Λ i and the unknown η will be defined thereafter..We can write

∆Ψ 1 = Ψ 1 (X) -Ψ 1 ( X) = -x 21 x41 0 0 -x 21 x42 0 0 0 0 ∆Ψ 2 = Ψ 2 (X) -Ψ 2 ( X) = -x 43 x 21 0 0 -x 44 x 21 0 0 0 0 we have Ψ(X) -Ψ( X) ≤ a x 4x + b x 4y
We can write

Ψ( X) θ -Ψ(X)θ = -Ψ( X) θ + ∆Ψθ ε = ∆Ψθ = (Ψ(X) -Ψ( X)θ = ∆Ψ 1 .Θ 1 ∆Ψ 2 .Θ 2 ε = -x 21     θ 1 0 0 0 0 θ 2 0 0 0 0 θ 5 0 0 0 0 θ 6     .x 4 = -x 21 . .x 4 ε ≤ (c. max(θ)) x 4          . x2 = Ω(u 1 ) x 4 -Λ 2 sign( x 2 ) . x3 = -H x 4x -Λ 3 sign( x 3 ) . x4 = -x 21 . .x 4 -Ψ( X) θ -Λ 4 sign( x 2 ) -Λ 5 sign( x 3 ) . θ = -η B.
Convergence analysis Define the state estimation errors x 2 , x 3 , x 1 = x 1 -x 1 , and θ = θ -θ, the parameters estimation errors can be written as:

         . x2 = Ω(u) x 4 -Λ 2 sign( x 2 ) . x3 = H x 4x -Λ 3 sign( x 3 ) . x4 = Ψ(X)θ -Ψ( X) θ -Λ 4 sign( x 2 ) -Λ 5 sign( x 3 ) . θ = -η (11)
In order to study the observer stability, we proceed, step by step. To proceed, let us consider the following Lyapunov function:

V 2 = 1 2 x T 2 x 2 (12)
The time derivative of this function is given by

V2 = x T 2 (Ω(u 1 ) x 4 -Λ 2 sign( x 2 )) (13)
By chosing Λ 2 > Ω(u 1 ) x 4 , then V2 < 0 therefore, from sliding mode theory, the surface defined by x 2 = 0 is attractive, leading x 2 to converge to x 2 in finite time t 0 . Moreover,we have .

x 2 = 0 ∀t ≥ t 0 and then for ∀t ≥ t 0 :

sign equ ( x 2 ) = Λ -1 2 (Ω(u 1 ) x 4 ) ( 14
)
where sign equ represents an equivalent form of the sign function on the sliding surface. Now let us consider a (second) Lyapunov V 3 and its time derivative .

V 3 :

V 3 = 1 2 x T 3 x 3 (15) V3 = x T 3 (-H x 4x -Λ 3 sign( x 3 )) ( 16
)
Also by considering Λ 3 > H x 4x , then .

V 3 < 0 therefore, the surface defined by x 3 = 0 is attractive and we obtain sign equ ( x 3 ) = Λ -1

3 (H x 4x ) . Then the system becomes (for ∀t ≥ t 0 ): . x2 = 0 and

. x3 = 0 ( 17
) . x4 = -A. x 4 -Ψ( X) θ and . θ = -η (18) A = Λ 4 Λ -1 2 Ω(u) + Λ 5 Λ -1 3 [H; 0 2 ] + x 21 . ( 19
)
Now consider the function of Lyapunov :

V 4 = 1 2 x T 4 G x 4 + 1 2 θ T P θ (20)
The time derivative of this function is given by

V4 = x T 4 G . x4 + θ T P -1 . θ (21) V4 = -x T 4 GA. x 4 -θ T P -1 η + Ψ( X) T G T x 4 (22)
An appropriate choice of the adaptation law would be η = P Ψ( X) T G T x 4 . Knowing that one will have (on the sliding surface) after convergence of x 2 and x 3 (in average): sign equ ( x 2 ) = Λ -1 2 (Ω(u 1 ) x 4 ) and sign equ ( x 3 ) = Λ -1 3 H x 4x , finally one can approach of this case (while using sign( x 2 ) and sign( x 3 )). One must choose G T to insure that A > 0 by the choice of Λ 2 , Λ 3 , Λ 4 and Λ 5 to assure a good convergence (see [START_REF] Canudas De Wit | Dynamic Friction Models for Road/Tire Longitudinal Interaction[END_REF])

IV. Estimation of the road profile

In the previous section we estimate the input forces, here we develop a method based on sliding mode to observe the road profile.

A. Observer design

The model considered for this objective is described in figure [START_REF] Clover | Longitudinal Tire Dynamics[END_REF]. When considering the vertical displacement along z axis, the model can be written as:

M q + C q + Kq = ζ (23)
where ( q, q) represent the velocities and accelerations vector respectively. G is related to the gravity effects. M is the inertia matrix, C is related to the damping effects, K is the springs stiffness vector, and q ∈ 8 is the coordinates vector defined by:

q = [z 1 , z 2 , z 3 , z 4 , z, θ, φ] T (24)
The variable ζ ∈ R 4×4 is defined by:

ζ = A 11 U + B 11 U (25)
where U = (u 1 , u 2 , u 3 , u 4 ) T is the vector of unknown inputs which characterize the road profile. The matrices A 11 and B 11 ∈ R 4×4 are:

A 11 =     k r1 0 0 0 0 k r2 0 0 0 0 k f 1 0 0 0 0 k f 2     (26)
B 11 =     B r1 0 0 0 0 B r2 0 0 0 0 B f 1 0 0 0 0 B f 2     (27)
(M i , C ij , i, j = 1, 2) are defined in R 4×4
, where C 11 is a positive diagonal matrix. The matrices M, C and K are defined by:

M = M 1 0 0 M 2 , C = C 11 C 12 C 21 C 22 (28)
Let us define a state vector:

x = (x 1 , x 2 ) T where x 1 = q = (x T 11 , x T 12 ) T , x 2 = q = (x T 21 , x T 22 ) T , x 11 = [z 1 z 2 z 3 z 4 ] T , x 12 = [z θ φ] T , x 21 = [ ż1 ż2 ż3 ż4 ] T and x 22 = [ ż θ φ] T .
We put then the model ( 23) in the state form as follows:

       ẋ1 = x 2 ; x 1 = q ẋ21 = -M -1 1 (C 11 x 21 + C 12 x 22 + K 11 x 11 + K 12 x 12 + ζ) ẋ22 = -M -1 2 (C 21 x 21 + C 22 x 22 + K 21 x 11 + K 22 x 12 ) y = [x T 12 , x T 22 , x T 31 ] T ; x 31 = ẋ21
.

(29) The positives matrices K 11 and C 11 are defined in R 4×4 , whereas the matrices K 12 and C 12 are in R 4×3 . The positives matrices K 22 and C 22 are in R 3×3 , whereas the matrices K 21 and C 21 are defined in R 3×4 . The system being posed, we can define an observer and study its convergence. In order to estimate the state vector x and to deduce the unknown input vector ζ, we propose the following sliding mode observer: In order to study the observer stability and to find the gain matrices, let us consider the following Lyapunov function:

                 . x11 = x12
V 1 = 1 2 xT 22 x22 (32)
The time derivative of this function is given by:

V1 = -x T 22 M -1 2 (C 21 x21 + C 22 x22 + K 21 x11 )(33) -x T
In order to estimate the elements of the unknown vector and according to (25), we can solve the following equation:

ζ = A 11 U + B 11 dU dt (37)
When we consider the initial conditions U (t = 0) = 0, we obtain from (37), the unknown input vector so that:

u i = ζ i a ii (1 -exp(- a ii b ii t)), i = 1..4 ( 38
)
where a ii and b ii are the elements of the matrices A 11 and B 11 respectively.

V. Simulation results

In this section, we give some results to validate our approach. In simulated models, forces are generated by the Magic formula tire model [ [START_REF] Pacejka | Magic Formula Tyre Model with Transient Properties[END_REF]]. The steering angle applied is shown in figure [START_REF] Canudas De Wit | Dynamic Friction Models for Road/Tire Longitudinal Interaction[END_REF]. The figures (4) and [START_REF] Clover | Longitudinal Tire Dynamics[END_REF] show the convergence of the estimated state vectors to the actual ones in finite time. In figure [START_REF] Delanne | Tire/Road Friction Perfor-mance Models from on-site Measurements[END_REF] we show the asymptotic convergence of the tire force to actual values. Then performance of sliding mode observer and adaptive estimation are satisfactory. The simulation results show that the adaptive observer is robust with respect to parameters uncertainties and the changes on the road conditions. After these simulations emphasizing features of the first method, let us consider some experimental results to validate road-profile estimation. The estimated road profile is compared to the profile measured by a longitudinal profile analyzer (LPA) developed at the LCPC Laboratory. It is equipped with laser sensor and accelerometer to measure the elevation of the road profile (see Figure [START_REF] Nacer | Observateurs robustes et estimateurs pour l'estimation de la dynamique des véhicules et du contact pneu -route[END_REF]). Through Figure [START_REF] Dugoff | An analysis of tire traction properties and their influence on vehicle dynamic performance[END_REF], we present the behavior of the front right wheel estimator and the equivalent velocities. In the first two plots on top, the vertical displacement of the wheels are presented. The bottom of this figure, represent the velocities. We can see that the estimated vertical velocities of the wheels are accurate compared to the measured signals. However, some chattering occurs in the estimation of positions. Figure [START_REF] Imine | Observateurs à entrées inconnues par mode glissant appliqués l'estimation du profil de route[END_REF] presents a comparison between the measured states, namely, vertical velocity of the body, roll velocity and the pitch velocity, with the estimated one. We notice that the estimated velocities converges well towards observed one. The Figure [START_REF] Lee | Adaptative vehicle traction force control for intelligent vehicle highway systems[END_REF], presents both the measured road profile and the estimated one. We can then observe, that the estimated values are quite close to the true ones (LPA measures).

VI. Conclusion

In this paper, we have developed a new estimation method for vehicle dynamics based sliding mode observer. Simulation results illustrate the ability of this approach to give well estimation of both vehicle states, road profile and tire forces. The robustness of the sliding mode observer versus uncertainties on the model parameters has also been shown in simulation. The future work will try to extend these observers and combine estimations.

Fig. 1 .

 1 Fig. 1. Half Vehicle Longitudinal and Vertical Model (x, z, ϕ)

Fig. 2 .

 2 Fig. 2. Half Vehicle Longitudinal and Lateral Model (x, y, ψ)

  ∂Fyi ∂λyi and the axis of F yi 2) Longitudinal dynamics and transient phenomenon: By analogy, the notion of relaxation length is used to describe the longitudinal dynamics. In [Clover 98], the authors present the variations of the slip rate by a first order differential equation. They use this representation and a longitudinal linear model to study the stability of the automotive dynamics in different rolling situations. In longitudinal case the variation can be represented by a first order model: a) during braking , i = f, r σ xi .

Fig. 3 .

 3 Fig. 3. Vehicle Model for estimation of road profile

  + H 11 sign(x 22 ) . x12 = x22 + H 12 sign(x 22 ) . x21 = -M -1 1 (C 11 x21 + C 12 x22 + K 11 x11 + K 12 x 12 ) +M -1 1 ζ + H 21 sign(x 31 ) . x22 = -M -1 2 (C 21 x21 + C 22 x22 + K 21 x11 + K 22 x 12 ) +H 22 sign(x 22 ) (30) H 11 , H 21 ∈ 4×4 and H 12 , H 22 ∈ 3×3 represent positive diagonal gain matrices. B. Convergence analysis Let us define the state estimation errors xi = x i -x i (i = 1..4), then we can write the dynamics estimation errors as follows: -H 11 sign(x 22 ) . x12 = x22 -H 12 sign(x 22 ) . x21 = -M -1 1 (C 11 x21 + C 12 x22 + K 11 x11 ) -M -1 1 ζ -H 21 sign(x 31 ) . x22 = -M -1 2 (C 21 x21 + C 22 x22 + K 21 x11 ) -H 22 sign(x 22 ) (31)

  22 H 22 sign(x 22 ) Since the states are bounded and while choosing the matrix H 22 components (H i2 , i = 1...3) such that H i2 > |C 21 x21 + K 21 x11 | , the equation (33) becomes: V1 = -x T 22 M -1 2 C 22 x22 -xT 22 H 22 sign(x 22 ) < 0 (34) Therefore, the surface x22 = 0 is attractive and we have the convergence of x22 towards x 22 in finite time t 0 . Then according to (31), we have . x12 = 0 and consequently . x11 = 0. Then we deduce the unknown vector ζ such as: ζ = ζ -ζ = C 11 x21 + K 11 x11 + M 1 H 21 sign(x 31 ) (35) Finally, we get the variable ζ: ζ = ζ + M 1 H 21 sign(x 31 )

Fig. 4 .

 4 Fig. 4. Estimated and Measured States

Fig. 6 .

 6 Fig. 6. Estimated and Measured Forces

Fig. 7 .Fig. 8 .

 78 Fig. 7. Longitudinal Profile Analyser (APL in french)

  x 11 x 31 -C x r f x 21 ; 0; x 11 ; 0 0; C x x 11 -x 11 x 32 -C x r r x 22 ; 0; x 11 ; Ψ 2 (X) = C y l f x 13 +C y x 12 -x 33 x 11 ; 0; x1 11 ; 0 0; C y x 12 -C y l r x 13 -x 34 x 11 ; 0; x 11


M ' Sirdi

THE VSAS APPROACH GIVES THE BEST MPPT FOR SOLAR ENERGY SOURCES (RES)

Keywords: Maximum Power Point Tracking, Renewable Energy Sources, Perturb and Observe, IncCond, Hill Climbing

More and more MPPT (Maximum Power Point Tracking) algorithms are in competition to maximize energy extracted from PV systems. This paper shows how to get the best algorithm (the most simple, fast and robust). VSAS (Variable Structure Automatic Systems) control methodology is applied to develop the control algorithm, to clarify the rationale behind and get the best optimization algorithm.

I. INTRODUCTION

The PV system operation principle needs a polarization depending on the weather to fix the operation point leading to extraction of the maximum power. The Renewable Energy Systems (RES) include commutations and discontinuities; this is one kind of Variable Structure System (VSS). However, the behaviour of the conversion systems of this kind of renewable energy is VSAS (Variable Structure Automatic Systems)[MsirdiEFEA14] and highly dependent on variations in climate parameters, such as temperature and irradiation.

The MPPT algorithms are necessary to maximize, at each time instant, the produced power. Several techniques have been designed to search this optimal Maximum Power Point (MPP). In the literature, more and more MPPT (Maximum Power Point Tracking) algorithms are in competition to maximize energy extracted from PV systems.

The maximum performance of a photovoltaic system depend, of course, on good weather conditions, but needs also appropriate MPPT algorithm [Mutoh, NianChun]. The great majority of MPPT control strategies are based on the (steady state) characteristics of PV panels, such as I-V or P-V plots, the duty cycle ratio control and sometimes using look-up tables [Amei]. A lot of MPPT techniques are well established in the literature. There are several methods: voltage feedback method, perturbation and observation method, linear approximation method, incremental conductance method, hill climbing method, actual measurement method, fuzzy control method and so on [Ting-Chung, Tavares, Hua Lin, Fangrui, ChihChuanHua, Weidong]. In general, there exist four types of MPPT techniques:

 the PV operation point perturbation and observation (PO) based algorithms in order to get the direction of tracking the MPP.

 the hill-climbing algorithm which makes a perturbation in duty cycle to reach the apex of the characteristics.

 the incremental conductance (InCond) algorithm which periodically checking heslope (conductance) of the P-V curve [Liu].

 the constant voltage algorithm based on

keeping constant the ratio between the PV voltage at the maximum power and the open circuit voltage (Voc) value; In this method the effect of solar irradiance variations is neglected [Hohm].

In[IbrahimHoussiny, Amei], Ibrahim and Houssing use a look-up table to track, when other author prefer the use of a dynamic MPP tracker to PV appliances [MidyaKerin]. A single-stage MPP controller using the slope the power versus voltage, like has been done in [KuoLiangChen]. In [Hua Lin] a DSP chips is used to implement the PO MPPT in order to get maximum power. They, also try to improve the efficiency of the PO and HC methods.

In [Fangrui], the time response of PO [ChihChuanHua]and HC methods [Fangrui, Weidong, Hohm]are compared for a grid connected system. The PO method can fail under rapidly changing atmospheric conditions. Several research activities have been carried out to improve the traditional Hillclimbing and P&O methods. Reference [Xiao] proposes to use three measurement points to compute MPP.

The method compares the obtained power measurement to the two preceding points before choosing the perturbation sign. In [KuoLiangChen] the authors propose a two stage algorithm that offers faster tracking in the first stage and more accurate tracking in the second stage. To prevent divergence from MPP, they use a modified adaptive algorithm. The paper is organized as follows. Section 1 gives the problem formulation and recalls most used algorithms used for MPPT. The second section presents PV system equations and features. Section three introduces our new algorithm design approach. We compare their results with the widely used MPPT algorithms; performance is evaluated considering different actual solar irradiation measured variations. The fourth section gives a conclusion and proposes perspectives.

The Solar irradiation Model

The electric power operating by a PV panel depends on the irradiation intensity of the sun shining, and the temperature of the PV cells. The simplest solar radiation, perturbation free, model that can be used is the following, from (roughly) 7h-19h, the period of sunshine. At the top of the sinusoid, the maximum power is assumed to be 1kw/m², at the time 13 hours. The sunless period is 19h-31h [schijndel1, Sallem]. Note that the wave is assumed piecewise sinusoidal. The solar radiation is represented by the following expression:

G 𝑆 = {G 𝑆-𝑝𝑖𝑐 sin ( 2πt 3600 ) if t ∈ [ 7h, 19h ] (1)
The electric diagram is equivalent to an average PV cell as shown in figure 1 (a)

Photovoltaic Systems Model and characteristics

A PhotoVoltaic Generator consists of a group of PV modules electrically connected in series or parallel or series-parallel combinations with each other in a RES to generate required currents and voltages

[MsirdiEFEA14, Liu].
The model is a photo-current source IL one diode with as reverse saturation current I0 , and a serial resistance RS , representing the PV cell resistance.

The circuit is connected to the load ( R L ) trough a converter in order to adjust (adapt) the operating voltage and current of the PV panel at optimal values to maximize the harnessed power and transmit it. The control has to tracks the Maximum Power Point.

The equations describing the I(V) relationship between the current and voltage of a solar cell are given by ( 2), with the parameters defined in Table1 [Shaefer, Hsiao, Walker].

𝐼 = 𝐼 𝐿 -𝐼 0 [exp ( 𝑞(𝑉-𝑅 𝑠 𝐼 𝑛𝐾 𝐵 𝑇 𝑠 ) -1] - 𝑉+𝑅 𝑠 𝐼 𝑅 𝑃 (2)
𝐼 0 = 𝐼 0 (𝑇 1 ) [exp ( 𝑉 𝑂𝐶 (𝑇 1,𝑟𝑒𝑓 ) 𝑛𝐾 𝐵 𝑇 1,𝑟𝑒𝑓 /𝑞 ) -1] (3)
The internal serial resistance RS relative to one cell in open circuit voltage VOC is:

𝑅 𝑆 = - 𝑑𝑉 𝑑𝐼 ⎥𝑉 𝑂𝐶 -( 𝑛𝐾 𝑜 𝑇 1,𝑟𝑒𝑓 𝑞 ) / 𝐼 0 (𝑇 1,𝑟𝑒𝑓 ) exp ( 𝑉 𝑂𝐶 (𝑇 1,𝑟𝑒𝑓 ) 𝑛𝐾 𝐵 𝑇 1,𝑟𝑒𝑓 /𝑞 ) (4)
The PV output voltage of the serial cells can be expressed as V:

𝑉 = 𝑁 𝑆 𝑛𝐾 𝐵 𝑇 𝑞 𝑙𝑛( 𝐼 𝑆𝐶 -𝐼 𝑞𝐼 0 ) (5)
IL varies as a function of the ambient temperature T and of the solar radiation G as:

𝐼L= 𝐼L(𝑇 1,𝑟𝑒𝑓 )+ 𝐾 𝑜 (𝑇 -𝑇 1,𝑟𝑒𝑓 ) (6)
With

𝐾 𝑜 = 𝐼 𝑆𝐶 (𝑇 2,𝑟𝑒𝑓 )-𝐼 𝑆𝐶 (𝑇 1,𝑟𝑒𝑓 ) 𝑇 2,𝑟𝑒𝑓 -𝑇 1,𝑟𝑒𝑓 (7)
𝐼L= 𝐼STC(𝑇 1,𝑟𝑒𝑓,𝑆𝑇𝐶 )

𝐺 𝑆 𝐺 𝑆𝑇𝐶 [START_REF] Liu | Comparison of P&O and hill climbing MPPT methods for grid-connected PV converter[END_REF] The PV system exhibits a nonlinear I(V) characteristic which depend on the temperature and the solar radiation which vary during a day, for exemple, as shown by the figure TempRadiation; where a) and b) show the profile of real irradiation measured and temperature for one day.

The photovoltaic Panel considered, in this paper, contains 36 photovoltaic cells connected in series. To track the Maximum Power Point (MPPT) and achieve the optimum matching, a good control of the DC-DC is necessary. As the temperature, insulation and load vary, an algorithm is used to ensures that the PV module always operates at its maximum power point. The boost converter operation principle is the one of Variable Structure System (VSS). The Switching frequency of the Boost IGBT (changing s) is in general around fsw = 1/Ts =20kHz . We will try to use, in the comparative analysis different frequency values.

B. The MPPT Controls for PV Systems

This part reviews 3 of the most frequently used control algorithms to get the maximum power. As presented previously, an MPPT controller is used to increase the PV system efficiency. The principle is to calculate the optimal reference output voltage and/or current which ensures that the PV system operate at its MPP. These techniques are different from each other in many aspects, including simplicity, convergence speed, hardware implementation, sensors required, cost, range of effectiveness and need for parametrization (see eg [Shaefer] for a survey).

Perturb and Observe Methods

The most commonly used MPPT algorithm is the Perturbation and Observation (PO) due to its easy implementation. It uses the P-V characteristics Ppv=f(V) of the PV module. Note that, for constant weather conditions, the operating power point P(n)=V(n)I(n) is obtained when the condition dP/dV=0 is accomplished. We can calculate the slope (dP/dV) using consecutive outputs measurements (voltages and output): The duty cycle of the Boost is changed and the process is repeated until the maximum power point has been reached. In actual experiments, the system oscillates around the MPP. To minimize the oscillations amplitude, we can reduce the perturbation step size. However, small step size slows down the convergence of the MPPT.

To solve this problem, we can use smaller perturbation size towards the MPP.

Perturb and Observe Methods

The incremental conductance (IncCond) [Femia], method is based on the fact that the slope (or the PV conductance G= dI/dV ) of the PV array, in the power curve is zero at the MPP and it is positive (constant) on the left of the MPP. The slope becomes negative on the right of th MPP. These relations can be rewritten in terms of the array current and voltage as: (this means that 𝑑𝑃 𝑑𝑉 > 0 ⁄ ), the operating point is on the left side of the MPP, V has to be raised.

 If G = 𝑑𝐼 𝑑𝑉 < - 𝐼 𝑉
(dP/dV<0), the operating voltage is on the right side of the MPP then V has to be reduced.

voltage is on the right side of the MPP then V has to be reduced.

Hill Climbing Method

In general the almost applications, use DC-DC converters and DC-AC inverters, as the power interface devices between PV modules and loads. The basic idea of the HC (Hill Climbing) method is the same as P&O method. It tests if P(n) is greater than P(n-1) or not, to reach MPP. The PO method uses instead a test on dP/dV to determine whether the maximum power point has been found or not. However, the HC method uses a test condition on P(n)-P(n-1) and uses the duty cycle (D) of these switching mode power interface devices as the decision action parameter for the maximum power point tracking. The figure7 shows a flow diagram of the hill climbing algorithm.

The Proposed VSAS-MPPT Algorithms

The desired objective to get is that the MPP reached when the maximum power is obtained (P(t)=Pmax and dP/dt=0). The power is function of the voltage V, the current I and the time t, then the required Maximum Power Point to Track is really defined by the following objective function:

𝐼 𝑃𝑉 ∆𝑃(𝑘) = 𝐼(𝑘). ∆𝑉(𝑘) + 𝑉(𝑘). ∆𝐼(𝑘) ∆𝑃(𝑘) = 𝐼(𝑘). 𝑢 1 (𝑘) + 𝑉(𝑘). 𝑢 2 (𝑘)

The two control variables are then u1(k)= ∆V(k) (the voltage variation) and u2(k)=∆I(k) (the current variation).

RUCA: Robust Unified Control Algorithm

For the proposed RUCA algorithm, both control inputs can be used if we look for adjusting both variables (V and I), either at each control step or alternatively. It can be noticed that the previous algorithms can be considered as particular cases of this one, when simplifying the proposed control method. P max is a constant chosen greater than the maximum power in any weather condition Pmax≥P(t) . Let us then consider the Lyapunov like function  Calculate the PV current and PV power increments (see equation eq:optim).

W(t) = ( 𝑃

MEPO: Modified Enhanced PO Algorithm

If we take as input u1(k) and put a constant current u2(k)=0 , which means that the voltage is perturbed and the current is fixed u2(k)= ∆I(k)=0 , we are in the same configuration as the PO algorithm. We propose, as a Modified Enanced PO Algorithm which will be more robust, the reference voltage is given by Vref=Vk+ Kα∆P(k)sign(∆V(k)) Knowing that we impose ∆I(k)=0 , we get for the proposed MEPO control algorithm we take u1(k)= Kα∆P(k)sign(∆V(k)) and u2(k)= ∆I(k)=0 . The reference voltage Vref(k) is calculated as below, where α is a gain weighting the perturbation variation step. Note that ΔVref= αsign(∆P(k)∆V(k)) produces exactly the same result as the classical PO algorithm with a much simpler implementation. This method gives an enhanced variable step size algorithm. The step size is adjusted in proportionally to the power variation produced in the previous step. The adaptive step adjustment gain K is used with ∆P(k), for weighting the variation (Voltage perturbation) step. It is useful for oscillation avoidance, fast convergence and noise sensitivity reduction. We can also use only the control variable u2 (k)= ∆I(k) (and put u1(k)=0 ), the current is perturbed and the voltage is fixed ∆V(k)=0.

C. Comparative Simulations and Tests

The simulation was performed under Psim software as shown in figure 3. The physical model of a PV panel is used. The panel is considered to have 36 cells. A boost converter is built using a MOSFET as a switch. The load is a 100 Ω resistor. The algorithms are implemented in a C block and the duty cycle is calculated from Vref using another C block. The actual, measured irradiation and panel temperature are read from a txt file as inpt to the simulation.

MPPT Alghorithms Comparison

In this study we compare four MPPT algorithms under the same experimental conditions. The first algorithm called P&O stands for perturb and observe. The second algorithm called IncCond stands for incremental conductance. The third algorithm called MEPO stands for Modified Enhanced Perturb and Observe (from VSAS-MPPT). And the fourth algorithm called RSMCA (from VSAS-MPPT approach) stands for Robust Sliding Mode Control Algorithm.

Note that RUCA algorithm is simply combining use of the two VSAS MPPTs alternatively.

The comparison is made under changing values of irradiation and temperature. Two cases are considered: first we update the irradiance and temperature each time without interpolation to create sudden variation in power output, second we interpolate the values of irradiance and temperature to create quasi-continuous irradiation and temperature. At the end we use real measured data for simulation. In this study we focus on the PWM commutation frequency of the MPPT algorithm.

Zero Order Hold Interpolation Analysis

We perform many simulations with varying irradiance and temperature. The values are not interpolated causing sudden variation in power output. For high MPPT algorithm frequencies (25Hz, 50Hz, 100Hz) no difference in power output is noted between the four algorithms. After a change in environmental conditions that makes an increase in power output, MEPO algorithm may present an oscillation because it is based on the value of DP and not on the sign of DP. At high irradiances, the RSMCA algorithm present oscillation caused by the gain in calculation of Vref.For low frequencies (5Hz, 10Hz, 12.5Hz), the MEPO and the RSMCA algorithms excel the PO and IncCond algorithms. We can see in figure 2 that the PO power and the IncCond power are too far from the maximum power that can be generated from the panel. This is due to the fixed or small step size of the Vref. We must recall that the MEPO and the RSMCA algorithms also converge rapidly to the maximum power after the start of the simulation.

First Order Interpolation Analysis

When interpolating data for irradiance and temperature (approaching real conditions), we can see that that the difference between PO and IncCond algorithms and MEPO and RSMCA algorithms are emphasized when we have high variation of environmental conditions and low algorithm frequency.

Real Data Analysis

In this simulation we take real data measured for one day. We choose 20 minutes that present high fluctuations of irradiance ( 878W/m2 -126W/m2 ) and temperature ( 26°C -14°C ). The data is available each 5 min which mitigates the high variation in environmental conditions. This is why the algorithms results are almost identical even for low frequency.

III. CONCLUSION

In this paper, a new Robust MPPT algorithm is proposed to optimize the production of a photovoltaic (PV) chain. We developed a new technique to design the maximum power point tracking (VSAS-MPPT) based on Variable Structure Automatic Systems approach. This approach generalizes the sliding mode control to systems with commutations in a simplified technique, easy to implement. New and very efficient algorithms are proposed from VSAS-MPPT approach:

 the MEPO (Modified Enhanced Perturb and Observe)

 the RUCA algorithm (Robust Unified Control Algorithm).  The MPP is computed online using a very simple algorithm which uses two control inputs (one in current and one in voltage). The other algorithms like Perturb and Observe (PO), Hill Climbing, Incremental Encoder (InCond) and SMC look as particular cases of the proposed Algorithm called RUCA (Robust Unified Control Algorithm).The proposed MPPT has several advantages: simplicity, high convergence speed, and is independent on PV array characteristics. We study and compare execution efficiency for the proposed VSAS-MPPT algorithms to the other methods, including Perturbation and Observation (P&O).

Realistic simulations are presented to show ease of implementation of our new algorithms, and compare its efficiency and accuracy to other MPPTs. The obtained results prove that the MPPT is tracked even under sudden change of irradiation level or temperature. In our comparative tests, a particular case (use only voltage input), called the Modified and Enhanced Perturb and Observe Algorithm (MEPO), is shown to be the most fast and efficient, despite using low frequency commutation and sudden fast temperature and irradiation changes. The proposed MPPT has several advantages: simplicity, high convergence speed, oscillation free and is independent on PV array characteristics. In summary the best algorithms are those designed using the SASV-MPPT approach and considering that the system can move from one characteristic to another. The proposed algorithms are the most efficient despite using low frequency commutation. They are the faster converging. The simulations was performed under PSIM software to use realistic physical models.

Fig . 1 .

 1 Fig .1. Equivalent circuit model of PV panel

Fig . 2 .

 2 Fig .2. a) Measurement of Solar irradiance taken for a day; b) Temperature taken for a day

Fig . 3

 3 Fig .3 Shows the electric circuit diagram of a PV panel connected to a DC/DC converter with the MPPT control circuit based on PO.



  If the operating voltage of the PV array is perturbed in a given direction and DP/DV>0, then the perturbation moves the array's operating point toward the MPP. The PO algorithm would then continue to change the PV array voltage in the same direction.  If DP/DV⟨0 then the change in operating point moves the PV array away from the MPP, and the PO algorithm reverses the direction of the perturbation [Femia]. Some limitations are encountered such as: oscillations around the MPP in steady state operation, slow response speed, and even tracking in wrong way when atmospheric condition instability [Kim, Liu, Femia].

  /DV = 0 we need: G=dI/dV=-I/V; DP/DV = 0 at the MPP  If G =

  Let us consider the control in case of discrete time, then the fetched MPPT may be defined by ΔP(k)=0, { ∆𝐼 = 𝐼 𝑃𝑉 (𝑘) -𝐼 𝑃𝑣 (𝑘 -1) ∆𝑉 = 𝑉 𝑃𝑉 (𝑘) -𝑉 𝑃𝑉 (𝑘 -1) ∆𝑃 = 𝑃 𝑃𝑉 (𝑘) -𝑃 𝑃𝑉 (𝑘 -1) 𝑃 𝑃𝑉 (𝑘) = 𝑉 𝑃𝑉 (𝑘).



  2 𝑚𝑎𝑥 -P(t) 2 ) > 0 which is strictly positive everywhere ( ∀t, ∀I, ∀V ). The derivative of the proposed Lyapunov function W(t)=( P^2〗_max-〖P(t)〗^2 )>0 , is 𝑊 ̇(t) = -P(t) 𝑑𝑃(made negative by choosing the appropriate control laws u1(k)= ∆V(k)and u2(k)=∆I(k) to get a decreasing Lyapunov function. The RUCA uses two control inputs which can be done either simultaneously or alternatively or one of them can be frozen depending on the needed voltage or current. Choosing u1(k)= Kα∆P(k)sign(∆V(k)) u2(k)=K_2 α_2 ∆P(k)sign(∆I(k)) ensures that W (t)<0 and then the convergence of the RUCA algorithm. This proves, theoretically also the convergence of the MEPO algorithm which follows. We can use both inputs or only one of them either, in a first case u1(k)= Kα∆P(k)sign(∆V(k)), or in the second case u2(k)=K_2 α_2 ∆P(k)sign(∆I(k)). In control context, the previously presented MPPT controllers use only one control variable u1 or u2 and impose the second to be zero. The implementation of the proposed enhanced MPPT controller RUCA can be summarised as follows: The reference voltage is set be equal to the double of the PV open circuit voltage.  Measurement of the of input signals (PV voltage, RESD © 2015 http://apc.aast.edu PV current and Load voltage).  Estimate the PV power at the sample time k : PPV (k)= Ipv(k).Vpv(k)

Fig . 4 .

 4 Fig .4. Shows the electric circuit diagram of the 4 MPPT algorithms

Fig . 5 .Fig . 7 .Fig . 8 .

 578 Fig .5. Sudden variation of irradiation and T with frequency (100Hz)

Table 1 .

 1 Parameters Definition

	Symbol	Units	Definition
	q	q=1.6.1 0 -19 C	The electron charge
	n	diode non-ideality factor
	K B	1.38065•1 0 -23 N.m/K	the Boltzmann's constant
	STC	at (1kW/m2, 25°C)	Standard Test condition
	I SC	A	Short Circuit Current
	I SC,STC	A	Short circuit current at STC
	V OC	V	Open Circuit Voltage
	K 0	% calculated for ISC	Temperature coefficient
	G	W/m2	Solar radiation
	T	in Kelvin °K	the cell temperature
	T 1,ref	at 1kW/m2, 298K (25°C)	Reference temperature
	T 2,ref	K	Temperature at 2nd STC
	V ref	Vref =5V	Reference voltage
	In [Tina]Tina et. al. proposed a mathematical model
	for the electrical-thermal coupling of a PV module
	with ambient temperature, wind speed, wind
	direction, relative humidity and electrical operating
	point (voltage and current values). An MPP tracking
	based on Dual boost converter is developed using
	fuzzy logic in [Veerachary]. Artificial neural network
	(ANN), trained offline with a gradient descent
	algorithm using a back-propagation have been used
	by Kaiser et. al.

June 2015 -ISSN 2356-8569

RESD © 2015 http://apc.aast.edu


Kamel Sia

A Namaane

N K M'sirdi

Aziz Naamane

N M'

Structuring residuals by faults structural decoupling from bicausal diagnostic bond graph

Keywords: Component Faults Detection and Isolation (CFDI), isolation enhancement, structured residuals, structural analysis, bond graph, diagnostic bond graph (DBG), fault decoupling, unknown input decoupling, analytical redundancy relations (ARRs)

This paper presents a new graphical algorithm to carry out the rule of structured residuals, in order to enhance the internal component faults isolation of linear and linearized dynamical systems. For this aim, it is shown that the diagnostic bond graph; which allows generating residuals by using the ARRs scheme also allows handling their structuration. The residuals structuration is achieved through the fault structural decoupling leading to selective sensitivity to subsets of faults. To deal with the fault structural decoupling, a new bond graph bicausal element is introduced; it is termed Effort and Flow Propagator (or join) element. This approach is within the frequency domain, so, the process of residuals generation is homogenized by using solely the rule of the transfer function in s-form.

INTRODUCTION

The problem addressed here is the isolation enhancement of the internal component (or plant) faults by performing the fault structural decoupling, from linear and linearized bond graphs [START_REF] Borutzky | Bond Graph Methodology -Development and Analysis of Multidisciplinary dynamic system models[END_REF][START_REF] Karnopp | System dynamics: A unified approach[END_REF][START_REF] Rimaux | Study of structural properties of specific classes of physical nonlinear systems modeled by bond graph[END_REF]. Regarding the output detectors (or sensors) and the power sources (or actuators), they are supposed either working normally, or they have their own FDI systems [START_REF] Chen | Robust Model Based Fault Diagnosis for Dynamic Systems[END_REF][START_REF] Staroswiecki | Analytical redundancy relations for fault detection and isolation in algebraic dynamic systems[END_REF]. In the considered problem, the complexity is due to the large number of internal components of dynamical systems w.r.t. that of output detectors. Face this complexity, thanks to the powerful rule of the structured residuals [START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF]; the fault isolation can be enhanced. It is then possible to generate the primary residuals, directly from the rearranged original (or initial) model. And in order to enhance the faults isolation, the primary residuals are transformed, by structuring them to be sensitive to specific faults, while being insensitive to others. So doing so, the effects of faults upon the residuals become distinguishable, and the faults can be isolated.

The task of residuals structuration is carried out off-line during the phase of residuals generation, it is performed by means of faults decoupling [START_REF] Chen | Robust Model Based Fault Diagnosis for Dynamic Systems[END_REF][START_REF] Ding | Model Based Fault Diagnosis Techniques: Design Schemes, Algorithms and Tools[END_REF][START_REF] Frank | Handling modelling uncertainty in fault detection and isolation systems[END_REF][START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF][START_REF] Ge | Detection of faulty components via robust observation[END_REF], and the resulted residuals are termed secondary residuals [START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF][START_REF] Korbicz | Fault diagnosis[END_REF].

For the aim of internal component FDI problem formulation, figure [START_REF] Benmoussa | Bond graph approach for algebraic component fault detection and isolation[END_REF] recalls the arrangement of the basic elements of the bond graph model, where the Simple Junction Structure (SJS) is chosen [START_REF] Thoma | Simulation by bond graph: Introduction to graphical method[END_REF]. It is easy to understand that the internal physical components of dynamical systems are mapped in the following discrete set { }; where each basic element has a specific treatment of the incoming energy. Indeed, these elements share a common property; their behavioural constitutive laws involve similar quantities; i.e. an output variable, an input variable and a physical parameter.

It is assumed that parametric (or multiplicative) as well as structural (or additive) component faults acting, in a similar way, as unknown inputs supplied by power sources, see figure 1. So then, a new fault structural decoupling algorithm is proposed to generate sufficient secondary residuals to enhance the component faults isolation. It is based on the annihilation of the influences of non-isolated faults on variables of the DBG that furnish ARRs, and its main tool is the bond graph bicausality concept [START_REF] Gawthrop | Bicausal bond graph[END_REF]. In graphical models based FDI literature, there are researches that have addressed the problem of component fault isolation enhancement by fault decoupling. In [START_REF] Benmoussa | Bond graph approach for algebraic component fault detection and isolation[END_REF], assuming the modeling of component faults as unknown inputs, thus the inverse bond graph model is used to allow output detectors to the specified component faults inversion, consequently the faults can be estimated and their effects can be decoupled. Noting that, here the computational derivative causality is preferred. In [START_REF] Gahlouz | Unknown inputs observers for component faults estimation: Bond graph approach[END_REF] and [START_REF] Tarasov | Unknown inputs estimation for control and diagnosis: Bond graph approach[END_REF], still with the same fault hypothesis modeling as previously, but by using the Unknown Input Observer (UIO) scheme from bond graphs; where for each component fault, a single UIO is dedicated to its estimation. Also, for the same issue of fault isolation enhancement, but this time by using the bipartite graph model [START_REF] Staroswiecki | Structural analysis[END_REF], deep knowledge about the time-dependency of faults (abrupt fault, incipient fault) are integrated in faults modeling to get further ARRs allowing decoupling the ffects of non-isolated faults [6] [7]. In this case, the isolability matrix is used to identify which faults require further modeling effort.

Moreover, the bond graph computational tools used in this paper are shared by some researches, but for the purpose of ARRs generation; so that in [START_REF] Samantaray | Model based supervision: A bond graph approach[END_REF], the bicausality benefit; which is separating the causality assignment to the effort and flow variables of the same power bond, is employed especially to handle nonlinear constitutive laws singularities. Likewise as in [START_REF] Borutzky | Bond Graph Methodology -Development and Analysis of Multidisciplinary dynamic system models[END_REF], the preferred integral causality is assigned to the diagnostic bond graph to avoid computing the unknown variables arising after the construction of the fault indicators (or residual sinks). This paper is organized as follows: Section II introduces the four hydraulic tanks system benchmark, through which the development of the CFDI algorithm is processed. It provides the bond graph model and it models the component faults to get the faulty behavior. Section III generates the primary residuals and particularly it points out their insufficiency to ensure the component faults isolation. Section IV presents the core of this paper, that is, the generation of secondary residuals by structural decoupling the primary residuals from non-isolated faults, where the bicausal diagnostic bond graph is the fundamental used tool. Section V runs the simulation example to prove the theoretical developments. Finally, section VI concludes with some benefits of the faults structural decoupling from the bond graph modeling.

BENCHMARK SYSTEM TO CFDI SYSTEM DESIGN

The nonlinear four hydraulic tanks system proposed for CFDI system design is shown in figure [START_REF] Borutzky | Bond Graph Methodology -Development and Analysis of Multidisciplinary dynamic system models[END_REF]. The tanks ; with , have the same section and they are coupled by nonlinear valves. The flow rate through the valve is given by Torricelli law:

( ); with is in ( ⁄ ), and are liquid levels in ( ) on either side of the valve and .

The input of the hydraulic system is the inflow :

( ) ( )
The linearized bond graph model around the equilibrium point, denoted , is shown by figure [START_REF] Brown | Engineering System Dynamics, a Unified Graph-Centered Approach[END_REF]. It should be noted that, the prime symbol standing for linearization is omitted for equations clarity. The tanks are represented by linear storage C-elements with parameters, ; with , but the nonlinear valves become linear and they are represented by R-elements with the same parameter ⁄ .

Let the complementary state variables at the storage Celements in integral causality be denoted ;

, and they represent liquid levels in the corresponding tanks:

( ) ( )
The output variables are furnished by effort detectors ; with , and they are defined as follows:

( ) ( ) ( )
While, the input variable is supplied by the unique flow power source .

Only for the clarity of the bond graph model, the state initial conditions, , are not represented. Still, it should be noted that, they can be modeled as power sources acting on the output variables of storage elements, and , in integral causality by means of adequate simple junctions [START_REF] Borutzky | Bond Graph Methodology -Development and Analysis of Multidisciplinary dynamic system models[END_REF][START_REF] Gawthrop | Bicausal bond graph[END_REF].

Internal component (or plant) faults modelling

In the presence of the internal component (or plant) faults, the bond graph model of the faulty behaviour is depicted by figure 4. It displays the parametric faults, and ; with , they are modelled by additive increments added to the nominal parameter, and , and they represent changes in physique properties of the components [START_REF] Samantaray | Model based supervision: A bond graph approach[END_REF]. Moreover, it is taken into account the structural faults that are modelled as flow sources, ; with , they acting as inputs to the storage elements and they represent leakage in different tanks [START_REF] Samantaray | Model based supervision: A bond graph approach[END_REF]. The individual treatment of each fault of each component makes the problem of faults isolation very complex to formulate [START_REF] Ding | Model Based Fault Diagnosis Techniques: Design Schemes, Algorithms and Tools[END_REF][START_REF] Ge | Detection of faulty components via robust observation[END_REF]. To overcome this complexity, the components are regrouped in a set of equivalence classes, where each equivalent class involves bond graph elements linked by causal loop of order 1 with a storage element in integral causality (in this case, elements), and it is recognizable by the label of this storage element. The table [START_REF] Benmoussa | Bond graph approach for algebraic component fault detection and isolation[END_REF] highlights the different found equivalence classes.

C C   : C 2 2 C C   : C 3 3 C C   : C 4 4 C C   : C 2 2 R R   : R 3 3 R R   : R 4 4 R R   : R 1 in,c f : Sf 2 in,c f : Sf 3 in,c f : Sf

Storage element in integral causality

Equivalence classes

Sets of basic elements

[ ] { } [ ] { } [ ] { } [ ] { }

Table 1. Equivalence classes of components

To each storage element in integral causality, ; with , a set of component faults is associated, it includes at once the parametric and structural faults of the components belonging to its own equivalence class and they are recognized as its own faults, denoted ; with . The table [START_REF] Borutzky | Bond Graph Methodology -Development and Analysis of Multidisciplinary dynamic system models[END_REF] highlights the component faults associated to the storage elements.

Faults of storage element in integral causality

Sets of component faults associated to storage element in integral causality

{ } { } { } { }

Table 2. Faults associated to storage elements in integral causality

In the end of the component faults modelling task, the simplified faulty linearized bond graph model is depicted by figure [START_REF] Ding | Model Based Fault Diagnosis Techniques: Design Schemes, Algorithms and Tools[END_REF], where the possible component faults associated to the storage elements in integral causality are modelled as unknown inputs supplied by flow power sources ; with , and they acting at junctions to which the storage elements in integral causality are attached.

Remark 1. For clarity, henceforth, the component faults associated to the storage elements in integral causality are simply called component faults and the flow power sources representing them are simply called component fault sources.

Since the linearized bond graph model, , contains four storage elements in integral causality; there are four component faults to be detected and isolated; which are regrouped and put in the following component faults vector:

( ) ( )
From now, the problem of component faults detection and isolation from linear and linearized bond graphs is posed and the rule of structured residuals is chosen to ensure the fault isolation.

PRIMARY RESIDUALS GENERATION

The Diagnostic Bond Graph (DBG), which is denoted here ̆, is introduced in [START_REF] Medjaher | Diagnostic bond graphs for direct residual evaluation[END_REF][START_REF] Samantaray | Model based supervision: A bond graph approach[END_REF]; it is the basis of residuals generation based on ARRs scheme. It is established by simple rearrangement of the causality in the original bond graph model , in order to generate straightly analytical constraints equal to zero [START_REF] Ould Bouamama | Derivation of constraint relations from bond graph models for fault detection and isolation[END_REF]. The causality rearrangement is obtained by substituting the output detectors ; with , by modulated power sources ̆ ; with , where the modulate signals are the measured outputs; see for this figure [START_REF] Düstegör | Structural analysis of fault isolability in the DAMADICS benchmark[END_REF].

Noting that the diagnostic bond graph, ̆, is assigned the preferred integral causality, in order to carry out the rule of the transfer function in s-form (see rule 1, in the appendix). For the sake of CFDI problem formulation [START_REF] Sia | A new structural approach for ARRs generation from linear & linearized bond graphs: Approach based on state structural observability[END_REF][START_REF] Sia | A new structural approach for ARRs generation from linear & linearized bond graphs: Approach based on the transfer function[END_REF], the new emerged modulated power sources are termed imposed output sources, and the junctions resulting after the redistribution of the causal strokes are termed null junctions and noted ; these latter are so called, because their net power are in principle identically zero. The response of the bond graph model to the solicitation of the imposed output sources are termed null signals and denoted .

Remark 2. The DBG ̆

is a generalized bond graph, due to the presence of storage elements ( and ) in derivative causality. Therefore, this property must be considered by the algorithms of the analytical models derivation.

The null junctions are labeled by the number of their original junctions from which they are issued:

{ } ( )
Whereas, the null signals are labeled by the same number of the null junctions to which they belong and they are put in a vector:

( ) ( )
Through the process of the null signals construction, the following important remark can be made.

Remark 3. It is possible to say that the null signals, in order,

, and are due respectively to the output detectors, , and .

In s-domain, by applying the rule of the transfer function in s-form, the null signals ( ) can be expressed from the input variables ̆( ) of the DBG ̆, which are specified in the equation [START_REF] Frisk | Improving fault isolability proprieties by structural analysis of faulty behavior models: Application to the damadics benchmark problem[END_REF] and that are; the input variables ( ), the imposed output variables ( ), the initial conditions ( ) and the component faults ( ), following rational transfer matrix accordingly to the equation ( 8):

̆ { } ( ) ( ) ( ) [ ( ) ( ) ( ) ( ) ( ) ( ) ( )] ( )
Where ( ), ( ), ( ) and ( ) are polynomial matrices in the Laplace variable ( ) with appropriate dimensions and ( ) is the characteristic polynomial of ̆ and it is a common denominator to all the elementary transfer functions.

To get the polynomial form of the ARRs, the two sides of the equation ( 8) are multiplied by the common denominator ( ), and after that by neglecting the effects of initial conditions , because they are transient in time [START_REF] Ding | Model Based Fault Diagnosis Techniques: Design Schemes, Algorithms and Tools[END_REF], it results then the equation ( 9):

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )
Where ( ) are the primary ARRs generated from the original DBG ̆, in the frequency domain [START_REF] Sia | A new structural approach for ARRs generation from linear & linearized bond graphs: Approach based on the transfer function[END_REF]. Their numerical values are supposed in the neighborhood of zero due to measurement noises, parameter uncertainties and linearization errors.

Let the primary residuals be defined as follows:

( ) ( ) ( ) ( ) ( ) ( )
Where, the part of the equation ( 10) involving known and measured variables is the residual computational form [START_REF] Staroswiecki | Analytical redundancy relations for fault detection and isolation in algebraic dynamic systems[END_REF], thus it can be computed on line for fault detection and isolation. While, the part of the equation ( 9) involving the component faults is termed residual evaluation form [START_REF] Staroswiecki | Analytical redundancy relations for fault detection and isolation in algebraic dynamic systems[END_REF]. , can influence an ARR that is generated from a null signal ; with , if its faults source is causally linked through at least one causal path to the null signal .

Corollary 3. Modeling the component faults as unknown inputs leads to direct separation between the computational and the evaluation forms.

In [START_REF] Sia | A new structural approach for ARRs generation from linear & linearized bond graphs: Approach based on the transfer function[END_REF], a graphical algorithm entitled the influences approach is provided to generate symbolically the polynomial form of the ARRs, by computing the influence (or contribution) due to each input variable of ̆; with ̆ { }, following the formulas [START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF][START_REF] Ge | Detection of faulty components via robust observation[END_REF]:

( ) ∑ ∑ ( → ) ( ) ( → ) ̃ ( ) ( ) ( ) ( )
Where is the primary ARR generated from the null signal , is the elementary influence due to the elementary causal path linking the input to , ̃ ( ) is the numerator of the gain ( ) of the causal path, ( ) is the characteristic polynomial of the reduced DBG to the causal path. The second sum is taken over all elementary causal paths and the first sum is taken over all the input variables, , of ̆; with ̆ { }. For the calculation of all these terms, see the appendix.

By analyzing the residual computation form, the following corollary (4) can be made: Corollary 4. A component fault is structurally detectable if, and only if, the component fault source is linked at least to one output detector.

The primary ARRs generated from the DBG ̆ of figure (6) are of number three; since there are three null signals due to three linearly independent detectors (see rule 3, in the appendix). Only their final expressions are given:

{ ( ) [ ( ) ( )] ( ) ( ) { [ ( ) ( ) ( )] ( ) ( ) ( ) { ( ) ( )
Clearly, the primary residuals , and are limited to parts of of equations [START_REF] Karnopp | System dynamics: A unified approach[END_REF][START_REF] Korbicz | Fault diagnosis[END_REF][START_REF] Lamb | Causal loops and Mason's rule for bond graphs[END_REF] involving the known input and measured output variables.

The Fault Signature Matrix (FSM) summarizing the sensitivity of the primary residuals to the component faults is illustrated by table [START_REF] Brown | Engineering System Dynamics, a Unified Graph-Centered Approach[END_REF]; where the binary value 1 (resp. 0) indicates that the considered component fault affects (resp. doesn't affect) the considered primary residual.

Faults

Residuals

1 1 0 0 0 1 1 0 0 0 0 1

Table 3. FSM of the primary residuals

The analysis of the FSM reveals that the effect of the component fault is not distinguishable from those of the faults and , since it affects both the residuals and . So, to overcome the insufficiency of the primary residuals to isolating all the component faults, the enhancement of fault isolation is required and it is ensured by the secondary residuals.

SECONDARY RESIDUALS GENERATION

There is only one component fault, , concerned by the enhancement of isolation. Thus, only one secondary residual, which is structurally insensitive w.r.t. the considered fault, is sufficient to ensure this task.

It is the fault structural decoupling that ensures the insensitivity condition of the secondary residual in respect to the non-isolated fault, and its principle from bond graphs is clarified by the definition (1). Definition 1. The principle of the fault structural decoupling is to perform causal transformations on the DBG ̆, in order to annihilate any causal link, in the sense to interdict any causal path that may cause a transfer (or influence), between the component fault source to be decoupled and the null signals that furnish ARRs.

For showing the causal transformations needed by the structural decoupling in order to make the null signals of the DBG ̆ insensitive to the component fault source , let the bicausal DBG devoted to generate the structurally decoupled secondary residual be denoted ̆ and it is written as the equation below:

̆ ̆ ( ) ̆ ( ) ( )
Where ̆ ( ) is the estimator bicausal DBG, ̆ ( ) is the decoupled bicausal DBG and ( ) is the binary causal joining operation that consists to attach, to each other, both the previous bicausal DBG. The subscript ( ) means that: (1) ̆ must be structurally decoupled from the component fault ; and (2) ̆ must furnish an estimation of the unknown complementary state variable ; whilst for a reason that will be argued in the subparagraph (4.1), the estimation must be structurally decoupled from the effect of . The subscript ( ) means that the joining operation is made through the variable .

The next subparagraphs outline the elements involved in the expression [START_REF] Medjaher | Diagnostic bond graphs for direct residual evaluation[END_REF] of the bicausal DBG for the secondary residual generation, ̆ , and through them, the causal transformations needed to the faults structural decoupling are highlighted.

Decoupled bicausal diagnostic bond graph

The decoupled bicausal DBG, ̆ , which is structurally decoupled from the component fault, , in the sense of the definition (1), is established by operating the following causal transformations on the original DBG ̆ of figure ( 6):

1. The component fault source,

, to be decoupled must be substituted by the bicausal double detectors element ; 2. The storage element in integral causality, , which is the origin of the component fault , must be substituted by the bicausal double sources element , and; 3. To satisfy the causal half-strokes distribution, the bicausality is propagated from until reaching .

At the end of the previous succession of causal transformations, it results the decoupled DBG ̆ depicted by the figure [START_REF] Frisk | Improving fault isolability proprieties by structural analysis of faulty behavior models: Application to the damadics benchmark problem[END_REF]. From the property (3), the structurally decoupled ARRs from the effect of the component fault must be generated by calculating the null signals belonging to the decoupled DBG ̆ . However, accordingly to the property (4), due to the unknown complementary state variable,

, it results that the null signals can't be calculated and they can't furnish ARRs. So, from now, the problem of unknown variables elimination is posed.

Here, the unknown variables elimination proceeds with the estimation of the unknown complementary state variable,

, by expressing it from known inputs, , and measured outputs, . But, due to the fact that is an input variable to the SJS of the decoupled bicausal DBG ̆ , thus the estimation must be structurally decoupled from the component fault , so that, its effect is not propagated through the estimated value of to the null signal that furnish the structurally decoupled ARRs.

Let the estimated value of the unknown complementary state variable be denoted and the following subparagraph is concerned by its graphical calculation.

Estimator bicausal diagnostic bond graph

The estimator bicausal DBG, ̆ ( ), which is devoted to furnish a structurally decoupled estimation of the unknown complementary state variable , is established by performing the following causal transformations on the original DBG, ̆, of the figure (6):

1. The storage element, , which is the origin of the component fault to be decoupled , is substituted by the bicausal double detectors element ; 2. Arbitrary, choosing one of the two following imposed output sources, ̆ or ̆ , which are causally linked to the complementary state variable , and substituting it by the bicausal double sources element, ̆ or ̆ . Let's choose the first one; 3. To satisfy the causal half-strokes distribution, the bicausality is propagated from ̆ until reaching . In the end, of the previous succession of causal transformations, the estimator bicausal DBG, ̆ ( ), is established and it is depicted in the figure [START_REF] Frank | Handling modelling uncertainty in fault detection and isolation systems[END_REF].

  c c f 1 : S f   c c f 2 : DeDf c c f 3 : S f   c c f 4 : S f   1 y : e S M       1   2   3   4   5   6   7   8   9   10   11   12   13   14   15   16   17   21   18   22   23   24   19   20   0 2 N s      ? : SeSf     out Q : Sf   3 y : e S M    0 3 N s   0 1 N s     1 z   2 z   3 z   4 z  
The following properties show what happened in the SJS of ̆ ( ), as a result of the structurally decoupled estimation of the unknown complementary state variable.

Property 5. The variables belonging to the estimator bicausal DBG, ̆ ( ), have the subscript ( ), and this does not change their meaning and their values. Property 6. The assigned value (!) to the bicausal double detectors element, , means that the values of its effort and flow variables depend from the internal environment of the estimator bicausal DBG ̆ ( ).

Property 7. The estimated value of the unknown complementary state variable, , is an output variable to the SJS of ̆ ( ). Property 8. The estimated value of the unknown complementary state variable, , is not causally linked to the component fault, , meaning that the estimation is structurally decoupled from the considered fault.

From the property (8), the alone condition on the estimation is hold. Indeed, when the estimation is structurally decoupled from the component fault in interest, this implies that the estimation error is consequently structurally decoupled from the same fault. Property 9. There is no null signal at the 0 1 -junction of the estimator bicausal DBG ̆ ( ).

From the property (9), the remark (4) is deduced. Remark 4. It is possible to say that the null signal belonging to the estimator bicausal DBG, ̆ ( ), is used to furnish the estimated value of the unknown complementary state variable, .

At this stage, there is one fundamental result of the fault structural decoupling from bond graphs that must be shown. For this aim, let us recall the following ideas, in the order:

1. Remark (3) stats that the output detector is used to construct the null signal ; 2. Remark (4) states that the null signal provides the estimated value ; 3. From the properties (1) and ( 5), it is possible to state that is numerically equivalent ; 4. In the end, the property (4) states that is used for the structural decoupling of the fault .

By relating the previous ideas, it is possible to conclude that, indeed, it is the output detector , which was used for the structural decoupling of the fault in interest.

From what has been said, the corollary (5) pointes out the first fundamental result of the fault structural decoupling from bond graphs.

Corollary 5. One output detector is needed to the structural decoupling of one additive component fault (or unknown input).

Corollary [START_REF] Ding | Model Based Fault Diagnosis Techniques: Design Schemes, Algorithms and Tools[END_REF] agrees with one of the fundamental result in additive faults (or unknown inputs) decoupling from analytical approaches [START_REF] Frank | Handling modelling uncertainty in fault detection and isolation systems[END_REF][START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF].

Remark 5.

It should be noted that the estimation procedure does not extend to numerically calculate the estimated value of the unknown complementary state variable, . Corollary 6. The construction of the estimator bicausal DBG requires at least one causal path linking the complementary state variable to be estimated to an output detector. And, the state structural observablity is a strong structural condition to ensure this property, since it necessitates causal link between the complementary state variables and the output detectors.

For state structural observability checking, see [START_REF] Rimaux | Study of structural properties of specific classes of physical nonlinear systems modeled by bond graph[END_REF][START_REF] Sia | Structural state observability & controllability of linear and linearized bond graphs: A Simple graph tools based reasoning approach[END_REF][START_REF] Sueur | Bond graph approach for structural analysis of MIMO linear systems[END_REF].

! : DeDf  in Q : Sf  0 2 N s  0 3 N s  1 z  2 z  3 z  4 z    
1   2   3   4   5   6   7   8   9   10   11   12   13   14   15   16   17   21   18   22   23   24   19   20   c c f 2 B    c c f 2 B    1 y : e S M   out Q : Sf   out Q : Sf ^03 N s ^02 N s   0 3 N s   0 1 N s   2

SECONDARY RESIDUALS GENERATION

The bicausal DBG for the secondary residual generation ̆ is depicted by figure [START_REF] Gahlouz | Unknown inputs observers for component faults estimation: Bond graph approach[END_REF], where at the top, it is the estimator bicausal DBG ̆ ( ), and at the bottom, it is the decoupled bicausal DBG ̆ and both they are joined by the Effort and Flow Propagator (EFP) element with gain .

The EFP element is the bond graph realization of the binary causal joining operation and it is defined below: Definition 2. An Effort and Flow Propagator element, noted EFP, is shown by figure [START_REF] Gawthrop | Bicausal bond graph[END_REF]. It is a two port and bicausal element, it has two inputs variables at its inward bond, two output variables at its outward bond and its gain is . Its constitutive law is: The EFB element is unidirectional; so it doesn't account for the feedback effect. Also, it is termed join element; this name is taken from the structural approach of ARRs generation based on relational algebra [START_REF] Ploix | A new efficient and flexible algorithm for the design of testable subsystems[END_REF], whose the joining relational operator is the fundamental tool. Indeed, the secondary residual, denoted , is the computational part of the secondary analytical redundancy relation, denoted , which must be generated from the null signal , but not from the null signal . This last one doesn't contain any term and it furnishes precisely 0; its terms reconcile with each other and they vanish all; because accordingly to remark (4), the null signal belonging to ̆ ( ) is used to furnish the estimated value of the unknown complementary state variable, .

{ ( ) k   EFP      f  e  
Here also, the secondary analytical redundancy relation, , is generated by the influences approach applied to the null signal . Consequently, the causal analysis of ̆ reveals six input variables causally linked to the null signal , that are:

̆ { ̆ ̆ ̆ } ( )
The first three input variables belong to ̆ , but the last three input variables belong to ̆ . Remark 6. Neither the component fault source , nor the component fault source , appear in the set ̆ of the input variables causally linked to the null signal . So, they have no influence on .

From the remark [START_REF] Düstegör | Structural analysis of fault isolability in the DAMADICS benchmark[END_REF], one has to agree the annihilation of the influence of the component fault on the null signal ; which furnishes in turn an insensitive secondary residual w.r.t. the same fault.

Let us show how to calculate the influence due to the input variable ̆ . In effect, there are two elementary causal paths from this input variable to the null signal , which are respectively determined by the equations ( 19) and [START_REF] Pichardo-Almarza | Bond Graph Approach to Build Reduced Order Observers in Linear Time Invariant Systems[END_REF], as successions of variables and edges and where each edge is augmented by its causal gain.

( ̆ → ) → → → → → → → → → → → → → → → ( ) ( ̆ → ) → → → → → → → → → → → ( )
The elementary influence due to the first elementary causal path is calculated as follows:

( → ) ̃ ( ) ( ) ( ) ( )
Where:

 ̃ ( ): is the part that does not contain integrators in the first causal path gain, ( ). Since, the causal path doesn't travers any storage element in integral causality, then ̃ ( ) is identical to ( ) and it is to ;  ( ): is the characteristic polynomial of the reduced bicausal DBG; which is obtained by supressing the first elementary causal path from ̆ , and it is equal to ; because the reduced bicausal DBG (not shown) doesn't contain any storage element in integral causality; Accordingly to the equation ( 21), the elementary influence due to the first elementary causal path is:

( ̆ → ) ( ) ( )
Following the same steps of symbolic calculus, the elementary influence due to the second elementary causal path is obtained and it is given below:

( ̆ → ) ( ) ( )
In the last step, by total summating the influences due to all the input variables ̆ of ̆ and after deleting the subscripts, ( ) and ( ), it results the final expression in the polynomial form of the and it is given hereafter:

{ ( ) ( ) ( )
Clearly, the secondary residual, , is reduced to the part of the , which involves only known input and measured output variables.

From the previous development, it has been stated that 1. The decoupled secondary residual is generated from the null signal , and;

2. The null signal is due the output detector ,

Then it is possible to conclude that the output detector , itself, which was used for generating the decoupled secondary residual.

At this stage, from what has been said and by taking into account the corollary (5), then the second fundamental result of the fault structural decoupling from bond graphs is pointed out by the corollary [START_REF] Frisk | Improving fault isolability proprieties by structural analysis of faulty behavior models: Application to the damadics benchmark problem[END_REF].

Corollary 7. To generate a secondary residual that is structurally decoupled from one additive component fault (or unknown input), two output detectors are necessary. The first one is used to structurally decoupling the fault (or unknown input), and the second one is used to generate the decoupled secondary residual itself.

The corollary [START_REF] Frisk | Improving fault isolability proprieties by structural analysis of faulty behavior models: Application to the damadics benchmark problem[END_REF] agrees well with the principle of additive faults isolation (or unknown inputs robustness) by means of perfect unknown inputs decoupling from the analytical approaches [START_REF] Frank | Handling modelling uncertainty in fault detection and isolation systems[END_REF][START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF].

Figure [START_REF] Gertler | Fault Detection and Diagnosis in Engineering Systems[END_REF] shows, by means of the secondary residual, another benefit of the influences approach for computer implementation of the residuals. It is the parallel realisation of the computational form, where the elementary influences due to the different elementary causal paths are taken such that as they are generated and summated to get the computer implemented form of the secondary residual; which means that there are no additional algebraic manipulations to be performed.

in Q 1 y 2 y 3 y s C R R 1 2 1  2 1 R  s C 3 2 1 R 3 1 R  s r 1 3 1 R  2 1 R R 

SIMULATION RESULTS

The simulation runs the residuals that are normalized in the interval [ ]. To handle measurement noises and input fluctuations, the residuals are followed by a State Variable Filter (SVF); with the cutoff frequency ⁄ , and where is the time response of the four hydraulic tanks system. The simulation is achieved with Matlab/Simulink software, version 7.10.0 (R2010a).

The validity of the proposed fault structural decoupling algorithm is tested through an abrupt structural fault at the junction (0 2 ); as a leakage in the tank 2, with magnitude 3% of and it is introduced at time .

The simulation results validate the developed theory. In effect, the figure [START_REF] Ge | Detection of faulty components via robust observation[END_REF] shows that after a transient behavior due to initial conditions, all the residuals tend to zero, indicating that the linearization has no effect. And, from the apparition of the considered component fault, the following observations can be deduced:

 The residuals and persist different from zero; since they are both structurally sensitive to the component fault .  The residual exhibits a transient behavior and it returns to its initial normal behavior; because it is structurally insensitive to the fault , and;  The secondary residual, , returns to its initial and normal behavior after a short transient; because it is made structurally insensitive to the considered component fault by structural decoupling. This confirms well the success of fault structural decoupling. The figure (13.a) illustrates both the structurally decoupled estimation and the simulated trajectories of the unknown complementary state variable, . While, the figure (13.b) shows the estimation error that is structurally decoupled from the considered component faults; because it is zero in steady state.

CONCLUSION

The proposed new structural algorithm of ARRs design, from linear and linearized bond graphs, in the frequency domain, for component FDI system, has achieved its objectives; so that it presents the following characteristics:

From FDI issues point of view: (1) it's structural; so it can be used at the design stage; (2) it is iterative; so it treats component by component; (3) it doesn't necessitate unknown variables elimination for the primary residuals generation, even if, for the secondary residuals it does. Fortunately, the estimation doesn't extend to numerical calculation; (4) it uses structural fault decoupling to carry out the rule of structured residuals for the aim of component faults isolation enhancement. Here, the fault structural decoupling, itself, requires structurally decoupled estimation from the component fault in interest; and (5) the fault structural decoupling is robust w.r.t. parameters uncertainties; since it is based on simple causality breaking rather than polynomial, rational or numerical matrices annihilator. From implementation issues point of view: (6) the residuals generation is homogenized by using the rule of the transfer function in s-form; and (7) the generated residuals can be implemented, such that, in parallel scheme without any additional algebraic manipulations.

Due to the advantages provided by the new structural algorithm for component fault detection and isolation, from one hand, and due to its automated calculus, from other hand, it is possible to conclude that the efficiency of the proposed structural approach is comparable to that of the analytical ones.

1. There exists storage elements ( and ) in integral causality that accept to be substituted by bicausal double detectors elements ; 2. The propagation of the bicausality, from until reaching , is satisfied, and; 3. The resulted bicausal bond graph is solvable.

The figures ( 14) and [START_REF] Lamb | Causal loops and Mason's rule for bond graphs[END_REF] show how operate the causal transformations to be performed to the output detectors and the storage elements ( and ) to carry out the rule (2). Indeed, the rule (3) point out the conditions that allow expressing the complementary state variables from the output variables. It should be noted that a quasi-similar graphical rule based on the bicausality is announced in [START_REF] Pichardo-Almarza | Bond Graph Approach to Build Reduced Order Observers in Linear Time Invariant Systems[END_REF] to check the invertible sub matrix of the observation matrix of the state space model, except that the solvability condition is not mentioned. Or, this condition is necessary to correctly performing numerical calculus.

Remark 9. If the number

of the output detectors is greater than the order of the bond graph model, then the output detectors are linearly structurally dependent. Also, there is a need to distinguish between the linearly independent output detectors from the linearly dependent ones. The following remark advises about this kind of output detectors and it is based on the violation of the conditions of the rule 2.

Remark 10. The output detectors, for which there are no state-to-output causal paths to propagate the bicausality half-strokes, they are said linearly dependents.

 q f   C q e     I p f  p e   !  f I !  e !  e !  f ! : DeDf     Figure 15.

Substitution of storage (I and C) elements by bicausal double detector elements

By applying the duality principle to the rule 3, it is possible to deduce the rule 4, which is dedicated to recognizing the linearly independent power sources. Its fundamental tools are the bicausality and the input-to-state causal paths concepts.

Rule 4. Linear independency of power sources

Let be a linear (or linearized) bond graph model, in integral causality, and it has storage elements in integral causality; i.e. its order is . It is supposed fitted by a set of power sources that is defined by the equation (35):

{ | ( )} ( )
Where stands for a power source with unspecified causality, (resp. ) stands for an effort (resp. flow) power source, is the input variable supplied by the power source and ( ) stands for equivalent to.

Then, the power sources are linearly independent if, and only if, when they are substituted by bicausal double detectors elements ; with , the following three properties hold simultaneously:

1. There exist storage elements ( and ) in integral causality that accept to be substituted by bicausal double sources elements ; 2. The propagation of the bicausality, from Until reaching , is satisfied, and; 3. The resulted bicausal bond graph is solvable.

The figures ( 16) and [START_REF] Merzouki | Intelligent Mechatronic Systems: Modeling, Control and Diagnosis[END_REF] show the causal transformations to be performed to the power sources and to the storage elements ( and ) to carry out the rule (4).

The remark (11) identifies a particular case for which the output detectors are structurally linearly dependent because of their number.

Remark 11. If the number of the power sources is greater than the order of the bond graph model, then the power sources are linearly structurally dependent.

There is a need to distinguish between the linearly independent power sources from the linearly dependent ones. The following remark advises this particular kind of power sources and it is based on the violation of the conditions of the rule 4.

Remark 12. The power sources, for which there are no input-to-state causal paths to propagate bicausality halfstrokes, they are said structurally linearly dependent.

In the end, the following concepts intervene in the rule of the transfer function in s-form, they are recalled succinctly. The interested reader may find more detail in [START_REF] Borutzky | Bond Graph Methodology -Development and Analysis of Multidisciplinary dynamic system models[END_REF][START_REF] Merzouki | Intelligent Mechatronic Systems: Modeling, Control and Diagnosis[END_REF][START_REF] Samantaray | Model based supervision: A bond graph approach[END_REF]. Definition 3. A causal path between two ports is an alternation of bonds and basic bond graph elements such that (i) all elements have a correct and complete causality, and (ii) two bonds of the path have opposite causal stroke direction. Definition 4. The causal path gain is the product of all the gains (or influence coefficients) of its edges. Where, the gain of each edge is equal to ( ) or ( ), when respectively there is or not a change of bond orientations at the (resp. ) junction, by following the flow (resp. effort) variable, but, if a bond graph element is traversed; the gain is equal to the element transmittance expressed in symbolic Laplace operator ( ).

Figure 1 .

 1 Figure 1. SJS of a faulty linear bond graph

Figure 2 .

 2 Figure 2. Four hydraulic tanks system.

4 zFigure 3 .

 43 Figure 3. Linearized bond graph model in integral causality of the four hydraulic tanks system.

Figure 4 .

 4 Figure 4. Faulty linearized bond graph in integral causality of the four hydraulic tanks system.

Figure 5 .

 5 Figure 5. Component faults modeling as unknown inputs supplied by flow power sources.

Corollary 1 .

 1 From the DBG ̆, by using the rule of the transfer function in s-form, three is a need only to compute the denominator of the null signals to get ARRs in polynomial form.

Figure 6 .

 6 Figure 6. Diagnostic bond graph ̆.

Property 1 .Property 2 .Property 3 .Property 4 .

 1234 The following properties show what happened in the SJS of ̆ , as a result of component fault structural decoupling. The variables belonging to the decoupled bicausal DBG, ̆ , have the subscript ( ), and this does not change their meaning and values. The assigned value (?) to the bicausal double sources element, , means that, from now, the values of its effort and flow variables depend from an external environment to the decoupled bicausal DBG, ̆ . The component fault, , to be decoupled becomes an output variable to the SJS of ̆ . Thus, all the null signals of ̆ are structurally decoupled from the considered fault, in the sense of the definition[START_REF] Benmoussa | Bond graph approach for algebraic component fault detection and isolation[END_REF]. The unknown complementary state variable, , becomes an input variable to the SJS of ̆ . As a result, the null signals belonging to ̆ , which are causally linked to become dependent on its value.

Figure 7 .

 7 Figure 7. Decoupled bicausal DBG ̆ from the component fault .

Figure 8 .

 8 Figure 8. Estimator bicausal DBG ̆ ( ) for the structurally decoupled estimation

Figure 9 .

 9 Figure 9. Bicausal DBG ̆ for secondary residual generation.

Figure 10 .

 10 Figure 10. Effort & Flow Propagator (or join) element.

Figure 11 .

 11 Figure 11. Parallel realization of the secondary residual.

Figure 12 .

 12 Figure 12. Trajectories of the residuals, in order, , ,

Figure 13 .

 13 Figure 13. a) Simulated and estimated trajectories of , b) Estimation error.

  Figure 14. Substitution of output detectors by bicausal double sources elements

Definition 5 .Definition 7 .

 57 A causal loop is a closed causal path starts and ends at the same element. It can exist between two elements of , , or type. Definition 6. A causal cycle is a closed causal path that can contain more than two storage ( ) elements (with any causality). Two closed causal paths are disjoint if they have neither junctions nor bonds in common while following the same type of variable. Definition 8. A causal cycle family is a set of disjoint causal cycles. The family is said to be of order , if it contains (resp. ) storage elements ( and ) in integral (resp. derivative) causality.The gains of causal loops and the causal cycles are calculated similarly as for the causal path.

Figure 16 .Figure 17 .

 1617 Figure 16. Substitution of power sources by bicausal double detectors elements.

Appendix

The rule of the transfer function in s-form [START_REF] Rahmani | Bond graph approach to structural analysis of linear systems[END_REF] uses the concept of causal cycle families, in order to write the numerator and the denominator as polynomials in the integer power of the Laplace variable ( ). It is then an alternative to the graphical rule furnished in [START_REF] Brown | Engineering System Dynamics, a Unified Graph-Centered Approach[END_REF]; which uses the Mason's loop rule to get the transfer function, but in -form. This last one writes the nominator and the denominator as sum of terms in the integer power of ( ), which needs consequently extra algebraic manipulations to get the transfer function in its classical s-form.

The following given rule is a general one, since it handles both the regular and the generalized bond graphs. It should be noted that the generalized bond graphs contain storage elements ( and ) in derivative causality. This new rule is similar to that given in [START_REF] Mouhri | Symbolic determination of generalized state equation for singular system modeled by bond graph[END_REF], but it rewrites the coefficients in a different manner.

Rule 1. Transfer function in s-form

Let be a linear (or linearized) bond graph in preferred integral causality, it may be regular or generalized and it has storage elements in integral causality; i.e. its order is . It is supposed fitted with the two sets of output detectors, , and power sources, , that are defined as follows:

Where (resp. ) stands for an output detector (resp. power source) with unspecified causality, (resp. ) stands for an effort output detector (resp. power source), (resp. ) stands for a flow output detector (resp. power source), (resp. ) is the output (resp. input) variable measured (resp. supplied) by the output detector (resp. power source ) and ( ) stands for "equivalent to".

Then, the output variable , with , can be expressed from the input variables , with , following the formula hereafter:

Where, the elementary influence, ( → ), due to the elementary causal path, ( → ), is:

With ̃ ( ) is the numerator of the gain ( ) of the elementary causal path that is determined as follows:

Where is the constant term of the gain and (resp. ) is the number of integrators (resp. derivatives) traversed by the causal path.

And ( ) is the characteristic polynomial of the whole bond graph and ( ) is the characteristic polynomial of the reduced (or complementary) bond graph ; it is obtained by suppressing the elementary causal path.

Rule 2. Characteristic polynomial of a bond graph model

The characteristic polynomial, of a linear (or linearized) bond graph model that may be regular or generalized, is given by the formula hereafter:

Where is the order of the bond graph model, is the symbolic Laplace variable. The coefficient and the coefficients ; with , are given by the equation (31):

Where:

: is the order of the causal cycle family; with . And ; with (resp. ) is the number of integrators (resp. derivatives) traversed by the causal cycles (or loops) involved in the causal cycle family.

: is the number of disjoint causal cycles (or loops) involved in the causal cycle family.

̃ ( ) : is the constant term of the gain ( ) ( ) of the causal cycle family of order .

The causal gain of the causal cycle family of order is:

Where ̃ is the constant term of the gain of the causal cycle (or loop) and the product is taken over all the disjoint causal cycles (or loops) involved in the causal cycle family of order .

Particularly, the coefficient is the graph (or Mason's) determinant of the bond graph model , and it is given by the equation (33):

Where is the constant term of the disjoint causal loops gain of order 0 taken one by one, is the product of the constant terms of the gains of disjoint loops of order 0 taken two by two and is the product of the constant terms of the gains of disjoint causal loops of order 0 taken three by three, and so on...

Remark 7.

If there is no causal loop of order 0 in the linear (or linearized) bond graph, then the coefficient is reduced to 1.

Remark 8. If in addition to the causal cycle family of order

, the linear (or linearized) bond graph model contains causal cycles (or loops) of order 0, then the causal cycle family is given by the equation (34).

( ) ( )

Where ( ) is the product of the gains of disjoint causal cycles (or loops) of order zero.

In the end, in the general case where the general junction structure (GJS) of the linear (or linearized) bond graph model involves bond loops; the gain of the causal cycle family must be multiplied by the Mason's determinant of the GJS, see for this [START_REF] Lamb | Causal loops and Mason's rule for bond graphs[END_REF].

The rule 3 is dedicated to recognizing the linearly independent output detectors, directly from the bond graph model in integral causality, by simple causal manipulations. It is fundamental tools are the bicausality and the state-tooutput causal paths concepts.

Rule 3. Linear independency of output detectors

Let be a linear (or linearized) bond graph model, in preferred integral causality and it has storage elements in integral causality; i.e. its order is . It is supposed fitted by a set of output detectors that is defined by the equation (35):

Where stands for an output detector with unspecified causality, (resp.

) stands for an effort (resp. flow) detector, is the output variable measured by the output detector and ( ) stands for equivalent.

Then, the output detectors are linearly independent if, and only if, when they are substituted by bicausal double sources elements ; with , the following three properties hold simultaneously: