Q Dollon

J Antoni

A Tahan

M Gagnon

C Monette

Operational Modal Analysis of Hydroelectric Turbines Using an Order Based Likelihood Approach

Keywords: Hydropower, Francis turbines, Structure Behavior, Ambient Modal Analysis, Fast-Bayesian Method

The purpose of this paper is to estimate hydroelectric turbine runner modal characteristics from experimental measurements coming from asynchronous regimes. This is achieved by investigating resonances generated by the interaction of a structural mode with harmonics of the rotating speed. Resonances are extracted using order tracking and processed with a Fast Bayesian algorithm in an ambient manner to estimate modal parameters and related uncertainties. Since data collection and effective processing of hydropower turbine field measurements is in its early stages, this paper lays some foundations in the treatment of transient regime measurements. A novelty in this approach lies in the use of a probabilistic identification tool in Order Based Modal Analysis (OBMA). A numerical experiment and a study from a hydroelectric Francis turbine field measurements are introduced to illustrate the method.

Scaled PSD of channel noise 1. Introduction A full grasp of turbine runner dynamic behavior is required in assessing the consequences of fatigue on life duration and crack propagation. Especially, an accurate determination of the actual behavior during transients is mandatory to determine the load levels [1, 2]. High amplitude stress cycles can lead to early damages and loss of reliability [3, 4]. Most of resonances observed in transient regimes emanate from the interaction between synchronous excitations coming from the harmonics of the 2 rotating speed, and a matching vibration mode of the turbine. Such phenomenon can generate periods of intense vibration, possibly leading to extensive damages or even failure of the structure. The prediction of structural modes allows defining optimal operating ranges in which resonances cannot develop, and predicting the dynamic stress amplitudes the structure withstands during operations. In some cases, operation during resonance may be allowed if excitation is sufficiently low and damping sufficiently high, and a proper design is even more important in such situations to ensure safety.

In a foreseeable future, experimental measurements on turbine prototypes are expected to become an integral part of this process [START_REF] Gagnon | Influence of load spectrum assumptions on the expected reliability of hydroelectric turbines: A case study[END_REF][START_REF] Liu | A review on fatigue damage mechanism in hydro turbines[END_REF]. For both designers and users, experimental tools can be used for validating and calibrating numerical models hence increasing structural reliability assessment accuracy. Such approaches are particularly welcomed in the field of hydropower, as operational conditions can strongly alter numerical predictions [START_REF] Valentin | Influence of the boundary conditions on the natural frequencies of a Francis turbine[END_REF][START_REF] Trivedi | Fluid-structure interactions in Francis turbines: A perspective review[END_REF][START_REF] Liu | Cavitation Effects on the Structural Resonance of Hydraulic Turbines: Failure Analysis in a Real Francis Turbine Runner[END_REF], and reduced model resonance studies cannot be transposed to turbine prototypes [START_REF] Valentin | Transposition of the mechanical behavior from model to prototype of Francis turbines[END_REF]. The most practical way to achieve modal analysis for operating devices is Operational Modal Analysis (OMA), where the structural behavior is investigated without prior knowledge of the excitation, which is in effect the ambient excitation occurring during operation. Experimental setup and cost are relatively low compared to other methods, so that OMA has been widely developed since the early 1990's. Usually, only minor assumptions are needed to characterize the response without knowing the excitation, but such analysis result in high uncertainties to be evaluated.

Analysis of turbine operational modal characteristics is still restricted to a few studies [START_REF] Gagnon | Operational Modal Analysis in presence of harmonic excitations: case study on hydroelectric turbine components[END_REF][START_REF] Valentin | Feasibility of Detecting Natural Frequencies of Hydraulic Turbines While in Operation[END_REF][START_REF] Dollon | Dynamic Characterization Of Hydroelectric Turbine In Transient Using OBMA And Phase-Shift Analysis[END_REF], but is gaining recognition with the demonstration of its benefits. Especially, it was shown that transient events are particularly useful for detecting modal signatures [START_REF] Valentin | Feasibility of Detecting Natural Frequencies of Hydraulic Turbines While in Operation[END_REF][START_REF] Dollon | Dynamic Characterization Of Hydroelectric Turbine In Transient Using OBMA And Phase-Shift Analysis[END_REF]. The purpose of this paper is to present an OMA procedure tailored for hydroelectric turbines, capable of extracting modal information and quantifying related uncertainties. The study focuses on asynchronous transient regimes, in which the harmonics usually lead to resonances. System response is processed with a Bayesian-based inference scheme which can return the modal statistical distributions with a low bias.

The paper is structured as follows. The first section introduces the main OMA tools available to estimate modal parameters, and presents the concept of order based modal analysis. the ML Order Based Modal Analysis (ML-OBMA) is presented in section 3, inspired by the works of S.-K. Au [START_REF] Au | Operational Modal Analysis, Modeling, Bayesian Inference, Uncertainty laws[END_REF]. In section 4, several numerical cases are studied to test the algorithm. The model is finally implemented on a hydroelectric case, where operational resonances during a coast-down are characterized.

2. Modal identification under random excitation OMA techniques are divided into four classes: they can process data in the time domain or in the frequency domain, and can be parametric or non-parametric. Non-parametric algorithms perform identification by splitting the signal into its modal contributions: see, e.g., Random Decrement Technique (RDT) [START_REF] Badaoui | The Random Decrement Technique applied to discrete-time structural vibration signals: closed-form solutions for the blind identification of modal parameters[END_REF], Frequency Domain Decomposition (FDD) [START_REF] Brincker | Modal identification of output-only systems using frequency domain decomposition[END_REF], Blind Source Separation (BSS) [START_REF] Sadhu | A review of output-only structural mode identification literature employing blind source separation methods[END_REF]. Parametric approaches, for their part, attempt to fit Multi-Degree-of-Freedom (MDoF) model parameters to experimental data. Time domain models are mostly based on signal autoregression or on the output correlation matrix: auto-regressive (ARMA) [START_REF] Poulimenos | Parametric time-domain methods for non-stationary random vibration modelling and analysis : A critical survey and comparison[END_REF], subspace methods (SSI) [START_REF] Ljung | System Identification[END_REF], eigenrealization algorithm (ERA) or Ibrahim algorithm [START_REF] Brincker | Introduction to Operational Modal Analysis[END_REF]. On the other hand, frequency domain models are mainly based on density matrices: poly-reference Least-Square Complex Exponential (p-LSCE) or polymax [START_REF] Peeters | The PolyMAX frequency-domain method: a new standard for modal parameter estimation?[END_REF], frequency polyreference [START_REF] Brincker | Introduction to Operational Modal Analysis[END_REF], maximum likelihood [START_REF] Yuen | Bayesian methods for Structural Dynamics and Civil Engineering[END_REF].

All these methods are deterministic and mostly rely on curve-fitting. Although ambient modal analysis offers a wide range of identification tools, assessment of uncertainties is still limited to specific procedures, and a general framework is yet to come. This is partly due to the many uncertainty sources and their propagation: data acquisition, identification model, computation process, etc. However, significant progress has recently been made in terms of taking stochastic uncertainties into account [START_REF] Pintelon | Uncertainty calculation in (operational) modal analysis[END_REF]. The determination of modal parameter uncertainties is a first step in estimating the uncertain vibrating response of structures [START_REF] Jian | Non-stationary random response analysis of structures with uncertain parameters[END_REF]. Most time domain modal parameter variance models are derived from the curve-fitting methods, using Monte Carlo algorithms or first-order developments. Essentially, Taylor series are expanded to the first order derivative, for models such as SSI [START_REF] Reynders | Uncertainty quantification in operational modal analysis with stochastic subspace identification: Validation and applications[END_REF][START_REF] Reynders | Uncertainty bounds of modal parameters obtained from stochastic subspace identification[END_REF] or autoregressive [START_REF] Vu | Uncertainties of Modal Parameters by Operational Modal Analysis[END_REF]. In the frequency domain, statistical approaches mostly rely on Maximum Likelihood Estimator (MLE). Such approach was combined with the polymax algorithm to reduce noise and obtain distribution intervals, producing the Polymax Plus algorithm [28, 29]. Some accuracy loss due to very noisy data were addressed more recently [START_REF] Amador | A new maximum likelihood estimator formulated in poleresidue modal model[END_REF]. These algorithms are efficient modal identifiers for long time-histories. MLE was finally extended into a Fast-Bayesian model through a Laplace approximation in 2017 [START_REF] Au | Operational Modal Analysis, Modeling, Bayesian Inference, Uncertainty laws[END_REF]. MLE has many benefits as it is asymptotically unbiased and efficient (with respect to the Cramer-Rao lower bound) and convergent.

OMA models generally consider linear time-invariant systems, which results in the fundamental frequency input-output relation suggested in eq. [START_REF] Presas | Fatigue life estimation of Francis turbines based on experimental strain measurements: Review of the actual data and future trends[END_REF].

N / × N / , m × m are respectively the discrete output and input spectral density matrices at radial frequency and 2 m × N / is the frequency response function matrix containing the system characteristics.

N / × N / is the noise spectral density matrix. The equivalent time domain model can be obtained through the Wiener-Khinchin theorem. The objective is to determine 2 for a given with unknown. This lack of knowledge is tackled by considering the ambient vibration as white noise having a constant power spectral density. Under this assumption, becomes a constant matrix, which is sufficient to estimate the poles of 2 .

= 2 * 2 " + , ∈ 1, (1) First introduced by K. Janssens et al. in 2006 [31], OBMA is a method combining advanced order tracking and an OMA algorithm to obtain modal characteristics. It has been extensively studied by E.D. Lorenzo [START_REF] Lorenzo | Industrial applications of advanced modal identification on operational rotating machineries[END_REF][START_REF] Lorenzo | Gear Dynamics Characterization by Using Order-Based Modal Analysis[END_REF]. By considering hydroelectric turbines as multi-sine sweep generators (especially, the Rotor-Stator Interactions, RSI), resonances are easily extracted from the signal using Order Tracking (OT) methods. OT tools are capable of extracting the harmonic content of a mixed signal. Once synchronous resonances are extracted, excited modes are identified with OMA. Traditionally, the p-LSCE (or Polymax in its industrial name) identification algorithm is used [START_REF] Peeters | The PolyMAX frequency-domain method: a new standard for modal parameter estimation?[END_REF]. This parametric, deterministic frequency domain procedure is a poly-reference extension of the Least Square Complex Frequency algorithm (LSCF). It is a two-stage least square procedure that first identifies frequencies and damping ratios, and then identifies mode shapes based on stable poles. The processing of synchronous resonances with OMA is possible because the spectral density matrix of chirps in the order domain turns out to be constant, thus making it equivalent to an 'order domain white noise' (see, e.g., Figure 1). The analysis of chirp excitations, rather than stochastic excitations, is more representative of transient regimes. Especially, they provide more intense excitation forces, and resonances with higher Signal-to-Noise ratio (SNR).

In the current context, the p-LSCE has several limitations. First, the classical formulation of the problem does not consider uncertainties, which is crucial for treating poorly conditioned signals. Furthermore, its accuracy in presence of small data is shown limited in section 4, making the algorithm improper for short-time signal processing. For these reasons, the Fast-Bayesian approach is chosen to experimentally identify turbine modal characteristics. Indeed, this framework is well suited for such identifications. Typically, the identification is performed in the frequency domain, as required for an order-based formulation. Moreover, it is well adapted to process short series, since likelihood is effective for making inference with small data sample available, and will be shown robust against short-time estimation biases. Short-time estimations result in non-averaged spectra and excitation variance, that can make identification fail as the modal force is not perfectly constant.

Bayesian order-based modal analysis

This section introduces a new OBMA formulation based on probabilistic modal analysis through Bayesian inferences. Each subsection explains one step of the process, successively:

• The resonance detection.

• The angular domain order tracking.

• The identification scheme.

Resonance Detection

The data consists of an asynchronous time series containing ( strain gauge measurement channels. The first step is to map all synchronous resonance regions contained in the dataset. This is done by investigating time-frequency spectra, using short-time discrete Fourier transforms. To ensure that the observed amplification regions truly are resonances, a dual-channel Phase-Shift Analysis (PSA) might be considered as it has been proven effective for diagnosis purposes [START_REF] Leonard | Phase spectrogram and frequency spectrogram as new diagnostic tools[END_REF][START_REF] Moisan | Self-excitation in Francis runner during load rejection[END_REF]. PSA examines the shorttime phase-shift diagram between two redundant sensors, as a complement to classical amplitude spectrum studies. This is relevant because modal phase shifting between sensors is constant. Indeed, rotational symmetric structure's modes are cylindrical or disk-like patterns, characterized by their number of nodal diameters 5. Nodal diameters are diameters along which modal displacement is null, thus separating isophase regions [START_REF] Dorfler | Flow-Induced Pulsation and Vibration in Hydroelectric Machinery[END_REF]. They can be positive or negative, the sign indicating the direction of rotation: forward for 5 6 0, backward for 5 8 0. The motion along a mode has a constant phase shift: between two consecutive blades separated by an angle of 9: = 2</> ? , the phase-shift is 59:. Then between N blades, the modal phase-shift is given by: 59: 2< , where . is the modulo operator, ensuring a result on the trigonometric circle.

Angular Domain Order Tracking

In this study, Angular Domain (AD) Order Tracking is used because it does not need phase reference and tackles the spectral leaking issue [START_REF] Lorenzo | Industrial applications of advanced modal identification on operational rotating machineries[END_REF][START_REF] Lorenzo | Gear Dynamics Characterization by Using Order-Based Modal Analysis[END_REF]. A Computed Order Tracking (COT) algorithm is performed to transform an asynchronous time series into a synchronous time series. This is based on a two-stage interpolation as follows:

1. Determine the time-angle relation by integrating the rotational speed. This result in a relation with constant 9A between each angular sample. 2. Interpolate this time-angle relation to obtain a new time reference with constant 9:.

3.

Interpolate the data with respect to the new time vector. The obtained data are synchronous, in the sense they are sampled at constant 9:. Then, a non-windowed short-angle Fourier transform is performed on the resampled signal to exhibit its frequency content. For the stake of convenience, each interval corresponds to one runner revolution such that the spectrum resolution coincides with orders. An order spectrum is obtained for each studied revolution. Each line of the spectrum carries the amplitude of a synchronous harmonic in the signal, that are called orders. For each revolution, an exciting radial frequency is defined by averaging the harmonic frequency over the lap.

Bias in AD comes from the interpolation method and synchronous interval sampling. The interpolation bias is induced by the interpolating method, and residual shaft torsion vibrations that can induce tachometer signal fluctuations. This results in a slight shift of the angular signal. Concerning the synchronous interval splits, each interval must represent exactly one revolution, which is not necessarily the case: the angular shift between the first and the last point of the Fourier interval can be slightly different than 2<, and depends on the angular resampling rate. Most such errors are asymptotically unbiased, and can be reduced if the data are recorded with a very high sampling frequency as compared to the structure's natural frequencies. Also, the harmonic frequency averaging over a revolution involves a quasi-static assumption, reasonable for low accelerating systems only.

SDoF Fast-Bayesian Inference

Hydroelectric runner modal frequencies are often close to each other, and might involve coupledmode responses. However, harmonic excitations are well shaped to separate modes, since they can excite only some specific runner modes (see, for instance, the RSI theory presented in section 5). The use of a SDoF model is then advocated as it reduces the model complexity and related uncertainties. The bandwidth selection can be based on modal coherence techniques [START_REF] Pastor | Modal Assurance Criterion[END_REF] or sensitivity to bandlimiting. The choice made here depends on the signal-to-noise ratio (SNR): Modal Assurance Criterion (MAC) should be preferred for high SNR, while sensitivity analysis is convenient for low SNR, since MAC becomes ineffective.

Once the SDoF assumption is assessed, the modal parameter vector estimate is determined using a Maximum Likelihood Estimator. This approach was initially developed in the early 2000's by Yuen [START_REF] Yuen | Bayesian fast Fourier transform approach for modal updating using ambient data[END_REF]. A comprehensive state of the art is given in [START_REF] Yuen | Bayesian methods for Structural Dynamics and Civil Engineering[END_REF]. The underlying assumption is that the posterior probability is proportional to the likelihood, which is an unscaled version of the Baye's theorem with uniform priors. The experimental frequency data B C N / × 1 is independent zeromean complex Gaussian, and the likelihood function reads as:

DB C| E = F < GH I | | -J- * GL M H N OL (2)
The likelihood is a function encapsulating the statistics of the data. For any parameter vector, it gives the conditional probability of observing the data given the physical model.

= P{ * } is the theoretical density matrix. It derives from the resonance model given in eq. [START_REF] Gagnon | Impact of startup scheme on Francis runner life expectancy[END_REF][START_REF] Goyal | Review of hydrodynamics instabilities in Francis turbine during off design and transient operations[END_REF], where = ) , * ) , + ) , +-) " is the modal parameter vector, including the natural radial frequency, the damping ratio, the modal excitation PSD and the noise PSD. The transition from eq. ( 1) to eq. ( 3) is reminded in Appendix A. ' N / × 1 is the partial mode shape and the amplification factor for displacement data. N / × N / is the identity matrix.

= + ) ' ' " + +-) H I (3) = 1 D ) S -S E S -2* ) ) S (4)
The main obstacle with the method is that the formulation given in eq. ( 2) is highly non-linear, and its maximization requires intensive and time-consuming computations. To reduce the complexity and some ill-conditioning problems, a modified version of the Negative Log-Likelihood Function (NLLF) was recently introduced [14, 39]. Based on an eigenvalue decomposition, this formulation allows to express ' as a function of the other modes, reducing the optimization problem to the minimization of eq. ( 5), in which only 4 parameters need to be estimated. For typical data sizes, this process is fast, as the NLLF has a unique minimum in the parameter space [START_REF] Au | Operational Modal Analysis, Modeling, Bayesian Inference, Uncertainty laws[END_REF].

, = ( ln < + V ln + ) + +-) H N OL + ( -1 ln +-) + 1 +-) WV * -' " ' H N OL X (5)
= V Y1 + +-) + ) Z GL ℛℯJ * M H N OL (6)
A good Gaussian approximation of the posterior densities relies on a second-order Taylor expansion of the NLLF. Since such computation requires intensive efforts, posterior distributions are approximated with the Laplace method, i.e., with normal distributions centered at the MPV with standard deviations deduced from the covariance matrix , obtained with eq. ( 7). This gives the Fast-Bayesian method which leads to conditional distributions, but is faster, and any loss of accuracy it provides as compared to the unconditional distribution are considered insignificant [START_REF] Au | Operational Modal Analysis, Modeling, Bayesian Inference, Uncertainty laws[END_REF]. In other terms, the likelihood around the MPV is assumed Gaussian and independent, which is apparently reasonable for large data processing.

= ]^_ D E. J] S ,D EM . ]^_ " D E (7
)

^_D E = ` , ' a b|' a|b c (8)
In eq. ( 7), is obtained as the inverse of the NLLF Hessian at the MPV, through a constraint function that takes into account the mode shape normalization, denoted ^_ and expressed in eq. ( 8). The solution is a "pseudo-inverse" computation in the sense it ignores the null eigenvalues arising from the constraint singularities.

Numerical Testing

This section verifies the model performance through several numerical studies. At first, the model is tested on a MDoF synthetic case with known masses, stiffnesses and damping, subjected to sinesweep excitations. The synthetic case is also identified using the p-LSCE, to demonstrate the benefits of the Bayesian framework. Then, a deeper study will highlight the influence of short-time data and non-averaged spectra on the identification accuracy. Subsequently, the modal force frequency modulation (FM) rate is studied. It is shown that high FM can distort the response and introduce epistemic uncertainties. Along this section, the synthetic data is obtained with a time domain statespace model solicited by chirp-like forces.

Illustrative computational example

The model is now used on a simulation example, and compared to the p-LSCE. Let consider a 3-DoF system with known modal properties, and d = 10 measured DOF. Three modes are expected at 24.7Hz, 31.2Hz and 39.7Hz, with damping ratios 2.7%, 5.3% and 4.5%, respectively. The system is first computed in the time domain using a state-space model. The sampling frequency is 2.5kHz, making the OT biases insignificant. The AD resampling is 2,500 samples per revolution. Two simulations are performed. On the first hand, the system is excited with a pure and low sine-sweep, capturing a well-shaped response. On the other hand, the system is excited by a mix of a sine-sweep and a stochastic source, which embodies the non-flat excitation behavior. The scanning harmonic is faster, delivering short-time spectra and ill-conditioned responses. The two simulations are presented in Figure 2a. and b. left, depicting the excitation input and response output time-frequency distributions. The responses are extracted using OT, and processed with the p-LSCE and with the Fast-Bayesian inference. The p-LSCE stabilization diagram and the Fast Bayesian bandwidth selection are presented on the right side of Figure 2. The p-LSCE converges fast for the well-shaped responses, and the three modes are stable from order 15. Notice that the stabilization diagrams are very clear and easy to interpret, since the spurious modes are computed with negative damping. However, the p-LSCE seems inadequate for processing ill-conditioned responses: the first and second modes are identified at low orders, but the third mode is very unstable, and the algorithm hardly identifies a stable pole around order 70. A stable spurious mode is present around 20Hz, between order 65 and 75. The Fast-Bayesian bandlimited inputs are obtained using the MAC criterion for the first simulation, and using a sensitivity analysis for the second simulation. The identification results are presented in Figure 3. When dealing with long time histories, the two approaches are consistent and identified modal parameters are close to the theoretical values. The Bayesian approach provides estimates with low uncertainties. The short time history case is properly identified by the Fast-Bayesian scheme, but the related uncertainties are very large, especially for the damping estimates. The p-LSCE finds accurate natural frequencies for modes 1 and 2, but fails at properly identifying the third natural frequency. The damping ratios are systematically underestimated, and the identified modes 2 and 3 are almost undamped. The processing of short-time spectra cannot be handled correctly with the p-LSCE. The stabilization diagram indicates over-fitting at high order that can be sufficient for capturing stable poles, but the results show that the damping ratios are not representative of the studied dynamical modes. Conversely, the Fast-Bayesian is a legitimate algorithm that can be used to study short-time responses.

As demonstrated in this section, the identification is still accurate with sparse data and the error between the MPV and the theoretical values is low. Nevertheless, the computed uncertainties grow when the data becomes limited. Although this is naturally expected, the posterior deviation can be wide and not representative of the real parameter variability. Consequently, these uncertainties can make the experimental characterization unusable in an industrial context.

Robustness against non-flat excitation spectra

Once the algorithm is shown adequate for processing short-time histories, it is necessary to assess in which extent the identification remains acceptable. For this, it is proposed to test the ML-OBMA on a linear sine sweep excitation passing through a resonance, with different levels of deviation from a constant excitation force. To be consistent with experimental data, a set of d = 10 channels with equally distributed mode shape ratios is considered. An ML-OBMA was performed on several datasets with different variance intensity levels. The variable control is the Stochastic Gain Ratio (SGR), defined as the power ratio of the white noise deviations (the 'stochastic' part) and the averaged white noise. The FM rate was set to 6 . 10 GS fg/h, and the modal frequency and damping are 50fg and 2.5%. It is observed that even for noise with a high SGR (near 500), the theoretical natural frequency still lies within a 95% credibility interval. When the stochastic gain ratio is lower than 100, the error between m n and m n drops below 1%. The damping ratio's MPV only converges toward the theoretical value for very low SGR, and tends to be underestimated for SGR greater than 20. This is explained by the modal force PSD standard deviation discrepancy, as shown in Figure 4.b. The underestimation of the damping ratio while the excitation variance is too high is in keeping with the well-known damping underestimation in OMA. The damping and radial frequency variance seem relatively insensitive to the stochastic excitation amplitude. This is because the modal force PSD and noise PSD uncertainties increase drastically with the SGR, dominating the other modal parameter variabilities. The MPVs are not displayed, but the noise PSD increases linearly with SGR, i.e., the model equates a fraction of stochastic excitation with channel noise. On the other hand, the modal force PSD decreases to conserve signal power, and is reduced by 50% for a SGR of 100.

The model appears to be robust against stochastic excitation processes, which is an important factor for ambient measurement processing. The identification fails as soon as the modal force PSD and noise PSD standard deviations begin to diverge. The SGR at which these deviations diverge is a function of the number of channels, the higher the number of channels the lower the critical SGR.

Sensitivity to the frequency modulation rate

Two aspects of the sine sweep excitation frequency modulation rate need to be considered. The first deals with model uncertainties, while the second undermines the structure of the non-stationary natural response. With a low frequency-modulated excitation o, uncertainties appear to be low and MPVs are well estimated. When the FM rate increases, uncertainties seem to increase linearly with it, as shown in Figure 5: when τ increases, resonances develop and vanish faster. This shrinks the time duration of the signal, as shown in the lower part of Figure 5 (N q decreases as the time duration decreases). The model accuracy is inevitably reduced because fewer observations are available. This remark is known as noise mitigation process. Damping estimations remain good, but the frequency MPV tends to be overestimated during the run-up (and underestimated during run-downs), the steeper the slope the higher the error.

For higher FM slopes, another issue may arise due to dynamic distortions of the classical steadystate response. These are due to the inability of the system to instantaneously match the excitation. Such distortions have been widely investigated for electronics, and more recently, in the field of mechanical engineering [START_REF] Markert | Analytically based estimation of the maximum amplitude during passage through resonance[END_REF]. The ISO-7626 standard defines some criteria under which the response distortion is low enough to be neglected. For example, the condition for linear FM excitation is h × τ stu 8 3.6 f n ξ S Hz/s . When these conditions are not fulfilled, the response may be distorted and cannot be modelled with a static transfer function. Figure 5 (right side) demonstrates an example of such a distorted response: epistemic uncertainties can be encountered when dealing with out-ofstandard systems, and natural frequencies deviate from their original value. At the same time, damping ratios are unpredictable if distortions are not modelled properly.

Dynamic distortions can also play a role in the uncertainty rise seen in Figure 5 (left), as it introduces a modeling error leading to the expansion of credibility intervals. Identification in the presence of such dynamic behaviors was not considered in this paper, but authors seemed to deem it worthwhile to mention its existence and potential effects on parameter identification.

Case study: an experimental analysis of a Francis turbine during coast-down 5.1. Turbine characteristics and experimental setup

The case study is based on a vertical medium head Francis hydroelectric turbine. It is part of a run-ofthe-river plant located in Quebec, Canada. The turbine has a specific speed of { | = 57, and a synchronous speed of 163.6 RPM. The flow in the penstock is directed into a spiral casing, and fed into the runner through a distributor of > ~= 20 guide vanes. The runner is a Francis runner of diameter 4.4m and height 2.7m, and is composed of > ? = 13 blades. For this kind of geometry, 7 cylindrical mode shapes are expected. Such patterns are cyclic-symmetric with nodal diameters |5| between 0 and 6.

In many medium-to-high head hydraulic turbines, the main harmonics of the rotating speed are the Rotor-Stator Interactions (RSI) [START_REF] Tanaka | Vibration Behavior and Dynamic Stress of Runners of Very High Head Reversible Pump-turbines[END_REF]. The RSI theory can be jointly used with Phase-Shift Analysis to determine the mode shapes. Indeed, RSI pressure waves have specific and predictable patterns and frequency, which will perform a selective excitation that only excite matching mode shapes. Typically, the RSI seen from the runner generate excitations at the wicket gate passing frequencies:

m •d€ • = {> ~mn .
The RSI is capable of exciting modes with 5 = '> ? -{> ~ nodal diameters on the runner. Table 1 summarizes the runner RSI properties, that will be used latter when identifying the mode shapes. The phase-shift is calculated using section 3.1. 2.42 Data were recorded during an asynchronous coast-down transient experiment in no-load conditions. Two blades separated by an angle of 111° were instrumented with strain gauges, as presented in Figure 6.a. The blade intrados were instrumented with three strain gauge rosettes located in the band junction to a blade leading edge and trailing edge, and in the middle crown-blade weld, as shown in Figure 6.b. The blade extrados were instrumented with two uniaxial gauges, one close to the crown, and the other close to the band. The locations were the same from one blade to the other. This enables the turbine runner mode shapes to be determined through a phase-shift analysis. The rosette and uniaxial gauges were oriented in agreement with the expected principal strain sensor.

Mode Shape Identification

The analysis of experimental amplitude and phase-shift spectra between redundant sensors was conducted first. Examples of time-frequency distributions are shown in Figure 7, as related to the crown extrados. Phase-shift spectra are amplitude-filtered and show only phases associated with a high enough amplitude. With a long enough Fourier windowing, amplitude spectra show five operating deflection shapes (ODS).

The ODS of interest are those triggered by synchronous harmonics. It should be recalled that a harmonic interaction requires both a matching frequency and a compatible nodal diameter. As expected by the RSI theory, a ND5 is excited by the 60-th harmonic around 60Hz when the rotating speed is around 1Hz (see Table 2). A ND0 is excited by the 13-th harmonic with a zero phase-shift. This correspond to the Blade Passing Frequency seen from the stationary cascade. Such frequency can be excited either by by an improper flow distribution, a local default on the structure or a casing asymmetry. Then, the ND0 can be excited every time a blade passes through the pressure fluctuation, i.e. 13 times per revolution.

In the studied turbine, the authors were surprised to find out that some resonances are not induced by the RSI harmonics, but by harmonics in their vicinity. As shown in Figure 8, two resonances are induced by the 41-th and 62-nd harmonics while the RSI theory predicts excitations with the 40-th and 60-th harmonics. While excited mode shape nodal diameters should be ND-1 and ND5, the observed phase-shift indicates that the excited nodal diameters are respectively ND-2 and ND3, as reported in Table 2. Without going further in details, these observations are the result of a stationary point force excitation. Readers should refer to the works of Wildheim for deeper explanations [START_REF] Wildheim | Excitation of rotating circumferentially periodic structures[END_REF].

Modal identification using ML-OBMA

Each resonance is processed with an ML-OBMA to extract both frequency and damping distributions. Bandwidth were selected using a MAC analysis for the ND0, while ND-2, 3 and 5 were bandlimited using a sensitivity analysis. An example of such an analysis for ND3 is provided in Figure 9. A band of 9m = 8fg is chosen as modal parameters stable and uncertainties low. Identification results are presented in Figure 10. Two types of result were investigated: a direct ML-OBMA inference as described in the previous sections and the distributions after taking into account short-time noise variance through Markov chains. 10 † process realizations were needed to ensure a good statistical convergence. Results show that the sole likelihood inference provides robust and accurate uncertainty bounds while Markov chains do not bring any significant change. This means that the Fast-Bayesian is not sensitive to non-flat noise PSD. ND0 distributions remain unchanged after experimental uncertainties are calculated. This is explained by the high quality of the signal, which causes the noise to have almost no impact. Frequency distributions are narrow distributions, almost Gaussian with small standard deviations, while damping distributions seem skewed, developing a wider tail for upper damping values. This is particularly observable for ND3. A high uncertainty exists for damping estimation, much like with many modal analysis extraction tools. Table 3 presents the parameter statistical distributions. The most probable values (MPV) and coefficients of variation (ratio of standard deviations over MPVs) are drawn.

The deviation of the estimated mode shape from the real mode shape can be estimated, without the latter being known. This is done using a stochastic representation of the mode shape: ' = ‡‖ ‡‖ G‰ , ‡~‹ ' a, ' , and using an eigen-representation of the MAC. For more details please refer to [START_REF] Au | Operational Modal Analysis, Modeling, Bayesian Inference, Uncertainty laws[END_REF]. ML-OBMA provides very close estimations of mode shapes, since the modal assurance criterion is high with a low c.v. As shown in Figure 11, the MAC may be close to 1, mode shape can strongly deviate from the expected target (ND5), or slightly (ND0). The MAC value does not reflect this inaccuracy. The MAC c.v. gives the global degree of adequacy, accounting for the norm constraint ‖'‖ = 1. This degree of uncertainty depends on the SNR, the higher the SNR, the lower the error: ND0 has an SNR of 7.4OE• and a MAC c.v. of 2.7%, while ND5 has an SNR of 0.95OE• with a MAC c.v. of 7.7%.

Comparison with numerical simulations

Once the experimental frequencies are characterized, it is possible to use them to validate the numerical predictions. In order to compare the experimental results to the numerical computations, several simulations were performed to obtain the runner modal analysis in different configurations. As it is well known that the modal behavior of a turbine is strongly influenced by the surrounding water, the fluid-structure interaction (FSI) added mass effect must be considered in the simulations when computing the natural frequencies. The full runner geometry is considered, and the meshing for one blade is presented in Figure 12, with tetrahedral elements. The shaft is not included in the simulation, which could have a significant influence on the prediction of nodal diameters 0 and 1. Instead, a zeromotion boundary condition is imposed at the shaft coupling ring. The FSI are treated like a structural-acoustical problem in Ansys. Three different water volumes were modelled, as shown in Figure 13. In the first simulation, a large volume of water is studied with a fluid domain of 0.5m around the runner. In the second simulation, a restricted volume of water is used with 0.1m distance between water external boundary and runner largest diameter, to take into account the runner confinement effect. In the third simulation, a more realistic geometry with a larger fluid domain is considered, including the upper and lower labyrinths real gaps (1-2mm), a portion of the distributor inflow, and the upper part of the draft tube. Obtained fluid volumes are meshed using quadratic elements with approximately 2M nodes. The external fluid bounds are treated as reflective walls. In the simulations, the system is considered undamped and no damping prediction is available. For this reason, the comparison with the experimental results is restricted to the natural frequencies. 4 presents the relative error between the identified frequencies and the numerical computations. The ND0 is a torsion mode, which is strongly influenced by the shaft. Since the finite element model did not incude the shaft, the large overestimation of the modal analysis is easily explainable. Theoretially, nodal diameters larger than 2 are fully balanced over the runner and they should not be influenced by the shaft. Consequectly, the modal prediction is very close to the experimental observations. From simulation 1 to 2, the error decreases by 10% for modes 2 and 3. The best results are naturally obtained from the complex fluid domain geometry used in simulation 3, with a remaining error of about 1%. This is mainly due to the confinement effect, which is known to reduce the natural frequencies.

Although experimental identification is incomplete and misses several modes, available results can be used to calibrate models and obtain more precise load levels throughout the turbine. This calibration will also have impact on non identified modes, what will result in an improvement of stress predictions, then turbine fatigue and structural health.

Conclusion

A new ambient modal analysis tool was implemented, based on the combination of order tracking and the so-called Fast-Bayesian algorithm. The ML-OBMA outperforms the OBMA in the sense that not only is able to identify modes, but it can also quantify related uncertainties, which can be relevant in many practical cases. Furthermore, it was demonstrated that the Fast-Bayesian is more adequate for analyzing limited data than the traditional Polymax. The model is unable to process synchronous resonances if spectra are influenced by modulation distortions, but it appears that such distortions involve FM rates that are very unlikely to occur in hydroelectric turbines. It is relatively robust against the flat spectrum assumption, and the inference still provides accurate estimates in presence of large residual excitation variance. This is an important characteristic because short-time responses systematically exhibit non-flat excitation behaviors. However, one should keep in mind that damping ratios tend to be underestimated in some cases where residual variance power is too high.

This study allows a statistical modal analysis through the investigation of synchronous harmonic resonances. Experimental parameter estimates can be used for many purposes in industry, including, for instance, numerical model validation, trade-off criteria to optimize model accuracy versus computing time, and model calibration. In this paper, the model was implemented on a prototype Francis turbine, and several resonances were featured, showing a close agreement with the numerical modal analysis. In future works, the modulated RSI phenomenon will be studied deeper. Such behavior is rare and absent from the literature, although it can produce unexpected resonances that may damage the structure.

In the Fast-Bayesian model, parameter distributions were assumed to be linearly dependent on the likelihood function, i.e., the priors were implicitly assumed to be uniform. The Bayes theorem was then used in its weakest version. Furthermore, the posterior distributions were approximated in a Taylor manner, using a Laplace approximation. Future works should transform this weak Bayesian theory into a fully Bayesian approach, by considering non-uniform priors and a sampling of the likelihood in the vicinity of the MPV.

7. Acknowledgement This study was funded by Mitacs Acceleration, grant IT10516, in collaboration with the research institute of HydroQuebec (IREQ) and Andritz Hydro Ltd.

Figure 1 :

 1 Figure 1: Time series, frequency content and order-tracked frequency response of a pure sine sweep and a harmonic signal buried in noise.

Figure 2 :

 2 Figure 2: 3-DOF system response. a) well-shaped responses. b) ill-conditioned responses.

Figure 3 :

 3 Figure 3: p-LSCE and Fast-Bayesian identification.

Figure 4 .

 4 a presents the sensitivity results. The markers draw the normalized most probable value (scaled on the theoretical value), and the bars represent k3 normalized standard deviations.

Figure 4 :

 4 Figure 4: a) Modal parameter most probable values as a function of SGR. Uncertainty bars cover 3 c.v. b) Modal force and PSD error standard deviations. The two distribution parameters evolve jointly and tend to when the noise intensity is too high.

Figure 5 :

 5 Figure 5: Left -Uncertainties increasing with FM slope. Right -The dynamic distortion of modal response when subjected to strong frequency modulated excitations.

Figure 6 :

 6 Figure 6: a) Instrumented blades. The two blades are separated by an angle of 111°. b) Strain gauge location on the instrumented blades. Circles represent rosette gauges, triangles, uni-axial gauges.

Figure 7 :

 7 Figure 7: ODS analysis of principal direction of extrados crown Rosette gauge.On the left, phase-shift spectrum of the redundant gauges. On the right, redundant amplitude spectra.

Figure 8 :

 8 Figure 8: Unexpected resonances with the 41-th and 62-nd harmonics.

Figure 9 :

 9 Figure 9: ND3 bandlimiting using a sensivity analysis.

Figure 10 :

 10 Figure 10: Frequency and damping distributions after ML-OBMA and after Markov chains ML-OBMA.

Figure 11 :

 11 Figure 11: Mode shapes for identified modes.

Figure 12 :

 12 Figure 12: Turbine mesh for one blade.

Figure 13 :

 13 Figure 13: FEM analysis for the simulations with acoustic fluid.

Table 1 :

 1 RSI characteristics in the runner rotating coordinate system.

	RSI harmonic { 1 2 3 4	Frequency m •d€ • 20 m n 40 m n 60 m n 80 m n	Excitation shape 5 6 -1 5 -2	Expected phase-shift [rad] -0.97 -1.93 -2.90

Table 2 :

 2 Excitation shape and frequency domain of the identified resonances.

	Resonance Band	Rotating speed	Phase-shift [rad] Mode shape nodal
	[15, 20] Hz [20, 40] Hz [40, 50] Hz [50, 60] Hz	[Hz]	1.3 0.7 0.8 1.0	0 2.42 -0.48 -2.90	diameters |5| 0 -2 3 5

Table 3 :

 3 Modal parameter distributions for identified modes.

		M.P.V. (SI)	C.V.	M.P.V. (SI)	C.V.
		ND0		ND3
	Frequency [Hz] Damping ratio [%] Excitation PSD [Ž+] Noise PSD [Ž+/fg] Shape (MAC)	17.43 1.26 2.2 × 10 • 1.7 × 10 GS 0.999 ND-2	3.8 × 10 G † 3.0 × 10 GL 1.6 × 10 G † 2.6 × 10 G' 2.7 × 10 GS	49.79 1.39 1.2 × 10 • 7.6 × 10 G' 0.999 ND5	5.2 × 10 G † 5.4 × 10 GL 7.3 × 10 G• 1.4 × 10 G• 4.5 × 10 GS
	Frequency Damping ratio Excitation PSD Noise PSD Shape (MAC)	28.72 2.75 3.3 × 10 • 1.2 × 10 G † 0.999	6.4 × 10 G † 2.8 × 10 GL 9.2 × 10 G• 1.6 × 10 G• 5.4 × 10 GS	59.25 2.01 4.3 × 10 • 2.0 × 10 G † 0.997	5.2 × 10 G † 2.4 × 10 GL 7.8 × 10 G• 5.9 × 10 G• 7.7 × 10 GS

Table 4 :

 4 Comparison of experimental and numerical natural frequencies of identified modes.

	Mode shape nodal	Experimental	Simulation	Simulation	Simulation
	diameters |5| 0 2 3 5	freq. [Hz] 17.43 28.72 49.79 59.25	1 232% 36% 16% 3%	2 211% 22% 9% 0.5%	3 180% -2% 1% -0.5%

=

L/S '〈"〉 " (10) where ' is the mode shape and 〈"〉 is the modal participation factor. Then, eq. ( 1) can be written as:

Without further development, it appears that 〈"〉 • 〈"〉 is a scalar independent of the frequency. Defining the scaled modal force as + = 〈"〉 • 〈"〉, and assuming to be constant diagonal (i.i.d white channel noise under homoscedastic assumption), then: = + '' " + + -H ˜, ∈ 1, (12) which is eq. ( 2).