Lghani Menhour
email: lghani.menhour@mines-paristech.frbrigitte.dandrea-novel@mines-paristech.fr

Brigitte D'andréa-Novel

Clément Boussard
email: clement.boussard@inria.fr

Michel Fliess
email: michel.fliess@polytechnique.edu

Hugues Mounier

Algebraic nonlinear estimation and flatness-based lateral/longitudinal control for automotive vehicles

A combined longitudinal and lateral vehicle control is presented. It employs flatness-based control and new algebraic estimation techniques for the numerical differentiation of noisy signals. This nonlinear control is designed for automatic path-tracking via vehicle steering angle and driving/braking wheel torque. It combines the control of the lateral and longitudinal movements in order to ensure an accurate tracking of straight or curved trajectories. It can also be used to perform a combined lane-keeping and steering control during critical driving situations such as obstacle avoidance, stop-and-go control, lane-change maneuvers or any other maneuvers. Promising results have been obtained using the noisy experimental data acquired by a laboratory vehicle under high dynamic loads and characterized by high lateral accelerations.

I. INTRODUCTION

The research on intelligent transportation systems, the aim of which is to improve road safety, is getting more and more active. Some assistance systems have been developed in order to help the driver in critical situations, and to improve the vehicle stability and steerability. Among those, we can find for example Anti-Lock Braking System, Electronic Stability program, Adaptive Cruise Control and recently four wheel steering systems.

To develop autonomous vehicles, the design of combined lateral and longitudinal control is required to perform some coupled maneuvers as stop-and-go control with obstacle avoidance, lane-change maneuvers, and the appropriate steering angle to ensure lane keeping. A majority of the strategies proposed in the previous published works deal with longitudinal or lateral control separately, and rarely with a combined implementation, for example to perform a braking action in curve. In [START_REF] Plochl | Driver models in automobile dynamics application[END_REF], a survey of driver models and their application in automotive dynamics is presented, and papers referred in this large inventory deal only with steering vehicle control. In [START_REF] Rajamani | Demonstration of integrated longitudinal and lateral control for the operation of automated vehicles in platoons[END_REF] an integration of longitudinal and lateral control for hugues.mounier@lss.supelec.fr platooning in highways is presented in a decoupled way. The steering control strategy for the autonomous lane change in highway or for urban steering driving are also presented in several works, like [START_REF] Ackermann | Linear and nonlinear controller design for robust automatic steering[END_REF], [START_REF] Cerone | Combined automatic lanekeeping and driver's steering through a 2-DOF control strategy[END_REF], [START_REF] Guldner | Robust automatic steering control for look-down reference systems with front and rear sensors[END_REF], [START_REF] Hatipoglu | Automated Lane Change Controller Design[END_REF], [START_REF] Marino | Input-output decoupling control by measurement feedback in four wheel-steering-vehicles[END_REF].

Longitudinal control is considered in some papers, for example in [START_REF] Martinez | A safe longitudinal control for adaptive cruise control and stop-and-go scenarios[END_REF], [START_REF] Nouvelière | Commandes robustes appliquées au contrôle assisté d'un véhicule à basse vitesse[END_REF], [START_REF] Tsugawa | Vision-based vehicles in Japan: Machine vision systems and driving control systems[END_REF], [START_REF] Villagra | Robust stop-and-go control strategy: an algebraic approach for nonlinear estimation and control[END_REF] with some stop-and-go control strategies to keep a safety distance and avoid collisions.

On the other hand, some works have dealt with the problem of comfort and improvements of a ground vehicle handling through the joint braking and suspension control, see for example [START_REF] Poussot-Vassal | Attitude and handling improvements through gain-scheduled suspensions and brakes control[END_REF]. In [START_REF] Poussot-Vassal | Vehicle dynamic stability improvements through gain-scheduled steering and braking control[END_REF] this approach is extended by the introduction of the steering control in order to improve the vehicle stability, the longitudinal motion is not considered as a state component but as a variable parameter. 1In [START_REF] Nobe | An overview of recent developments in automated lateral and longitudinal vehicle controls[END_REF], [START_REF] Khodayari | A historical review on lateral and longitudinal control of autonomous vehicle motions[END_REF], the authors propose an inventory concerning longitudinal and lateral control, and point out that in the literature the longitudinal control and the steering control are treated separately in most of the cases. However, in [START_REF] Andréa-Novel | Tracking with stability for a vehicle braking in a corner[END_REF], [START_REF] Chou | Global vehicle control using differential braking torques and active suspension forces[END_REF], a nonlinear global chassis control law is proposed to follow desired trajectories in yaw rate and in longitudinal acceleration via braking torques only without active steering control.

The present paper deals with the problem of a nonlinear coupled longitudinal and lateral control based on algebraic nonlinear estimation techniques [START_REF] Fliess | Non-linear estimation is easy[END_REF], [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF] and differential flatness [START_REF] Fliess | Flatness and defect of nonlinear systems: introductory theory and examples[END_REF], [START_REF] Fliess | A Lie-Bäcklund approach to equivalence and flatness of nonlinear systems[END_REF] (see also [START_REF] Fuchshumer | Nonlinear vehicle dynamics control -A flatness based approach[END_REF] for a related approach). To deal with this method, a three degrees-of-freedom two wheels vehicle model is considered. The estimation approach allows to obtain a very good estimation of the derivatives of noisy signals. Indeed the sensors used in a real automotive application are generally low cost and their measurements are very noisy. Hence it is not easy to use these measurements directly to obtain the reference signals for the control applications. For this reason, algebraic estimation is useful.

The remainder of this paper is organized as follows. Section II presents the 3DoF two wheels vehicle model used for the control law design. Section III describes the design steps of a flatness-based robust combined longitudinal and lateral non linear control and algebraic nonlinear estimation. Section IV presents detailed experimental results using the noisy measurements acquired on real roads with high dynamic loads. Some concluding remarks and future work are presented in Section V.

II. VEHICLE MODELS AND COUPLED DYNAMIC

Two nonlinear vehicle models are used: the first one is a 3DoF two wheels nonlinear model (3DoF-NLTWVM) which will be used to design a combined nonlinear control law. The second is a 10DoF four wheels nonlinear model (10DoF-NLFWVM) which will be used as a complete vehicle simulator to test the proposed control law under high loads with curved trajectories. In the following subsections we present the reasons for which the coupled nonlinear vehicle models are used.

A. Coupled dynamics

The aim of designing a combined lateral and longitudinal controller is to enhance and ensure a certain level of control performance, and realize some combined maneuvers. Therefore, to ensure a combined robust control, coupled vehicle dynamics must be taken into account. Clearly, to satisfy these requirements, the identification of some coupling effects is an important step to design a combined control law:

• The first coupling that can be cited is the one related to the kinetics vehicle. The longitudinal motion is affected by the steering angle; in fact when the vehicle wheels are steered, the longitudinal forces are functions the wheel forces and steering angle. The lateral deviation is also a function of the longitudinal speed and yaw angle. The same reasoning is done for the lateral forces and the longitudinal displacement.

B. 3DoF Nonlinear Two Wheels Vehicle Control Model

The 3DoF-NLTWVM shown in Figure 1 used to design the combined control law provides a sufficient approximation of the longitudinal and lateral dynamics of the vehicle in normal driving situations. This model is composed from the longitudinal V x , lateral V y and yaw ψ motions (see table I for notations).

The coupled equations of this model are: The forces and moments of Eq. ( 1) in the vehicle coordinates taking into account the kinetics coupling are:

   ma x = m( Vx -ψV y ) = (F x1 + F x2 ) ma y = m( Vy + ψV x ) = (F y1 + F y2 ) I z ψ = M z1 + M z2 (1)
               F x1 = F x f cos(δ ) -F y f sin(δ ) F x2 = F xr F y1 = F x f sin(δ ) + F y f cos(δ ) F y2 = F yr M z1 = L f (F y f cos(δ ) + F x f sin(δ )) M z2 = -L r F yr (2)
Considering the small angles assumption, the forces and moments, in Eq. ( 2) become as follows: Eq. ( 1) the longitudinal forces are calculated using the dynamical model of the wheels. For the front wheel:

               F x1 = F x f -F y f δ F x2 = F xr F y1 = F x f δ + F y f F y2 = F yr M z1 = L f (F y f + F x f δ ) M z2 = -L r F yr (3) f L r L 2 x xr F F  1 y F  v x 2 y yr F F  xf F 1 x F yf F x V y V  v y V CoG
F x f = (1/R)(-I r ω f + T m -T b f ) (4)
Considering that the vehicle is only propelled by the front wheel, the dynamical model of the rear wheel is as follows:

F xr = -(1/R)(T br + I r ωr ) (5)
The wheels rotations ω f and ω r are supposed to be measured through odometers. They are introduced as external inputs in the equations of the 3DoF-NLTWVM. We assume that the rear braking torque T br is equal to the front one T b f . The wheels rotation accelerations will be estimated using the algebraic nonlinear estimation and the measured rotation speeds which are available on the vehicle's CAN bus. Note that the lateral tire forces, we have, modeled as proportional to the slip angles of each axle. The front and rear forces (F y f = C f α f , F yr = C r α r ) are defined as follows:

   F y f = C f δ - V y + ψL f V x F yr = -C r V y -ψL r V x (6)
where C f and C r are the cornering stiffness for the front and rear tires respectively. The longitudinal speed is considered different from zero (V x > ε). Replacing Eq. ( 3), ( 4), ( 5) and ( 6) in Eq. ( 1) yields:

                         m Vx = m ψV y -I r R ( ωr + ω f ) + 1 R (T m -T b f -T br ) +C f V y +L f ψ V x δ -C f δ 2 m Vy = -m ψV x -C f V y +L f ψ V x -C r V y -ψL r x +(1/R)(T m -T b f )δ + C f -I r R ω f δ I z ψ = -L f C f V y +L f ψ V x + L r C r V y -L r ψ V x +(L f /R)(T m -T b f )δ + L f (T m -I r R ω f )δ (7)
The longitudinal movement is controlled via the traction/braking wheel torque

T ω = T m -T b with T b = T b f + T br ,
and the lateral movement via the steering angle δ . Define the two control variables:

u 1 = T ω u 2 = δ
The model given by Eq. ( 7) may now be written as follows:

ẋ = f (x,t) + g(x)u + g 1 u 1 u 2 + g 2 u 2 2 ( 8
)
where

f (x,t) =      ψV y -I r mR ( ωr + ω f ) -ψV x + 1 m -C f V y +L f ψ V x -C r V y -L r ψ V x 1 I z -L f C f V y +L f ψ V x + L r C r V y -L r ψ V x      , g(x,t) =     1 mR C f m V y +L f ψ V x 0 (C f R -I r ω f )/mR 0 (L f C f R -L f I r ω f )/I z R     , g 1 =    0 1 mR L f I z R    , g 2 =   -C f m 0 0   T , x =   V x V y ψ   T , u = u 1 u 2 T

III. FLATNESS-BASED CONTROL

In order to reduce the complexity of the nonlinear model in Eq. ( 8), nonlinear terms such as u 1 u 2 and u 2 2 are neglected. 2 Despite these simplifications some coupled behaviors are kept as shown by the functions f (x,t) and g(x,t). Eq. ( 8) becomes:

ẋ = f (x,t) + g(x,t)u (9)
A. Flatness property

Consider the system ẋ = f (x, u) where x = (x, • • • , x n ) ∈ R n and u = (u, • • • , u m ) ∈ R m .
It is said to be differentially flat (see [START_REF] Fliess | Flatness and defect of nonlinear systems: introductory theory and examples[END_REF], [START_REF] Fliess | A Lie-Bäcklund approach to equivalence and flatness of nonlinear systems[END_REF] and [START_REF] Lévine | Analysis and Control of Nonlinear Systems -A Flatnessbased Approach[END_REF], [START_REF] Sira-Ramírez | Differentially Flat Systems[END_REF]) if, and only if,

• there exists a vector-valued function h such that

y = h(x, u, u, • • • , u (r) ) (10)
where y = (y,

• • • , y m ) ∈ R m , r ∈ N; • the components of x = (x, • • • , x n ) and u = (u, • • • , u m )
may be expressed as

x = A(y, ẏ, • • • , y (r x ) ), r x ∈ N u = B(y, ẏ, • • • , y (r u ) ), r u ∈ N
Remember that y in Eq. ( 10) is called a flat output.

B. Flatness-based longitudinal and lateral control

1) A proof of flatness: We want to show that y 1 and y 2 given by y 1 = V x y 2 = L f mV y -I z ψ define a flat output. Some algebraic manipulations yield 2 The results of Section IV fully justify those approximations. and

x = V x V y ψ T = A(y 1 , y 2 , ẏ2 ) =       y 1 y 2 L f m -I z L f m L f my 1 ẏ2 +C r (L f +L r )y 2 C r (L f +L r )(I z -L r L f m)+(L f my 1 ) 2 - L f my 1 ẏ2 +C r (L f +L r )y 2 C r (L f +L r )(I z -L r L f m)+(L f my 1 ) 2       (11)
ẏ1 ÿ2 = ∆(y 1 y 2 , ẏ2 ) u 1 u 2 + Φ(y 1 , y 2 , ẏ2 ) ⇒ u 1 u 2 = ∆ -1 (y 1 , y 2 , ẏ2 ) ẏ1 ÿ2 -Φ(y 1 , y 2 , ẏ2 ) (12)
where

∆ 11 (y 1 , y 2 , ẏ2 ) = 1 mR ∆ 12 (y 1 , y 2 , ẏ2 ) = C f m V y + L f ψ y 1 , ∆ 22 (y 1 , y 2 , ẏ2 ) = + -L f my 1 + L r C r (L f +L r ) y 1 (L f C f R-L f I ω ω f ) I z R + (C r (L f +L r )(V y -L r ψ)-L f m ψy 2 1 ) y 2 1 C f (V y +L f ψ) my 1 - C r (L f +L r ) y 1 RC f -I ω ω f mR ∆ 21 (y 1 , y 2 , ẏ2 ) = (C r (L f + L r )(V y -L r ψ) -L f m ψy 2 1 ) mRy 2 1 Φ 1 (y 1 , y 2 , ẏ2 ) = ψV y - I ω mR ( ωr + ω f ) Φ 2 (y 1 , y 2 , ẏ2 ) = -L f my 1 f 3 (x,t) - C r (L f +L r ) y 1 f 2 (x,t) + C r (L f +L r y -L r ψ)-L f m ψy 2 1 y 2 1 f 1 (x,t) + L r C r (L f +L r ) y 1 f 3 (x,t)
The flatness property holds therefore if the matrix

∆(y 1 , y 2 , ẏ2 ) is invertible. It reads det(∆(y 1 , y 2 , ẏ2 )) = ∆ 11 ∆ 22 -∆ 21 ∆ 12 = (Iω ω f -C f R) L 2 f y 2 1 m 2 -Cr(L f +L r )L r L f m+C r I z L I z R 2 y 1 m 2 = 0
This determinant, which depends only on the longitudinal speed y 1 = V x , is indeed nonzero:

• The wheel rotation acceleration is less than RC f /I ω : RC f /I ω is around 10 4 , then

I ω ω f -C f R = 0. • Notice that I z > L f m, then: C r (L r + L f )(I z -L f m) + L 2 f m 2 y 2 1 = 0. Thus u = T ω δ = B(y 1 , ẏ1 , y 2 , ẏ2 , ÿ2 ) = ∆ -1 (y 1 , y 2 , ẏ2 ) ẏ1 ÿ2 -Φ(y 1 , y 2 , ẏ2 ) (13)
with r x = 1 and r u = 2.

2) A tracking feedback control: In order to track the desired output trajectories y re f 1 and y re f 2 , set

ẏ1 ÿ2 = ẏre f 1 + K 1 1 e y 1 + K 2 1 e y 1 dt ÿre f 2 + K 1 2 ėy 2 + K 2 2 e y 2 + K 3 2 e y 2 dt (14)
where, e y 1 = y re f 1y 1 = V re f x -V x and e y 2 = y re f 2y 2 . The choice of the gain parameters K 1 1 , K 2 1 , K 1 2 , K 2 2 and K3 2 is straightforward.

C. Algebraic nonlinear estimation

It should be pointed out that the control law contains derivatives of reference signals which are estimated from measurements such as V re f

x , V re f y , ψre f and the derivatives of the measured front and rear rotation speed wheels ω f and ω r . In order to perform this task and avoid the effects of the noisy measurements we propose to use the numerical differentiation based on an algebraic nonlinear estimation. 3 This estimation is performed using the recent advances in [START_REF] Fliess | Non-linear estimation is easy[END_REF], [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF], which yield efficient real-time filters. The following formulae (see, e.g., [START_REF] García Collado | Analyse fréquentielle des dérivateurs algébriques[END_REF]) may be used:

• Denoising: ŷ(t) = 2! T 2 t t-T (3(t -τ) -T )y(τ)dτ (15)
• Numerical differentiation of a noisy signal:

ẏ(t) = - 3! T 3 t t-T (2T (t -τ) -T )y(τ)dτ (16)
Note that the sliding time window [t -T, t] may be quite short.

Remark 3.1: The above estimation methods are not of asymptotic type and do not require any statistical knowledge of the corrupting noises (see [START_REF] Fliess | Analyse non standard du bruit[END_REF] for details).

IV. COMPARISON WITH EXPERIMENTAL TESTS

Several trials in presence of coupled traction/braking and steering maneuvers have been realized under high dynamic loads using the laboratory vehicle. During these tests a large set of dynamic parameters has been acquired such as longitudinal and lateral speeds, lateral and longitudinal accelerations, yaw rate, wheel rotation speeds, all forces and moments on four wheels, steering angle, etc. For these tests the acquisition device operates at frequency 400 Hz.

These data have been used in our simulation process as reference trajectories, and compared to the obtained simulation results in closed loop with the experimental ones. More precisely, the combined control approach presented here is tested using several trials data performed on real race track as shown by curvature of Figure 4. This section presents a test conducted on a race track which is performed under high dynamic loads as shown in Figures 2, 3, 4 and5 which summarize some dynamic parameters of longitudinal and lateral motions. This test presents a coupled steering, accelerating, a sudden and sharp braking maneuver.

Figures 2 and3 show the dynamic parameters of the two controlled vehicle models (3DoF-NLTWVM and 10DoF-NLFWVM) which are close to those measured such as yaw rate, longitudinal speed, longitudinal and lateral accelerations. Notice also that the steering angle and braking/traction wheel torque provided by the combined control law are similar to the measured ones, this, for any given vehicle model (3DoF or 10DoF) and the performed maneuvers as illustrated in Figures 2 and3 between the positions from 600m to 750m and from 1000m to 1080m.

these two phases as mentioned on Figures 2 and3, the combined steering and braking maneuvers are performed simultaneously in the most dangerous bends as illustrated in Figure 4. The performances of the control law in terms of tracking trajectories errors are depicted in Figure 5. Let us emphasize that these errors are quite small and satisfactory whatever the vehicle model used, the nature of the trajectory and the dynamic loads (lateral acceleration -5 m/s 2 a y 5 m/s 2 ).

During this test, the control law is efficient in simulation which is emphasized by the behavior of the controlled models being close to the measured one. This remark is verified especially when the 10DoF vehicle model is coupled with the nonlinear flat control law, indeed, as shown by Figure 4, the curvature of 3DoF-NLTWVM + nonlinear flat control is slightly different from the reference road curvature. From the results obtained with the two tests, the control law provided more realistic behaviors in terms of the tracking trajectories, when coupling the braking/traction torque and steering angle on a realistic simulator 10DoF-NLFWVM. This shows the effectiveness and the ability of the combined control law to operate under extreme and coupled nonlinear vehicle behaviors (for lateral acceleration -5 m/s 2 a y 5 m/s 2 and longitudinal acceleration -5 m/s 2 a x 3.5 m/s 2 ) in curved trajectories (20 m < R < 35 m). It should be pointed out that the emergency situations, such as collisions and obstacle avoidance, lane-change maneuvers, can easily cause critical dynamical driving situations. For this reason the performance of the control law must be tested under high dynamics loads.

V. CONCLUSIONS AND FUTURE WORK

The comparison with experimental results shows that the proposed control law is able to track a coupled longitudinal and lateral desired dynamics even under high dynamics loads and sudden and sharp maneuvers.

Future works will be devoted to the study of some combined tasks, like obstacle avoidance or lane change maneuvers during stop-and-go control to maintain a constant spacing between the cars using the inter-distance model [START_REF] Martinez | A safe longitudinal control for adaptive cruise control and stop-and-go scenarios[END_REF].

•

  The load transfer is an important coupling of the vehicle dynamics. This coupling is the consequence of the braking/traction (change the weight distribution between the front and rear tires), lateral acceleration (change the weight distribution between the right and left tires), road bank angle, road slope angle, roll and suspension motions. All these coupling are considered in the computation of the tire vertical forces.• Another coupling effect comes from the behavior of the tires; this coupling is expressed by the longitudinal and lateral forces. In fact, to obtain a realistic representation of the vehicle behaviors in extreme driving situations, the coupling of longitudinal and lateral forces is required. To achieve this goal, the coupled formula of Pacejka[START_REF] Pacejka | Tire and Vehicle Dynamics[END_REF] model is used. In this model the coupling of longitudinal slip, lateral slip, vertical forces and camber angle is taken into account.

Fig. 1 .

 1 Fig. 1. Nonlinear two wheeled vehicle model

Fig. 2 .

 2 Fig. 2.Results of longitudinal behaviors: longitudinal acceleration, longitudinal speed and traction/braking torque control

Fig. 3 .Fig. 4 . 4 . 1 :

 3441 Fig. 3. Results of lateral behaviors: lateral acceleration, yaw rate and steering angle control

Fig. 5 .

 5 Fig. 5. Tracking errors on first and second outputs

See also[START_REF] Zheng | Controller design for vehicle stability enhancement[END_REF].

See[START_REF] Abouaïssa | Vers une caractérisation non linéaire d'un réseau autoroutier[END_REF],[START_REF] Villagra | Robust stop-and-go control strategy: an algebraic approach for nonlinear estimation and control[END_REF],[START_REF] Villagra | A diagnosisbased approach for tireroad forces and maximum friction estimation[END_REF] for previous successful applications to intelligent transportation systems.

This work was supported by the French national project INOVE/ANR 2010 BLANC 308. L. Menhour and B.


Pauline ; Krug

Vincent Moriniere

Sandrine Marlin

Gabriel Heinz

Valérie Koubi

Dominique Bonneau

Estelle Colin

Remi Salomon

Corinne Antignac

M.D Laurence Heidet
email: laurence.heidet@nck.aphp.fr

Ap-

Vincent Morinière

Heinz D Gabriel

Rémi Salomon

Mutation screening

Keywords: BOR syndrome, EYA1, SIX1, SIX5, developmental defects BOR syndrome, EYA1, SIX1, SIX5, develomental defect

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Mutations in SIX1, a DNA binding protein that associates with EYA1, have been reported less frequently. One group has recently described 4 missense mutations in SIX5 in 5 unrelated patients with BOR.

Here, we report a screening of these three genes in a cohort of 140 patients from 124 families with BOR. We identified 36 EYA1 mutations in 42 unrelated patients, 2 mutations and one change of unknown significance in SIX1 in 3 unrelated patients, but no mutation in SIX5. We did not find correlation between genotype and phenotype, and observed a high phenotypic variability between and within BOR families. We show the difficulty in establishing a molecular diagnosis strategy in BOR syndrome: the screening focusing on patients with typical BOR would detect a mutation rate of 76%, but would also miss mutations in 9% of patients with atypical BOR. We detected a deletion removing three EYA1 exons in a patient who was previously reported to carry the SIX5 Thr552Met mutation. This led us to reconsider the role of SIX5 in the development of BOR.

Formatted: Don't adjust space between Latin and Asian text, Don't adjust space between Asian text and numbers Deleted: the screening of only patients with typical BOR would lead to a 76% rate of mutation detection but would also miss mutations in 9% of atypical BOR Deleted: .

INTRODUCTION

Branchio-oto-renal (BOR) syndrome is an autosomal-dominant developmental disorder which is characterized by hearing loss, branchial arch defects and various renal anomalies. The prevalence of BOR syndrome is estimated to be 1 case per 40 000 [START_REF] Chen | Phenotypic manifestations of branchio-oto-renal syndrome[END_REF][START_REF] Fraser | Genetic aspects of the BOR syndrome-branchial fistulas, ear pits, hearing loss, and renal anomalies[END_REF][START_REF] Fraser | Frequency of the branchio-oto-renal (BOR) syndrome in children with profound hearing loss[END_REF][START_REF] Melnick | Autosomal dominant branchiootorenal dysplasia[END_REF][START_REF] Melnick | Branchio-oto-renal dysplasia and branchio-oto dysplasia: two distinct autosomal dominant disorders[END_REF]. The syndrome is clinically heterogeneous and has a high penetrance with variable expressivity [START_REF] Fraser | Genetic aspects of the BOR syndrome-branchial fistulas, ear pits, hearing loss, and renal anomalies[END_REF][START_REF] Fraser | Frequency of the branchio-oto-renal (BOR) syndrome in children with profound hearing loss[END_REF], Chen et al., 2004). BOR syndrome is also genetically heterogeneous.

Over 80 mutations in EYA1 (MIM ID 601653), the human homolog of the Drosophila eyes absent gene, encoding a transcriptional regulator, have been identified. These include large and small heterozygous deletions, frameshift, stop, splice-site and missense heterozygous mutations (Abdelhak et al., 1997b[START_REF] Ni | Refined localization of the branchiootorenal syndrome gene by linkage and haplotype analysis[END_REF], Vincent et al., 1997). The rate of detection of EYA1 mutations varies from 7% to 40% of patients tested according to the clinical criteria required for molecular testing (Abdelhak et al., 1997a, Abdelhak et al., 1997b, Chang et al., 2004, Orten et al., 2008). Mutations in SIX1 (MIM ID 601205) (mainly missense mutations and small deletions), the human homolog of sine oculis encoding a DNA binding protein that associates with EYA1, have also been associated with BOR syndrome [START_REF] Kochhar | SIX1 mutation screening in 247 branchio-oto-renal syndrome families: a recurrent missense mutation associated with BOR[END_REF][START_REF] Ruf | A gene locus for branchio-otic syndrome maps to chromosome 14q21.3-q24.3[END_REF][START_REF] Ruf | Formatted: German (Germany) 17 SIX1 mutations cause branchio-oto-renal syndrome by disruption of EYA1-SIX1-DNA complexes[END_REF][START_REF] Sanggaard | Branchio-oto-renal syndrome: detection of EYA1 and SIX1 mutations in five out of six Danish families by combining linkage, MLPA and sequencing analyses[END_REF], though much less frequently than EYA1 mutations. More recently, missense mutations in another SIX family member, SIX5 (MIM ID 600963), have been reported by one group in patients with BOR syndrome [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF]. SIX5 homologous is known to interact with eya-1 in C. elegans. In vitro functional analyses of the BOR-associated SIX5 variants showed that some of these variants modified EYA1-SIX5 binding and the ability of the EYA1-SIX5 complex to transactivate a [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF]. However, the association of SIX5 mutations with BOR syndrome has not been confirmed by other groups.

In the present study, we screen for EYA1, SIX1 and SIX5 mutations a large cohort of patients with BOR syndrome. We describe the clinical features associated with the mutations and the rate of mutations identified, according to the clinical phenotypes. We also show that one family previously reported as carrying a SIX5 missense mutation harbors a heterozygous deletion of three EYA1 exons, which therefore questions the role of the reported SIX5 change.

PATIENTS AND METHODS

Patients

A total of 140 patients from 124 families with a diagnosis of BOR syndrome were included in the study. Subjects were classified according to the criteria defined by Chang [START_REF] Chang | Branchio-oto-renal syndrome: the mutation spectrum in EYA1 and its phenotypic consequences[END_REF] on the basis of clinical history, audiometry and renal ultrasonography. Patients were considered as typical BOR when they had at least three major criteria (branchial anomalies, deafness, preauricular pits or renal anomalies), or two major and two minor criteria (internal, middle and/or external ear anomalies, preauricular tags, facial asymmetry or palatine anomalies) or one major criterion and an affected first-degree relative meeting the above criteria for typical BOR. Other were considered as atypical BOR and were tested only when they demonstrated at least two features of the syndrome.

Patient 1062 was previously reported as carrying a heterozygous SIX5 c.1655C>T (p.Thr552Met) mutation (patient A500 in Hoskins et al.). His DNA had been tested for EYA1 mutations by direct sequencing, but not for abnormal copy number [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF].

This patient was having assisted reproduction, and thus was making inquiries regarding the possibility of preimplantation genetic testing.

Mutation analysis

Genomic DNA was isolated from peripheral blood using standard methods. The 16 exons of EYA1 were screened for mutations by direct sequencing. When no mutation was found, quantitative multiplex PCR amplification of short fluorescence fragments [START_REF] Charbonnier | Detection of exon deletions and duplications of the mismatch repair genes in hereditary nonpolyposis colorectal cancer families using multiplex polymerase chain reaction of short fluorescent fragments[END_REF] was performed for EYA1 exons 1, 5, 10, 15 and 16. When a deletion was found to remove some but not all of the exons tested, the DNA sample was analysed by multiplex ligation dependent probes amplification (Salsa MLPA kit P153 EYA1, MRC-Holland, Amsterdam, Netherlands). When neither mutation nor deletion was found in EYA1, the 2 exons of SIX1 and the 3 exons of SIX5 were screened by direct sequencing. For previously unreported missense mutations, 92 control individuals were tested by direct sequencing. The cDNAs NM_172060.2 for EYA1, NM_005982.3 for SIX1 and NM_175875.4 for SIX5 were used for numbering, with nt +1 corresponding to the A of the ATG translation initiation codon. One DNA sample (patient 608), was shown to have large 8q13.3 deletion by FISH analysis and was used as a positive control for deletion screening. Missense mutations were evaluated using the softwares PolyPhen (http://genetics.bwh.harvard.edu/pph/) and ConSurf (http://consurftest.tau.ac.il/).

Statistical tests

Testing for difference in proportions was performed using either the X 2 or Fisher's exact test.

All tests were two sided. P-values <0.05 were considered significant.

RESULTS

Mutations

The mutations identified and the associated phenotypes are shown in table 1. In the entire cohort (140 patients from 124 families) we identified 36 EYA1 mutations in 42 families (55 patients), two SIX1 mutations and one SIX1 variant of unknown significance in 3 families (4 patients), but no SIX5 mutation. We identified two SIX5 variants which were not considered to be responsible for the phenotype: one was found in controls, and the other, previously reported as a disease causing mutation in two families [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF], was associated with a partial EYA1 deletion in one of these two families included in our cohort. Therefore, mutations were identified in 36% (45/124) of the tested families.

EYA1 gene analysis

Thirty-six EYA1 heterozygous mutations, spread over the entire length of the gene (figure 1), were identified in 42 probands (table 1). Thirty-three were small mutations (8 missense including a mutation of the stop codon, 14 frameshift, 6 stop, and 5 splice-site mutations), and 24 of these were novel. All missense mutations but one [c.319G>A (p.Gly107Ser)] were considered as possibly or probably damaging by the PolyPhen software (http://genetics.bwh.harvard.edu/pph/). That mutation c.319G>A, which appeared de novo in the patient, was scored as benign by Polyphen but was in the last base of exon 4 and thus was expected to modify the splicing of intron 4 (GeneSplicer score changes from 6.97 to 2.49) .

Amino acid conservation scores according to ConSurf (varying from 1 to 9) for previously unpublished missense mutations are shown in table 1.Two previously reported mutations, c.982C>T and c.1220G>A, were respectively found in three and two unrelated patients. Three different deletions were identified in 5 unrelated patients by quantitative multiplex PCR amplification of short fluorescence fragments. In 3 of these probands (patients 608, 821, 991) all tested exons were missing and the deletion was considered to remove the entire gene. In the other two cases, the deletion was partial and its precise extent was determined by MLPA analysis: one removed exons 8 to 16 (patient 1216) and one removed exons 3, 4 and 5 (patients 1062, 1063 and 1064). Patient 1062 was one of the 2 probands previously reported to carry the SIX5 p.Thr552Met mutation (patient A500 in [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF]. Parent status was tested for 26 probands with an identified EYA1 mutation : 8 mutations out of 26 were de novo and 18 were inherited.

SIX1 gene analysis

Three different, potentially pathogenic variations in SIX1 were identified in 3 families (table 1). The SIX1 mutation (p.Tyr129Cys) has already been described [START_REF] Ruf | Formatted: German (Germany) 17 SIX1 mutations cause branchio-oto-renal syndrome by disruption of EYA1-SIX1-DNA complexes[END_REF], and affects a conserved tyrosine in the homeodomain. It is predicted to be probably damaging by Polyphen (score 2.945), and inhibits the transcription activation in vitro [START_REF] Patrick | Biochemical and functional characterization of six SIX1 Branchio-oto-renal syndrome mutations[END_REF]. This mutation was also present in the affected father of the proband. The mutation c.560+3A>T, probably leading to aberrant splicing (GeneSplicer score changes from 9.74 to 2.89), has never been previously described. Finally, the c.746C>T change in exon 2 affects a strongly conserved amino acid (p.Pro249Leu) and was considered as possibly damaging by PolyPhen (score 1.806). Although this change was not present in 92 healthy control chromosomes it is located in a region of unknown function and replaces a non polar side chain amino acid with another amino acid of the same family. The segregation of these two last changes could not be tested because DNA samples from family members were not available.

SIX5 gene analysis

We did not identify any novel SIX5 mutations in our entire cohort. We confirmed the finding of the SIX5 p.Thr552Met heterozygous variant in patient 1062, and also found the same SIX5 ). However, we identified an EYA1 partial deletion in the three affected members of this family (see above). We found another SIX5 variant, c.156_161dup (p.Gly55Ala56dup), in a patient from Guadeloupe. That variant, which introduces two amino acids in the N-terminus of the protein, was also found in 3/86 controls from the West Indies and thus was considered as non pathogenic.

Phenotypes (tables 1 and 2)

According to previously described criteria [START_REF] Chang | Branchio-oto-renal syndrome: the mutation spectrum in EYA1 and its phenotypic consequences[END_REF], our population included 67 patients with typical BOR and 55 patients with atypical BOR. In 18 cases the information we had was insufficient to classify patients as typical or atypical. Sixty eight probands had a family history of BOR (38 typical cases of BOR and 30 cases of renal anomaly, branchial arches defects and/or deafness, which did not fulfill the typical BOR diagnosis criteria). The phenotypic features of patients with an identified mutation are detailed in table 1.

The frequency of each symptom in our entire population as well as in patients with mutation and in patients without mutation is summarized in table 2. The various renal phenotypes observed in patients with an identified mutation are described figure 2

Genotype-phenotype correlation

Because the type and severity of the symptoms were very variable, we searched whether there was a correlation between the phenotype and either the mutated gene (EYA1 or SIX1) or the type of mutation (missense mutation, truncating mutation because of stop, frameshift, or splice-site mutation), or deletion. Of the 67 patients with typical BOR syndrome, 50 (75%), patients with atypical BOR syndrome, 5 (9%) had a mutation (in EYA1 in 4 cases and in SIX1 in 1 case). The proportion of patients affected with typical or atypical BOR syndrome was not significantly different (p=0.38) between patients carrying a missense mutation (9/11 typical), a truncating mutation (34/37 typical) or a deletion (8/8 typical). The same mutation (p.Arg328X) was identified in a patient with atypical BOR and in 2 patients with typical BOR. We did not observed any particular phenotypic features associated with SIX1 mutations. Rare features were observed both in patients carrying EYA1 or SIX1 mutations (assuming that the SIX1 p.Pro249Leu is responsible for the phenotype). Among the 55 patients with EYA1 mutations, the proportion of deletion (n=9), missense (n=10), and truncating (n=36) mutations was not significantly different between patients with (n=40) or without (n=15) deafness (p=0,46), with (n=30) or without (n=25) kidney involvement (p=0.35), with (n=32) or without (n=23) pits (p=1), or with (n=34) or without (n=11) branchial defect (p=0.66) (Table 3).

Regarding the renal disease, cases with prenatal renal failure associated with oligoamnios (n=4) or with severe renal failure leading to renal transplant (n=5) were associated with SIX1 mutations in one case, and with EYA1 mutations in 7 cases. These proportions were not significantly different from that observed in all patients. However, in the 7 cases with EYA1 mutations, none of these mutations were a missense mutations (3 were frameshift, 1 splicesite, 1 stop and 2 were entire gene deletions).

DISCUSSION

To our knowledge, this study is the first one to analyze all of the genes currently known to be implicated in BOR syndrome in a large cohort of patients. Overall, we detected a mutation in 45/124 (36%) probands. Forty two probands were carrying an EYA1 mutation thought to be pathogenic. However, the impact of the c.867+5G>A change on mRNA splicing has not yet . Three probands were carrying mutations in SIX1. Of these, one (patient 714) was carrying the c.806C>T, p.Pro249Leu change, which has never been described before, and is predicted to be possibly damaging by Polyphen. This mutation was not found in 92 control chromosomes, but does not affect a protein domain with a known function. All other SIX1 mutations reported so far affect either the SIX or the homeodomain encoding nucleotides [START_REF] Kochhar | SIX1 mutation screening in 247 branchio-oto-renal syndrome families: a recurrent missense mutation associated with BOR[END_REF][START_REF] Ruf | Formatted: German (Germany) 17 SIX1 mutations cause branchio-oto-renal syndrome by disruption of EYA1-SIX1-DNA complexes[END_REF]. We were unable to test other family members of patients 714 and 715, so it is difficult to conclude whether this change is or not a disease-causing mutation.

In our entire cohort, we did not identify any pathogenic mutation in the SIX5 gene. In one patient previously reported to carry a SIX5 missense variant (case 1062), we found a partial (exons 3-5) EYA1 deletion. The three affected patients in this family were carrying both the EYA1 deletion and the SIX5 variant. We believe that the EYA1 deletion is responsible for the phenotype in this family, though we cannot rule out the hypothesis that the SIX5 variant may modify the EYA1-associated phenotype. However, whereas the three patients had deafness, the renal disease was more severe in the two siblings (undergoing renal transplantation at 22 and 23 years of age) than in their father (who had not reached end-stage renal failure at 58 years). This was despite the fact that all three carried the SIX5 variant. The finding of an EYA1 mutation in that family made us reconsider the role of SIX5 in the development of BOR syndrome. Among the 5 index cases reported by Hoskins et al. as carrying a SIX5 mutation [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF], all carried a missense variant, including two cases with the c.1655C>T p.Thr552Met variant (patient 1062 and another patient). The segregation of the variants with the phenotype had not been studied. These variants modestly (20 to 48%) although significantly decreased the ability of SIX5/EYA1 to activate gene transcription in vitro [START_REF] Hoskins | Transcription factor SIX5 is mutated in patients with branchio-oto-renal syndrome[END_REF]. No other SIX5 mutation (whether missense or other type of mutation/rearrangment) has been reported since this initial report. In addition, whereas mice Six5-/-mice is limited to abnormalities in the eye (cataract), and does not affect the ear or the kidney. Taken together, these findings strongly suggest that SIX5 mutations are not directly responsible for BOR syndrome.

Among the patients that were classified with either typical or atypical BOR (122 patients), we identified a mutation in 75% of cases with typical BOR syndrome and in 9% of cases with atypical BOR syndrome. These results are different from those reported recently in a smaller cohort in which no mutations were detected in any subject with atypical BOR (Rickard et al., 2008). This highlights the difficulty in reconciling the need for performing molecular testing in a consequential and cost effective manner, and the fact that a screening limited to typical BOR syndrome will miss few mutations and prevent accurate genetic counselling in these few families. The rate of mutation that we report here is not different from that recently reported in another large cohort [START_REF] Orten | Branchio-oto-renal syndrome (BOR): novel mutations in the EYA1 gene, and a review of the mutational genetics of BOR[END_REF]. However, the rate of EYA1 deletion in the present study is lower than that (18%) reported by [START_REF] Chang | Branchio-oto-renal syndrome: the mutation spectrum in EYA1 and its phenotypic consequences[END_REF]. Although we used the same semiquantitative fluorescence multiplex PCR approach for tracking EYA1 deletions, we only tested exons 1, 5, 10, 15, and 16 in a first attempt. We may have thus missed small or complex deletions involving other exons.

In our series as in others [START_REF] Chang | Branchio-oto-renal syndrome: the mutation spectrum in EYA1 and its phenotypic consequences[END_REF][START_REF] Ruf | Formatted: German (Germany) 17 SIX1 mutations cause branchio-oto-renal syndrome by disruption of EYA1-SIX1-DNA complexes[END_REF], Saanggard et al., 2007, Okada et al., Orten et al., 2008), the type and severity of the phenotype does not seem to correlate with the type of mutation and is very variable, even within a given family. Only the severity of renal failure may correlate to some extent with the type of EYA1 mutation, as none of the 7 patients with the most severe renal insufficiency were carrying a missense mutation. However the small number of patients does not allow any conclusion to be made, and it would be interesting to analyze the severity of the renal failure in a larger number of cases carrying an EYA1 mutation. The high frequency of renal anomalies in our series may be due to the fact F o r P e e r R e v i e w 13 that our laboratory is focused on renal diseases. The important phenotypic variability observed in our cohort as in others may be due to genetic and/or environmental factors (in particular maternal factors during embryonic and fetal development) that modify the phenotype. This makes genetic counselling particularly difficult for these families.

We report some interesting clinical features associated with EYA1 mutations. Although already reported by others [START_REF] Chen | Phenotypic manifestations of branchio-oto-renal syndrome[END_REF][START_REF] Shimasaki | EYA1 mutation in a newborn female presenting with cardiofacial syndrome[END_REF] the association with hypothyroidism or with persistent ductus arteriosus may be fortuitous, as may be the association with a conotruncal cardiopathy. Two patients suffered from cataract, which could be associated with a defect of early expression of EYA1 in the anterior ocular segment, and has already been described [START_REF] Azuma | Mutations of a human homologue of the Drosophila eyes absent gene (EYA1) detected in patients with congenital cataracts and ocular anterior segment anomalies[END_REF]. A facial nerve impairment was present in six patients, which may be explained by anomalies of inner ear, affecting the nerve trajectory.

In conclusion, our results confirm that EYA1 is by far the most frequently mutated gene in BOR syndrome today and seriously question the role of SIX5 variants in the pathophysiology of BOR syndrome. They confirm the lack of genotype-phenotype correlation and illustrate the difficulty in establishing an algorithm for molecular diagnosis in BOR syndrome. The screening of patients with typical BOR only would greatly increase the rate of identified mutations but would also lead to missed mutations in a few families. Our data suggest that testing patients with atypical BOR still results in the identification of a few mutations, but screening of the SIX5 gene can be given up. Human Mutation

  -renal (BOR) syndrome is an autosomal dominant disorder characterized by branchial, ear and renal anomalies. Over 80 mutations in EYA1 have been reported in BOR.

  , medical records, genealogy and written informed consent from patient and/or parents were sent from Paediatric, Paediatric Nephrology, Nephrology, or Genetics departments between August 2004 and December 2009.

Formatted

Formatted

  twin brother (patient 1063) and in his father (patient 1064), who were both affected. This variant was predicted as possibly damaging by the Polyphen program (score 1.711

  . The renal function of these patients varied greatly, from a normal glomerular filtration rate to end-stage renal failure. In four cases prenatal renal failure diagnosed by oligoamnios during the second trimester of pregnancy was observed in fetuses displaying kidney hypoplasia(patients 175, 700, 991 and 1126), and this led to termination of pregnancy in three of these cases. Five patients received a renal transplant: two reached end stage renal failure during childhood (at 14 and 16 years), one at 26 years, and two others at unknown age (but one received a transplant at 22). In some families (see cases 700 and 1126 as examples), although a fetus or a Deleted: (figure2) Deleted: c.161_162insCGGGGC (p.54_55 ins GA) renal disease, the affected parent had normal renal function and normal kidneys on ultrasound examination. After renal disease, the most frequent symptom was deafness (conductive and/or sensorineural), followed by pits, branchial defect and then by ear tags. Some patients presented with rare phenotypic features. In patients with an identified mutation, two (patients 311 and 1291) presented with cataract (unilateral in one case and bilateral in the other). Five patients had facial nerve impairment: unilateral facial paralysis, crocodile tears syndrome and facial asymmetry. Palate anomalies were observed in 5 patients: short palate responsible for hypernasal speech, arched palate, palatine cleft, bifid uvula and posterior agenesia of uvula. We observe two heart defects: a persistent ductus arteriosus in patient 708 and a conotroncal cardiopathy (pulmonary atresia with interventricular communication) in patient 700. Finally, patient 229 suffered from hypothyroidism of unknown etiology.In patients without any identified mutation, rare clinical features were also observed : palate or laryngeal anomalies (n=4), facial asymmetry (n=2), bifid uterus (n=2), interventricular communication (n=2), clinodactyly (n=1), aortic arch hypoplasia (n=1), cutaneous pigmentation anomalies (n=2), Malherbe's calcifying epithelioma (n=1), dorsal medullar atrophia responsible for pyramidal syndrome (n=1), pulmonary hypertension of unknown etiology (n=1), cerebral ventricular dilation (n=1), thyroglossal duct cyst (n=1), jejuno-ileal atresia (n=1), and bone anomalies (radial aplasia, mandibular hypoplasia, dental agenesis, postaxial polydactyly) in four patients.

  in EYA1, and one had a change in SIX1 of unknown significance. Of the 55

with

  Six1 or Eya1 gene inactivation develop ear and kidney abnormalities, the phenotype in

  Tours, France), Pr Cochat (Hôpital de Lyon, France), Dr Demeer (Hôpital d'Amiens, France), Pr Dollfus (Hôpital de Strasbourg, France), Dr Drouin (Hôpital de Rouen, France), Dr Faivre (Hôpital de Dijon, France), Dr Duban-Bedu (Hôpital de Lille, France), Pr Grunfeld (AP-HP, Hôpital Necker, Paris, France), Pr Lacombe (Hôpital de Bordeaux, France), Pr Manach (AP-HP, Hôpital Necker, Paris, France), Pr Cormier-Daire (AP-HP, Hôpital Necker, Paris, France), Dr Mathieu (Hôpital d'Amiens, France), Dr Francannet (Hôpital de Clermont-Ferrand, France), Pr Manouvrier (Hôpital de Lille, France), Dr Baumann (AP-HP, Hôpital Robert Debré, Paris, France), Dr Isidor (Hôpital de Nantes, France), Pr Deschênes (AP-HP, Hôpital Robert Debré, Paris, France), Dr Dupin-Deguine (Hôpital de Toulouse, France).

Figure 1 :

 1 Figure legendFigure 1 : Schematic representation of EYA1 gene. Boxes represent EYA1 exons. Mutations: frameshift , nonsense , missense , splice site , deletion EyaHR: eyes absent homologous region.

Figure 2 :

 2 Figure 2: Renal phenotypes in patients with an identified mutation.

Formatted

  type of EYA1 mutation according to the type of symptom.

Table 1 : mutations and phenotypes in patients with BOR syndrome

 1 . : T : typical BOR, A : atypical BOR, b : branchial anomalies, d : deafness, p : pits, t : tags, ie : inner ear, me : middle ear, ee : extern ear anomalies, r : renal anomalies, f : familial, U : unknown, TOP : termination of pregnancy, I : inherited, N : de novo. BTNRH and MORL : Boys Town National Research Hospital and Molecular Otolaryngology Research Laboratory (UI, USA), see http://www.healthcare.uiowa.edu/labs/pendredandbor/. * : not found in 92 healthy individuals. ** :according to ConSurf (http://consurftest.tau.ac.il/). The reference sequences for numbering are NM_172060.2 for EYA1, NM_005982.3 for SIX1 and NM_175875.4 for SIX5 .

	1

Formatted: Font: Not Bold Formatted: Font: (Default) Arial, 8 pt, Not Bold, English (U.K.) Formatted: Font: (Default) Arial, 8 pt, Not Bold, English (U.K.) Formatted: Font: (Default) Arial, 8 pt, Not Bold, English (U.K.) Page 22 of 29 John Wiley & Sons, Inc.

Table 2 :

 2 clinical features of the whole population (n=140). Columns with and without mutation represent the ratio of patients with each symptoms reported on the total number of patients with (n=59) or without (n= 81) mutation, respectively

John Wiley & Sons, Inc.Human Mutation

Acknowledgments

We thank the patients and their family for their participation. We are grateful to the following physicians for contribution of material and clinical data from patients : Pr Melki (AP-HP Hôpital de Bicetre, France), Dr Toutain (Hôpital de Tours, France), Pr Bensman (AP-HP, Hôpital Trousseau, Paris, France), Dr Cartault (Hôpital de St Denis, France), Dr Blesson Formatted: English (U.S.)


Farouk Benmeddour

Treyss Fabien

Laurent Ède

Laguerre

Numerical modeling of guided wave interaction with non-axisymmetric cracks in elastic cylinders

Keywords: Guided waves, Cylinders, Cracks, 3D Hybrid FE-SAFE

come

Numerical modeling of guided wave interaction with non-axisymmetric cracks in elastic cylinders

Farouk Benmeddour, Fabien Treyssede,

INTRODUCTION

Non-destructive testing (NDT) techniques using guided ultrasonic waves constitute viable inspection means for tubes, pipes and plates owing to their potential to carry out energy over long distances and their sensitivity to internal damages. Nevertheless, the interaction of such waves with damages is complex since elastic guided waves are multi-modal and dispersive. This highlights the essential role of modeling techniques in providing a better understanding of and a deeper insight into the nature of phenomena arising from scattering problems.

Various analytical approaches have been proposed in the literature. For instance, one can cite quasi-static approximations based on the S -parameter formalism for plates [START_REF] Auld | Acoustic fields and waves in solids[END_REF] and tubes [START_REF] Ditri | Utilization of guided elastic waves for the characterization of circumferential cracks in hollow cylinders[END_REF], and modal decomposition methods recently developed for plates [START_REF] Clézio | The interaction of the S0 Lamb mode with vertical cracks in an aluminum plate[END_REF][START_REF] Castaings | Modal decomposition method for modeling the interaction of Lamb waves with cracks[END_REF][START_REF] Demma | Scattering of the fundamental shear horizontal mode from steps and notches in plates[END_REF]Shkerdin andGlorieux, 2004, 2005). Such methods are fast from a computational point of view but generally limited to simple geometries with horizontal or vertical discontinuities.

Fully numerical approaches based on transient finite element (FE) analyses have been conducted for plates [START_REF] Alleyne | The interaction of lamb waves with defects[END_REF]Lowe and Dili-gent, 2002;Benmeddour et al., 2008b,a) and pipes [START_REF] Lowe | Defect detection in pipes using guided waves[END_REF][START_REF] Demma | The reflection of guided waves from notches in pipes: a guide for interpreting corrosion measurements[END_REF][START_REF] Ma | Scattering of the fundamental torsional mode by an axisymmetric layer inside a pipe[END_REF]. The main advantage of these approaches is that complex-shaped waveguides or damages can be handled with standard FE codes. However, transient FE models are time consuming, which often limits their practical use to two-dimensional problems and short propagation distances. Besides, a modal post-processing step in the frequency domain must be performed, which requires an a priori knowledge of guided modes. Recently, [START_REF] Gunawan | Mode-exciting method for lamb wavescattering analysis[END_REF] have developed a mode-exciting method to analyse the interaction of Lamb waves with defects. The originality of their work is in the analysis of a finite plate model using a standard FE or boundary element (BE) method software with a time-harmonic regime.

However, in addition to the a priori knowledge of guided modes, their analysis requires a great number of computations at high frequencies. [START_REF] Moreau | An orthogonality relation-based technique for post-processing finite element predictions of waves scattering in solid waveguides[END_REF] and [START_REF] Predoi | Influence of material viscoelasticity on the scattering of guided waves by defects[END_REF] have used an existing FE software with absorbing regions at the edge of plates. Harmonic FE solutions have been post-processed with the help of a general orthogonality relation.

Hybrid approaches, combining a FE method and the so-called normal mode expansion technique have been developed for studying wave scattering by damages in plates [START_REF] Karunasena | Plane-straine-wave scattering by cracks in laminated composite plates[END_REF][START_REF] Karunasena | Guided waves in a jointed composite plate[END_REF][START_REF] Al-Nassar | Scattering of lamb waves by a normal rectangular strip weldment[END_REF][START_REF] Mal | A semi-numerical method for elastic wave scattering calculations[END_REF]. They consist in performing a modal decomposition both at the inlet and outlet of a waveguide FE model. It can then be limited to a small region surrounding the damage. The modal decomposition is performed by virtue of analytical solutions for the guided modes of the considered structure. The scattered solution is computed for each frequency and directly yields the coefficients of scattered modes without any post-processing. The same procedure based on a BE method has also been applied by [START_REF] Zhao | Boundary element mode for defect characteriza[END_REF].

In the literature, analytical solutions for guided modes, limited to simple cross-section geometries, have been replaced in several occasions by the semianalytical finite element (SAFE) technique (see e.g. [START_REF] Lagasse | Higher-order finite-element analysis of topographic guides supporting elastic surface waves[END_REF]; [START_REF] Gavrić | Computation of propagative waves in free rail using a finite element technique[END_REF]; [START_REF] Hladky-Hennion | Finite element analysis of the propagation of acoustic waves in waveguides[END_REF]; [START_REF] Hladky-Hennion | Finite element analysis of the propagation of acoustic waves along waveguides immersed in water[END_REF]; [START_REF] Damljanović | Propagating and evanescent elastic waves in cylindrical waveguides of arbitrary cross section[END_REF]; [START_REF] Hayashi | Guided wave dispersion curves for a bar with an arbitrary cross-section, a rod and rail example[END_REF][START_REF] Hayashi | Wave structure analysis of guided waves in a bar with arbitrary cross-section[END_REF]). This has given rise to the so-called hybrid FE-SAFE and BE-SAFE methods capable of handling arbitrary cross-sections with complex-shaped damages. The former has been successfully applied to two-dimensional geometries such as plates [START_REF] Karunasena | Hybrid analysis of lamb wave reflection by a crack at the fixed edge of a composite plate[END_REF][START_REF] Datta | Ultrasonic waves and material and defect characterization in composite plates[END_REF] and cylinders [START_REF] Rattanawangchroen | Axisymmetric guided waves in jointed laminated cylinders[END_REF][START_REF] Rattanawangchroen | Axisymmetric guided waves in jointed laminated cylinders[END_REF], and so has been the latter for plates (Galán andAbascal, 2003, 2005). It is noteworthy that the periodic FE method can be seen as an alternative for the SAFE method when dealing with arbitrary cross-sections, as proposed by Zhou et al. (2009) for the study of thin pipes with non-axisymmetric damages. In a mathematical framework, a Dirichlet-to-Neumann approach has also recently been proposed for isotropic plates [START_REF] Baronian | Transparent boundary conditions for the harmonic diffraction problem in an elastic waveguide[END_REF] based on a mixed formulation and biorthogonality relation.

The aim of this work is to study the interaction of guided waves with non-axisymmetric cracks in infinite circular cylinders. Such a study requires a three dimensional (3D) formulation. To this end, a general 3D hybrid FE-SAFE procedure is developed in Sec. 2. The 3D variational formulation of elastodynamics is first recalled. A modal expansion is applied at the crosssection boundaries of the waveguide, which enables the separation between ingoing and outgoing modes. Eigenmodes are computed using the SAFE technique. A formula is proposed for the direct calculation of eigenforces from SAFE matrices. For hybrid FE-SAFE methods, this may greatly simplify the post-process of consistent load eigenvectors, usually performed from displacement derivatives. In the same way, expressions for power reflection and transmission coefficients are also presented. This paves the way for the implementation of hybrid FE-SAFE methods inside standard FE softwares.

Section 3 deals with the interaction of the fundamental Pochhammer-Chree modes with a free-end and non-axisymmetric vertical cracks inside circular isotropic cylinders. For the free-end cylinder, the obtained results are in good agreement with those already published in the literature. For non-axisymmetric vertical cracks, the results of their interaction with the fundamental guided modes are obtained and discussed. It should be noted that some of the analytical approaches mentioned earlier could also be used for this kind of discontinuities. However, to the authors' knowledge, vertical cracks inside cylinders have surprisingly not received a great attention in the literature.

Finally, the case of oblique discontinuities is briefly studied in Sec. 4 where different inclination angles are analysed. Results are compared with those obtained for vertical discontinuities. These waveguides are situated on the left (L) and right (R) of the interior region (I) of the volume V . Each cross-section has the same uniform geometry along the propagation direction. Incident waves are launched from the left waveguide in the positive direction (+). Then, reflected waves from the damage travel in the negative direction (-). Transmitted waves in the right waveguide propagate in the positive direction. The interior region is connected to the left and right waveguides through boundaries denoted S L and S R , respectively. The whole volume V is bounded by S L and S R and by a traction free boundary, denoted S I .

Hybrid FE-SAFE method

General description of the problem

Finite element formulation for the interior region

For the volume V , a standard FE method with a time harmonic regime is used. The application of virtual works to the equation of motion governing elastodynamics gives, after dropping the time dependence e -jωt , the following variational formulation:

V δǫ T σdV -ω 2 V ρδu T udV - S δu T tdS = 0, (1)
where S = S L ∪ S I ∪ S R is the boundary surrounding the volume V , δ(•) denotes virtual fields, ω is the angular frequency and ρ is the material

density. u = [u x u y u z ] T is
t i = σ ij n j (i, j = x, y, z), where n = [n x n y n z ]
T is the outward normal to the surface S.

The FE discretisation of the volume V leads to an approximation of the displacement field u based on the use of interpolating functions. For a given element e, the approximation of displacement is u = N e U e , where N e is the matrix of nodal interpolating polynomials and U e is the vector of nodal displacements of the element. Assembling elementary matrices into global ones, Eq. (1) finally yields:

δU T K -ω 2 M U -δU T f = 0, ( 2
)
where U is the vector of degrees of freedom (dofs) of displacement, K and M are the global stiffness and mass matrices, respectively, and f represents a vector of external time harmonic forces.

In this study, interior forces are not considered. The above discretised variational formulation can be partitioned as follows:

         δU L δU I δU R          T      D LL D LI D LR D IL D II D IR D RL D RI D RR               U L U I U R          -          δU L δU I δU R          T          f L 0 f R          = 0, (3)
where U L , U I and U R respectively denote the left, interior and right displacement dofs.

D ij = K ij -ω 2 M ij (i, j = L, I, R
) are the associated partitions.

Mode expansions on waveguide cross-sections

On S L and S R , displacements and forces are expanded into sums of modes.

On S L , the sum is decomposed into ingoing (+) and outgoing (-) modes. On S R , only the outgoing (+) modes are taken into account which enforces the coefficients of ingoing (-) modes to zero giving rise to a transparent boundary condition. For propagating modes, the traveling direction is determined from the sign of the energy velocity, whereas for non-propagating modes it is determined from the sign of the imaginary part of the wave number.

Displacements and forces on the left and right sections are then given by:

U L = N + L n=1 α + Ln U + Ln + N - L n=1 α - Ln U - Ln , f L = N + L n=1 α + Ln f + Ln + N - L n=1 α - Ln f - Ln , (4)
and:

U R = N + R n=1 α + Rn U + Rn , f R = N + R n=1 α + Rn f + Rn , (5)
where n denotes the mode number, N ± L and N + R are the number of modes considered in the expansions. α ± Ln (n = 1, ..., N ± L ) and α + Rn (n = 1, ..., N + R ) are the left and right displacement modal amplitudes, respectively. From a theoretical point of view, N ± L and N + R should tend to infinity. However, in practice, series are truncated by retaining the less attenuated modes. All the propagating and only some of the non-propagating modes are taken into account. The choice of the number of these non-propagating modes is based on their attenuations, which are related to the imaginary part of the wavenumber.

Eqs. ( 4) and ( 5) can be transformed into products of matrices. For the left section, coefficients of positive-going modes are imposed so that δα + Ln = 0. Displacements, forces and virtual displacements can then be expressed as:

U L = B + L α + L + B - L α - L , U R = B + R α + R , f L = T + L α + L + T - L α - L , f R = T + R α + R , δU L = δα - L T B - L T , δU R = δα + R T B + R T , (6)
with the following notations:

α ± L = α ± L1 α ± L2 ... α ± LN ± L T , B ± L = U ± L1 U ± L2 ... U ± LN ± L , T ± L = f ± L1 f ± L2 ... f ± LN ± L , (7)
α ± L denotes column vectors of the amplitudes of positive and negative going modes of the left section. B ± L defines displacement modal bases. T ± L denote modal bases of forces f ± Ln . B + R , T + R and α + R are defined in the same way by replacing the subscripts L with R.

With the help of Eqs. ( 6), δU can be rewritten as:

δU L δU I δU R = δα - L δU I δα + R      B - L T 0 0 0 I 0 0 0 B + R T      , ( 8
)
where I is the identity matrix.

For any trial field δα - L , δU I and δα + R , the substitution of Eqs. ( 6) and ( 8) into Eq. ( 3) yields after rearrangements the following linear global system:

     B - L T 0 0 0 I 0 0 0 B + R T                D LL D LI D LR D IL D II D IR D RL D RI D RR           B - L 0 0 0 I 0 0 0 B + R      -      T - L 0 0 0 0 0 0 0 T + R                    α - L U I α + R          = -      B - L T 0 0 0 I 0 0 0 B + R T                D LL D LI D LR D IL D II D IR D RL D RI D RR           B + L 0 0      -      T + L 0 0           α + L . (9)
The unknown vectors α - L , U I and α + R are found by solving the above linear system at each frequency. If only the left semi-infinite waveguide is present, the system given by Eq. ( 9) is modified by removing the third row and column from matrices and the third element from vectors. The right cross-section surface, in this case, is traction free.

Modal bases of displacements and forces are needed. In this paper, they are obtained by using a SAFE technique, which is explained in the next subsection.

Semi-analytical finite element formulation

The SAFE method is used in order to determine the elastic guided modes for the left and right cross-sections of the interior region. This technique reduces the analysis of three dimensional waveguides to two dimensions by using a spatial Fourier transform along the propagation direction (u(x, y, z) = u(x, y)e +jkz , where k is the wavenumber and z is the propagation direction).

Therefore, only the cross-sections have to be meshed, which allows a fast and accurate computation of eigenmodes (wavenumbers and modeshapes) for any arbitrary cross-section.

The SAFE variational formulation can be written as:

δU K 1 -ω 2 M S -jk K 2 -K T 2 + k 2 K 3 U = 0 , (10)
obtained from the following elementary matrices:

K e 1 = S e B T 1 CB 1 dS , K e 3 = S e B T 2 CB 2 dS , K e 2 = S e B T 2 CB 1 dS , M e S = S e ρN eT N e dS ,
where B 1 = L x N e ,x + L y N e ,y and B 2 = L z N e (N e ,x and N e ,y are the differentiation of the interpolation functions with respect to x and y, respectively).

In the above equation, when the left (right) section is considered, U and S would designate U L (U R ) and S L (S R ). L x , L y and L z are defined as follows:

L x =               1 0 0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 0               , L y =               0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 1               , L z =               0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 0 1 0               . ( 11
)
For more details, the reader can refer to [START_REF] Hayashi | Guided wave dispersion curves for a bar with an arbitrary cross-section, a rod and rail example[END_REF] for instance.

By fixing the angular frequency ω and finding the wavenumber k, Eq. ( 10) leads to a quadratic eigenvalue problem. This eigenproblem can be linearized (see e.g. [START_REF] Tisseur | The quadratic eigenvalue problem[END_REF]). The resolution of this system for the left (right) cross-section gives wavenumbers k Ln (k Rn ) and modeshapes U Ln (U Rn ). The separation into ingoing and outgoing modes then yields the displacement modal bases B ± L,R used in the previous subsection. The next step is to determine the modal bases of forces T ± L,R . The idea is to derive forces directly from SAFE matrices. For clarity, let us first consider the right cross-section. From the third left integral of Eq. ( 1) and the application of normal mode expansions, forces on the right section can be rewritten as:

δU T R f R = S R δu T N + R n=1 α + Rn t + Rn dS. ( 12
)
For one element e of surface S e R , δu T = δU e T N e T . Besides, following the same approach as in [START_REF] Treyssède | Elastic waves in helical waveguides[END_REF] for computing the energy velocity, it can be

verified that t + Rn = L T z σ + Rn = L T z C B 1 + jk + Rn B 2 U +e
Rn . The substitution of this expression into Eq. ( 12) yields:

δU T R f R = δU T R N + R n=1 α + Rn K 2 + jk + Rn K 3 U + Rn . (13)
The above result defines explicitly the forces f + Rn , involved inside the modal basis T + R of the previous subsection, as:

f + Rn = K 2 + jk + Rn K 3 U + Rn . ( 14
)
where K 2 and K 3 are the SAFE matrices associated with the right section.

Following the same approach for the left section, the same kind of expression is obtained for the left forces f ± Ln :

f ± Ln = K 2 + jk ± Ln K 3 U ± Ln , (15)
where K 2 and K 3 are the SAFE matrices now associated with the left section.

Once the displacement modal bases are obtained, the calculation of modal bases T ± L and T + R is then direct from Eqs. ( 14) and ( 15) and consistent with the used FE approximation. For hybrid FE-SAFE methods, this may simplify greatly the tedious post-process of consistent loads based on displacement derivatives associated with each mode.

Power coefficients

The definition of the power flow through a surface S (S = S L or S R ) is given by:

P = S P • n dS, (16)
where P i = -σ ij uj is the Poynting vector and the bar denotes time averaging.

This yields:

P = -j ω 4 S (u * T t -u T t * )dS. ( 17
)
The power of one single mode n, denoted P nn , can then be directly obtained from:

P nn = |α n | 2 ω 2 ℑ U * T n f n , (18)
ℑ designates the imaginary part.

For conciseness of notations, the subscripts L and R are omitted as well as superscripts ±. However, one should be careful with such a definition when checking energy conservation in computations involving flexural modes.

Auld's orthogonality relationship [START_REF] Auld | Acoustic fields and waves in solids[END_REF] states that:

(k m -k * n )P mn = 0, (19)
with:

P mn = -j ω 4 S (u * T m t n -u T n t * m )dS. (20)
Eq. ( 19) shows that mode orthogonality P mn = 0 holds for k m = k * n . Double roots k m = k * n typically occur for flexural modes in cylinders. Thus, the power of a superposition of two flexural modes m and n might not be equal to the sum of their individual powers [START_REF] Gunawan | Mode-exciting method for lamb wavescattering analysis[END_REF]). In such a case, the power flow P ′ mn of such modes ought to be calculated as follows:

P ′ mn = ω 2 ℑ       α m α n    * T U m U n * T f m f n    α m α n       . (21)
Hence, P ′ mn represents the total power associated with the couple of modes (m,n).

Finally, the power reflection and transmission coefficients are computed by dividing the reflected and transmitted powers by the incident power, re- spectively:

R mn = P - Ln P + Lm , T mn = P + Rn P + Lm , ( 22
)
where m is the incident mode, n is the reflected or transmitted mode, and P + Lm , P - Ln and P + Rn are computed with Eq. ( 18) (or Eq. ( 21) for flexural modes).

Numerical results for vertical discontinuities

This section gives a detailed study of the interaction of fundamental guided modes (L(0,1), T(0,1) and F(1,1)) with vertical discontinuities.

Generalities

In the following, a steel circular cylinder with radius r is considered in simulations. The values of the Poisson coefficient, the density and the Young modulus are: ν = 0.25, ρ = 7800 kg/m 3 and E = 2 × 10 11 P a, respectively. The dimensionless frequency is Ω = ω r c L where c L is the longitudinal velocity. The crack geometry considered in this section is depicted in Fig. 2. It can be referred to as a transverse notch [START_REF] Mohr | On inspection of thin-walled tubes for transverse and longitudinal flaws by guided ultrasonic waves[END_REF] characterised by an arbitrary maximum depth d and an area located in the vertical transverse plane at a given z. The opening widths (w) of cracks is small with respect to the wavelength considered so that its effect can be neglected. The only crack parameter is hence the ratio d r . The selected incident mode is launched at the left section of the cylinder with a unit amplitude. Structures are meshed with quadratic hexahedrons with 20 nodes per element. The maximum element length is h = 0.2r yielding min(λ) h > 10, where λ = 2πr Ω . The rough criterion of h < λ 4 for quadratic elements is hence quite satisfied.

It should be emphasized that the modal expansions used at the left and right sections also account for non-propagating modes. This allows to set both sections relatively close to the crack thus reducing the FE system. In SAFE computations, the searching interval for the eigenmodes is ℑ(kr) ∈ [-5; +5], ∀Ω. These computed modes are all retained in the modal expansions of the hybrid method. A good accuracy is expected because the amplitudes of non-propagating modes scattered from cracks, but not retained in the expansions, would then be divided by at least e 5 ≃ 148 at cross-section boundaries, which are located at z = ∓r (see Fig. 3). This has been further confirmed by numerical tests (not shown here) performed by changing the number of modes in the truncation.

Figure 4 shows the dispersion curves of dimensionless energy velocity versus dimensionless frequency obtained with the SAFE method on one of the two identical sections. Mode labels are in accordance with [START_REF] Meitzler | Mode coupling occurring in the propagation of elastic pulses in wires[END_REF].

Free-end cylinder: numerical validation

This section deals with the scattering of the L(0,1) and F(1,1) modes from the free-end cylinder. The frequency ranges of both modes are chosen to enable the comparison with results found in the literature [START_REF] Gregory | Axisymmetric waves in a semi-infinite elastic rod[END_REF][START_REF] Rattanawangcharoen | Reflection of waves at the free edge of a laminated circular cylinder[END_REF][START_REF] Taweel | Wave reflection from the free end of a cylinder with an arbitrary cross-section[END_REF]. Modal expansions are applied to the left section while the right section and the surrounding surface are traction free. and L(0,3) (-•) when the L(0,1) mode encounters a free-end cylinder. The bullets (•) are results found in [START_REF] Gregory | Axisymmetric waves in a semi-infinite elastic rod[END_REF].

Ω R F(1,1) F(1,2) F(1,1) F(1,1) F(1,3) F(1,2) F(1,4) F(1,1)
Figure 6: Power reflection coefficients of flexural modes when the F(1,1) mode encounters a free-end cylinder. The bullets (•) are results found in [START_REF] Taweel | Wave reflection from the free end of a cylinder with an arbitrary cross-section[END_REF].

encounters a free-end cylinder. The range of the dimensionless frequency Ω is restricted to vary between 1.8 and 2.4. As expected, no mode conversions to flexural modes are observed due to the axisymmetric geometry. These results are compared with those obtained by [START_REF] Gregory | Axisymmetric waves in a semi-infinite elastic rod[END_REF].

A very good agreement is found.

Figure 6 shows the power reflection coefficients when the F(1,1) mode interacts with a free-end cylinder. The range of the dimensionless frequency Ω, in this case, varies between 0.6 and 3.3. As before, no mode conversions to compressional modes are observed. These results compare successfully with those obtained by [START_REF] Taweel | Wave reflection from the free end of a cylinder with an arbitrary cross-section[END_REF].

The obtained results validate the modal expansion procedure at crosssection boundaries and its implementation.

Infinite cylinder: non-axisymmetric crack analysis

In this section, the scattering by non-axisymmetric damages is analysed for vertical cracks having three different depths and for the incident L(0,1), F(1,1) and T(0,1) modes. Modal expansions are now applied on both the left and right sections. The surfaces of cracks and the surrounding surface are traction free.

First, a numerical test has been performed on the cylinder without cracks (d = 0). Though not presented for conciseness, results have shown that incident modes on one side are fully transmitted to the other side with negligible reflection.

It is emphasised that singularities at crack corners are not taken into account in the FE models. However, as already found by [START_REF] Bai | Scattering of guided waves by circumferential cracks in steel pipes[END_REF] and [START_REF] Alleyne | The reflection of guided waves from circumferential notches in pipes[END_REF], it is expected that inaccuracy remains only very close to the singularities without affecting the global behaviour in the far field and hereby, the power reflection and transmission coefficients. In this paper, it has been also checked that refining the mesh near the crack corner does not affect the computation of reflection and transmission coefficients.

Compressional mode L(0,1)

Figures 7a and 7b illustrate the power reflection (R) and transmission (T ) coefficients of the compressional L(0,1) mode when it encounters a crack with the following depths: d = 0.5r, 1.0r and 1.5r. As can be seen from Fig. 7a, R is significant except for the case of the crack with d = 0.5r. Fig. 7b shows that T is more sensitive to all depths including d = 0.5r. Note that the sharp changes of both coefficients coincide with the cut-off frequencies of different compressional and flexural modes (see Fig. 4). Moreover, as opposed to the free-end case, the axisymmetric L(0,1) mode is converted to flexural modes due to the non-axisymmetric nature of vertical cracks. These mode conversions are discussed in the following.

As shown in Fig. 4, the L(0,1) mode is less dispersive below the first cutoff frequency of the F(1,2) mode (Ω c = 1.06). R and T reveal a monotonic variation with the crack depth over the range defined by the first cut-off frequency (Figs. 7a andb). According to this and due to the fact that the single L(0,1) mode can be easily launched in this frequency range, this fundamental compressional mode is expected to be a good candidate for experiments.

Figures 8a and 8b depict the power reflection coefficients of the converted L(0,3) and F(1,1) modes from the incident L(0,1) mode. The L(0,3) mode presents a significant mode conversion above its cut-off frequency except for d = 0.5r. In the very low frequency range (Ω < 0.2), it can be observed that the sensitivity of the converted F(1,1) mode is of the same order of magnitude as the reflected L(0,1) mode (Fig. 8b). R and T of the converted mode L(0,2) and T of L(0,3) are found to be insignificant in value and are not shown for conciseness. R and T of the converted flexural modes lower than 0.12 for the whole frequency range are not shown either. As in the free-end cylinder case, note that no mode conversion to the torsional T(0,1) mode has been found. This might be explained by the specificity of the T(0,1) mode, which is non-dispersive and whose displacement field remains in the cross-section plane.

For these computations, it has been verified that the sum of the power reflection and transmission coefficients at each frequency equals unity with a maximum error of 0.2%.

Torsional mode T(0,1)

Figures 9a and9b show the power reflection and transmission coefficients of the torsional mode T(0,1) when it encounters the vertical cracks.

In the low-frequency range, the T(0,1) mode exhibits a monotonic increase in the amplitude of the power reflection coefficient (and a decrease in the amplitude of the power transmission coefficient) with increasing crack depth.

However, this behaviour is limited to a shorter frequency range in the case of the T(0,1) mode, up to approximately Ω = 0.2 against Ω = 1 for the L(0,1) Mode. From Ω=0.2 to 1.0, the reflected T(0,1) mode appears to be less attractive than the L(0,1) mode because of the non-monotonic variation of its amplitude with crack depth. Under the cut-off frequency of the F(1,2) mode, we note that the sensitivity of the T(0,1) mode is very low for a crack depth d = 0.5r, as observed earlier for the L(0,1) mode case.

In the frequency range 1 < Ω < 3, the reflected T(0,1) mode turns out to be more sensitive to a crack of depth 0.5r than the reflection of the L(0,1) mode (compare Figs. 7a and9a). As for the incident L(0,1) mode, T is globally more sensitive than R. Although strong local changes may occur at cut-off frequencies, a monotonic increase of R with crack depth is found at certain frequencies in [2.0, 3.0].

According to the non dispersive nature of the T(0,1) mode and since this is the only torsional mode to propagate in this frequency range, this confers to the fundamental T(0,1) mode potential interest for long-range non-destructive testing of low vertical crack depths. encounters the cracks. As can be observed, R remains small except below the cut-off frequency of the F(1,2) mode. T , being of the same order of magnitude as R, is omitted for brevity and so are the converted flexural modes having insignificant values. No mode conversion to the compressional modes L(0,1), L(0,2) or L(0,3) has been found in this frequency range. For these computations, the maximum error on the power balance is 0.07%.

Flexural mode F(1,1):

The flexural mode F(1,1) corresponds to double eigensolutions having the same wavenumber but modeshapes of different orientations. Analytically, one mode is in cos(θ), the other is in sin(θ) (θ being the azimuthal coordinate in the cylindrical system). Consequently, the scattering of both F(1,1) modes will be different when non-axisymmetric cracks are considered. This can be overcome by launching "rotating" F(1,1) modes, written in terms of e ±iθ (instead of cos(θ) or sin(θ)): at a given axial position z, the displacement field of such modes is indeed time-rotating in the cross-section so that the whole section is "scanned". Hence, the incident mode F(1,1) in e ±iθ will be scattered identically by non-axisymmetric vertical cracks. In the following results, the SAFE method gives modes in cos(θ) and sin(θ) but the incident mode F(1,1) is chosen as a linear combination cos(θ) + isin(θ) = e iθ , yielding a rotating mode.

Figures 11a and11b show R and T of the flexural mode F(1,1) when it encounters the cracks. Their abrupt changes coincide with cut-off frequencies as for the modes mentioned earlier. R is not very sensitive to cracks except for the case d = 1.5r at some frequencies. In addition, similar to the incident L(0,1) and T(0,1) modes, T is more sensitive than R. For these computations, the maximum error on the power balance is less than 0.13%.

Application to oblique free-ends and cracks

In order to highlight the potential of the hybrid 3D FE-SAFE method, this section is devoted to the study of the fundamental L(0,1) mode interaction with oblique discontinuities. The sensitivity of this mode to different inclination angles is illustrated. cylinder radius are the same as those in Sec. 3.1.

Figure 14 shows the power reflection coefficients (R) when the L(0,1) mode encounters free-ends with an inclination angle φ of 0 • , 10 • and 20 • .

As observed in this figure, R is sensitive to the free-end inclination. For oblique free-ends, R starts to decrease before the cut-off frequency of L(0,2) and L(0,3) modes, which does not occur for the vertical free-end cylinder (φ = 0 • ). This is explained by mode conversions to flexural modes due to the non-axisymmetry of oblique free-ends. In this figure, it is observed that R decreases as φ increases in the frequency range Ω < 1.9.

The maximum error on the power balance computations is less than 0.38%. vertical crack (φ = 0 • ) and two oblique cracks (φ = 10 • and 20 • ) are studied.

It can be observed that the maximum difference is located around frequencies Ω = 1.0 and 2.7. In the low frequency range (Ω < 0.8), L(0,1) has a low sensitivity to crack inclination.

Figure 16 represents the power reflection coefficients (R) of the converted F(1,1) mode from the incident L(0,1) mode. A significant mode conversion to the F(1,1) mode is observed for the three inclination angles. The reflection coefficients are very close to each other for φ = 10 • and 20 • . However, it can be noticed that a stronger mode conversion occurs for φ = 10 • and 20 • around the frequency Ω = 1.0.

A maximum error of 0.14% on the power balance is found for these computations. From a NDT point of view, it could be concluded that a good sensitivity of the L(0,1) mode to the considered inclination angles is observed when it interacts with free-ends (full cracks). However, this sensitivity is poor for oblique cracks whose depth is lower than the cylinder radius except at particular frequencies.

Conclusions

In this paper, a 3D hybrid numerical method has been developed for the investigation of wave scattering in elastic waveguides. This method combines the classical FE method and the so-called SAFE technique. It has the main advantage of being able to handle complex-shape inhomogeneities and waveguides of arbitrary cross-sections. Moreover, it gives a formula that renders straightforward the calculation of eigenforces and modal power flows.

The validation of the hybrid 3D method has been established for the case of a free-end cylinder by comparing results with those of the literature. Then, the analysis of the scattering of the fundamental longitudinal, torsional and flexural modes in circular cylinders with non-axisymmetric vertical cracks of different depths has been carried out. Cylinders with oblique free-ends and oblique cracks with different inclination angles have also been considered. The interactions of the L(0,1) mode with these discontinuities have been briefly studied. As pointed out throughout the analysis, the proposed method is likely to provide valuable information to further improve experimental investigation on cracks detection and characterisation. It has to be emphasised that the method is quite versatile and that complex-shaped damages can be readily handled.

Figure 1

 1 Figure1depicts a damage located inside an arbitrary volume V . This volume is connected to two semi-infinite waveguides of arbitrary cross-sections.

Figure 1 :

 1 Figure 1: Description of the problem in the case of a damaged volume located between two semi-infinite waveguides (2D representation).

Figure 2 :

 2 Figure 2: Description of a vertical non-axisymmetric crack: side view (a) and cross-section view (b).

Figure 3 :

 3 Figure 3: Meshing of a cylinder with a free-end of d = 2.0r (a) and vertical cracks of depths of 0.5r (b), 1.0r (c) and 1.5r (d).

Figure

  Figure3adepicts the FE mesh used for a free-end cylinder equivalent to a vertical crack of depth d = 2.0r which leads to 4128 dofs. Figs.3b-drepresent the meshing of an infinite cylinder with a depth-variant vertical crack.

FigureFigure 4 :Figure 5 :

 45 Figure 5 depicts the power reflection coefficients when the L(0,1) mode

Figure 7 :

 7 Figure 7: Power reflection (a) and transmission (b) coefficients of the compressional mode L(0,1) when the L(0,1) mode encounters a crack with depths d = 0.5r (--), d = 1.0r (--) and d = 1.5r (-•).

Figure 8 :

 8 Figure 8: Power reflection coefficient of the L(0,3) (a) and F(1,1) (b) modes when the L(0,1) mode encounters a crack with depths d = 0.5r (--), d = 1.0r (--) and d = 1.5r (-•).

Figure 10 Figure 9 :

 109 Figure 10 depicts R of the converted flexural mode F(1,1) when T(0,1)

Figure 10 :

 10 Figure 10: Power reflection coefficient of the flexural mode F(1,1) when the T(0,1) mode encounters a crack with depths d = 0.5r (--), d = 1.0r (--) and d = 1.5r (-•).

Figure 11 :

 11 Figure 11: Power reflection (a) and transmission (b) coefficients of the flexural mode F(1,1) when the F(1,1) mode encounters a crack with depths d = 0.5r (--), d = 1.0r (--) and d = 1.5r (-•).

Figure 12 :

 12 Figure 12: Power reflection coefficient of the flexural mode F(1,2) when the F(1,1) mode encounters a crack with depths d = 0.5r (--), d = 1.0r (--) and d = 1.5r (-•).

R

  of the converted flexural mode F(1,2) is depicted in Fig. 12. It is worth noting the strong reflection of the incident F(1,1) mode into the F(1,2) flexural one near Ω=1.25. Due to the low dispersion of both flexural modes around this frequency (the F(1,2) mode is located near a maximum of energy velocity) and the slower velocity of the F(1,2) mode, flexural polarized waves in this frequency region could provide valuable information to experimental crack detection. Though results are not shown here, significant mode conversions to the reflected F(2,2) mode and to the transmitted F(2,1) are also found. The reflected mode F(1,3) is only sensitive to the case of d = 1.5r. The other converted flexural modes are not sensitive to these crack depths (less than 0.06).

FigureFigure 13 :

 13 Figure 13a sketches an oblique free-end cylinder of angle φ with respect to the vertical plane. The inclined crack is shown in Fig. 13b. The range of the dimensionless frequency Ω is chosen to vary between 0 and 3 as in the previous subsections (see e.g. Sec.3.3.1). The material characteristics and

FiguresFigure 14 :

 14 Figures 15a and b depict the power reflection (R) and transmission (T )coefficients when the L(0,1) mode encounters cracks of depth d = 1.0r. One

Figure 15 :

 15 Figure 15: Power reflection (a) and transmission (b) coefficients of the compressional mode L(0,1) when it encounters a d = 1.0r crack depth of angle 0 • (--), 10 • (--) and 20 • (-•).

Figure 16 :

 16 Figure 16: Power reflection coefficients of the flexural mode F(1,1) when the L(0,1) mode encounters a d = 1.0r crack depth of angle 0 • (--), 10 • (--) and 20 • (-•).

  tion potential in a wave guide. International Journal of Solids and Structures 40, 2645 -2658. Zhou, W., Ichchou, M., Mencik, J., 2009. Analysis of wave propagation in cylindrical pipes with local inhomogeneities. Journal of Sound and Vibration 319, 335-354. Zhuang, W., Shah, A. H., Datta, S. K., 1997. Axisymmetric guided wave scattering by cracks in welded steel pipes. ASME Journal of Pressure Vessel Technology 119 (4), 401-406.

  the vector of time harmonic displacements, where x, y and z are cartesian coordinates and the superscript T denotes transpose, ǫ = [ǫ xx ǫ yy ǫ zz 2ǫ xy 2ǫ xz 2ǫ yz ] T is the strain vector and σ = [σ xx σ yy σ zz σ xy σ xz σ yz ] T is the stress vector. The stress-strain relationship is given by σ = Cǫ, C is the elastic coefficient matrix. t = [t x t y t z ] T represents the external traction defined from:

Acknowledgements

This research work is part of the project ACTENA (Auscultation des Câbles Tendus Non Accessibles) funded by the French ANR (Agence National de Recherche) and EDF ( Électricité De France) and piloted by LCPC (Laboratoire Central des Ponts et Chaussées).