A Cellier
email: cellier@cerfacs.fr

C J Lapeyre

G Öztarlik

T Poinsot

T Schuller

L Selle

Detection of Precursors of Combustion Instability using Convolutional Recurrent Neural Networks

Keywords: Thermoacoustic Instability, Instability Precursors, Deep Learning, Convolutional Recurrent Neural Networks 1

Many combustors are prone to Thermoacoustic Instabilities (TAI). Being able to avoid TAI is mandatory to efficiently operate a system without sacrificing neither performance nor safety. Based on Deep Learning techniques, and more specifically Convolutional Recurrent Neural Networks (CRNN) 1 , this study presents a tool able to detect and translate precursors of TAI in a swirled combustor for different fuel injection strategies. The tool is trained to use only time-series of unsteady sensors in stable conditions to predict the proximity of unstable operating points on a mass flow rate / equivalence ratio operating map, offering a real-time information on the margin of the system versus TAI. This allows to change operating conditions, and detect the directions to avoid in order to remain in the stable domain.

Introduction

The way towards leaner and fuel flexible combustion is paved with critical issues [1,2]. Thermoacoustic Instabilities (TAI) threaten the integrity of the system and the safety of its user. Their prediction is crucial to further improve combustion devices. This work highlights a step toward the implementation of a Deep Learning tool for the prognosis of instability in a ducted lean multi-fuel swirled injector. Based on a real-time analysis of combustion noise, it predicts the proximity of unstable operating points in a bulk velocity and equivalence ratio diagram.

The framework of this task is a well known field of research in thermoacoustics known as "Precursor Detection". Several groups have highlighted that precursors of TAI are observed as the system gets closer to bifurcation [3][4][5][6]. However identifying a robust definition of precursors remains difficult. As pointed out by Juniper and Sujith [START_REF] Juniper | Sensitivity and Nonlinearity of Thermoacoustic Oscillations[END_REF], strategies to extract precursors for experimental systems are still challenging to set but robust stable signal analysis for instability prediction should be found in the near future. In fact, the existence of tangible information in combustion noise pressure has been discussed by Lieuwen [8] on a laboratory-scale combustor. An et al. [START_REF] An | Early Warning Signs of Imminent Thermoacoustic Oscillations Through Critical Slowing Down[END_REF] observed that pressure signals can warn about impending instability in a turbine engine. Rouwenhorst et al. [START_REF] Rouwenhorst | Online Monitoring of Thermoacoustic Eigenmodes in Annular Combustion Systems Based on a State-Space Model[END_REF] successfully produced an online strategy to obtain stability margins based on modal analysis of time-series for an annular combustor. Murayama et al. [START_REF] Murayama | Characterization of dynamic behavior of combustion noise and detection of blowout in a laboratory-scale gas-turbine model combustor[END_REF] and Gotoda et al. [START_REF] Gotoda | Detection and control of combustion instability based on the concept of dynamical system theory[END_REF] have also shown that specific dynamic behaviors can be detected in noise, preceding a critical event such as lean blow-off, and proposed methods to identify precursors for such events. Yet, finding the right precursors to focus on for one specific device can be a thorny problem. Coupling Machine Learning techniques to pre-existing methods eases the identification and detection of precursors.

Sarkar et al. [START_REF] Sarkar | Early Detection of Combustion Instability from Hi-speed Flame Images via Deep Learning and Symbolic Time Series Analysis[END_REF] have implemented real-time diagnostics based on Deep Learning applied to symbolic time series analysis of high-speed camera imaging. Furthermore, Kobayashi et. al. [START_REF] Kobayashi | Early Detection of Thermoacoustic Combustion Instability Using a Methodology Combining Complex Networks and Machine Learning[END_REF] highlighted then the possibility to extract meaningful information from pressure and heat release rate measurements for critical event prediction. Evidences of the benefits of Machine Learning techniques are also observed when affordable sensors are used. Sengupta et al. [START_REF] Sengupta | Bayesian Machine Learning for the Prognosis of Combustion Instabilities from Noise[END_REF], through the analysis of a single pressure sensor, obtained insights on the stability margins of the system. Finally, McCartney et al. could couple Detrended Fluctuation Analysis to Machine Learning and outperformed former methods on the detection of precursors of instability from dynamic pressure measurements [START_REF] Mccartney | Online detection of combustion instabilities using supervised machine learning[END_REF].

These results suggest that the time signals of pressure, velocity and heat release rate can contain sufficient information to predict that the instability is close, prior to bifurcation, and that a machine learning strategy can leverage data to automatically extract it. An additional step toward a user-centered precursor detection system is taken here by proposing a Convolutional Recurrent Neural Network (CRNN). The objective of this study is to use real-time pressure, velocity and heat release rate signals of a combustor to be able to change its operating conditions, here the equivalence ratio φ and total mass flow rate, controlled by the bulk velocity U sw , without ever entering a domain of combustion instability. The tool must have been trained to recognize if a stable operating point is close in parameters (φ, U sw ) to unstable operating points, and predict in which directions the unstable points exist, offering the possibility to avoid them.

A concrete application of this tool could be explained as follows. Maps of operating conditions (U sw , φ), tested a priori, could be used to know which regions are unstable and to be avoided. However, these maps can be subject to modifications due to changes in fuel, atmospheric conditions, hardware wears, ... It is thus interesting to use real time analysis of the combustor activity to predict whether combustion instability is close or not, without relying on initial maps anymore. The practical challenge of this study is that two operating parameters are studied simultaneously (φ, U sw ), implying that precursors of instability should have different signatures regarding a notion of direction on this bi-parameter plane. In addition, different injection strategies are tested in order to enhance the diversity of the dataset. In this context, Deep Learning, a sub-field of Machine Learning, can be a precious ally as it gives the possibility to analyze high-dimensional data such as high sampling rate time signals. Data-based methods applied to physics are currently in development in many fields, following progress of software and hardware. Accessible tools for Deep Learning [START_REF] Abadi | TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems[END_REF] open the scope of possibilities. This study investigates the ability of neural networks to recognize instability precursors from time signals with good confidence and high accuracy for various fuel injection strategies, making a real-time evaluation of the stability of the burner possible.

First a description of the experimental setup and dataset is given. The learning task is then exposed with explanations about the learning model, its training and use. Finally, comments on the interpretability of the approach are presented, along with a discussion about choices made during the implementation of the method.

Experimental setup

The MIRADAS experimental setup used in this study is presented in Fig. 1 [START_REF] Oztarlik | Suppression of Instabilities of Swirled Premixed Flames with Minimal Secondary Hydrogen Injection[END_REF]. The instrumented system is composed of a swirled injector confined in a cylindrical quartz combustion chamber with an internal diameter of 46 mm and a length of 297 mm, opened at the top-end. The main mixture is injected in the lower part of the plenum of circular cross-section with a diameter of 65 mm, and homogenized by a honeycomb panel. The flow passes through a convergent and a swirler with feed channels oriented 15°off the swirler axis. Flame piloting is offered by the central injection tube dedicated to pure fuels. In this study, the configurations involve methane (CH 4 ) and hydrogen (H 2 ). Methane, hydrogen and air mass flow rates are controlled separately (Brooks SLA 585x series), making this system a modular burner for the study of multi-fuel/multi-injection control strategies.

Combustion stability is diagnosed with one Brüel & Kjaer-TYPE 4954 microphone for the measurement of acoustic pressure, located at the exit section of the convergent (M1). It is coupled with a constant temperature hot-wire (Dantec 55P16 miniature hot wire probe and a Dantec 54T42 MiniCTA) for velocity measurement (HW). The unsteady heat release rate is deduced from the light emission intensity recorded by a THORLABS FB430-10 photo-multiplier and a narrow band-pass filter centered on the emission spectrum peak of CH * radicals (PM).

The experimental dataset

The stability of the system under different operating conditions, controlled by the bulk velocity U sw and the equivalence ratio φ, has been assessed in [START_REF] Oztarlik | Suppression of Instabilities of Swirled Premixed Flames with Minimal Secondary Hydrogen Injection[END_REF] where seven different fuel injection configurations have been tested, for premixed and piloted injections of methane (CH 4 ) and hydrogen (H 2 ). The fuel injection configurations are defined in Tab. 1 and represented in Fig. 2. The influence of each injection strategy is measured in percentages of power and presented in Tab. 1. Most of the power is provided by premixed CH 4 /air combustion.

For each fuel injection configuration, the bulk velocity U sw is varied from 14 m.s -1 to 40 m.s -1 with a 2 m.s -1 step, and the equivalence ratio φ is varied from 0.70 to 0.85 with a 0.05 step. A maximum of 56 points per fuel injection configuration are acquired at a sampling frequency of 10 kHz for a 10.3 s duration. For several configurations, operating points at high bulk velocity and high equivalence ratio were not acquired as strong TAI threatened the integrity of the test bench. Furthermore, at low bulk velocity or low equivalence ratio, blow-off or flash-back phenomena were observed. This explains the variations in number of operating points acquired N op in Tab. 1.

A simple definition of the proximity of combustion instability

Stability Analysis

The database is composed of 366 different operating points, indexed by the corresponding injection configuration, the bulk velocity U sw and the equivalence ratio φ. For each operating point, the Root Mean Square (RMS) value of the acoustic pressure recorded by the microphone M1 is computed. Unstable operating points display high pressure oscillations, a main characteristic of the limit cycle during TAI [START_REF] Poinsot | Theoretical and Numerical Combustion[END_REF]. The empirical limit of 250 Pa (141.9 dB) of RMS is observed to be the threshold of TAI for this combustor. Operating points with RMS pressure fluctuations levels above the threshold are labeled unstable. Among the 366 operating points obtained, 187 are evaluated stable, i.e. with a RMS pressure level lower than 250 Pa (141.9 dB).

Specificity of the database

In the MIRADAS burner, TAI does not monotonically increase with the bulk velocity and the equivalence ratio. For example in the Ref configuration in Tab. 1, the Root Mean Square (RMS) values of the acoustic pressure recorded by the microphone M1 are represented in Fig. 3, as a (U sw , φ) map. Darker regions are subject to strong instability. It is thus possible to see stable regions framed by unstable regions in the way that it is difficult to properly define a simple rule associating a bulk velocity and an equivalence ratio to the stability of the combustor. Furthermore, when considering a different injection configuration where 1 % of the power is produced by the injection of H 2 in the pilot lance (Tab. 1, PH1), regions of instability are modified in Fig. 4 and cannot be associated precisely to a given bulk velocity and equivalence ratio. It is impossible to predict the position of unstable regions with respect to the operating parameters (U sw , φ) without running unstable.

The following paragraphs present an approach based on the sole exploration of stable operating points in order to get information on the position of these unstable operating regions.

Label attribution

In a first step, a label is attributed to each stable operating point of the operating map. This index contains information on surrounding unstable points on a (U sw , φ) plane. The label contains 4 binary digits, 1 digit per direction, indicating whether this direction points toward instability or not. An example of label attribution is given in Fig. 5. The boundary conditions for this attribution are treated such that a stable operating point located at a border is given a label 0 in the direction of outside the map. Examples are shown in Fig. 5 for U sw ∈ [START_REF] Kobayashi | Early Detection of Thermoacoustic Combustion Instability Using a Methodology Combining Complex Networks and Machine Learning[END_REF][START_REF] Oztarlik | Suppression of Instabilities of Swirled Premixed Flames with Minimal Secondary Hydrogen Injection[END_REF] m.s -1 , φ = 0.7 or U sw = 14 m.s -1 , φ ∈ [0.7, 0.75]. This procedure is done manually and produces a map of labels for the entire operating region. This label is the simplest indicator of the proximity of unstable points on a bi-parameter plane.

Definition of the objective

In order to build the labelling presented here, it is necessary to operate the experimental bench at every operating point, including unstable ones. The practical objective of this study is to see if it is possible to predict this labelling without ever entering unstable zones. Every stable operating point is associated with two sorts of information which are its position on a bi-parameter plane (U sw , φ), and the signal acquired by the sensors (M1, HW, PM).

The following study offers a comparison between two different approaches to predict the label. First of all, a baseline is built to predict the label from the treatment of the position (U sw , φ) of a stable operating point. Then, an analysis of the three-channels signal using Deep Learning techniques is performed on stable operating points in order to predict the label.

Presentation of the learning environment

State prediction in the machine learning context

Predicting the proximity of a system to instability can be translated to a Machine Learning task by defining the process to learn the association of an input X to a corresponding output y such that X ∈ IR N T ×N C is a sample of time signal provided by the sensors, for a stable operating condition, where N T is the number of points in the sampling window and N C is the number of channels considered, and y is a label translating the proximity of unstable operating points.

A supervised learning of this task is a two step procedure. First, the training dataset composed of samples of signal X i and their corresponding output label y i known exactly is built.

It contains N s samples (X i , y i ) i∈ 1,N s . Let f such that ∀i ∈ 1, N s , y i = f (X i ).
A function fθ is searched under the form of optimizable models in order to approximate f , where θ is the set of parameters of the model. A metric called the loss function is used to assess how close fθ approximates f , and the process of optimizing fθ to minimize the resulting loss function on the dataset (X i , y i ) i∈ 1,N s is called training the model. At the end of the training process, the best trained model defined by f * θ , and its parameters θ * is frozen and saved [START_REF] Goodfellow | Machine Learning Basics[END_REF]. It can then be used to associate new input samples X to predicted labels ŷ , where the exact label y is not known a priori, because the operating point has not been tested before and its proximity to instability is thus unknown. This second phase is separated from the training, and called the inference. The predicted label is computed as ŷ = f * θ (X ) which approximates y = f (X ). The learning environment can be summarized as follows. The signal (X i ) i∈ 1,N s composes the Input of the model, and (y i ) i∈ 1,N s represents the corresponding Output that the model is expected to produce. The optimizable Model f is searched in order to link efficiently the Input to the Output for further use. Eventually, a Baseline is defined in order to assert the performance of the Model. It is thus necessary to define the Input, the Output, the Model and its Baseline.

Choice of Input and corresponding Output

The 10.3 s long signal acquired at a given operating point is sampled by randomly cropping 300 windows of 0.3 s. This strategy is preferred to a classical regularly spaced sampling procedure to avoid injecting spurious patterns based on the choice of number of windows extracted. The choice of the duration 0.3 s of the sample has been influenced by two factors, the necessity to have a sufficient number of samples for the training in order to reach convergence, and the will to capture enough periods of oscillation given that the TAI frequency of the system always lies close to 590 Hz. The samples are built with three channels provided by the three sensors mounted on the experimental setup:

1. The acoustic pressure M1 upstream the swirler 2. The acoustic flux M1.HW, which corresponds to the product of the unsteady pressure and the hot wire velocity at the same location 3. The unsteady heat release rate PM Each channel of the sample is normalized by its own maximum absolute value such that the signal lies in [-1, 1]. It insures that the most obvious traces of the bulk velocity and the equivalence ratio are impossible to retrieve in the signal. Besides, it is usual practice to feed a learning model with normalized signals as amplitude can vary from a configuration to another making generalization harder. One sample of an operating point is attributed the label corresponding to this operating point. The Input is thus a three-channels 0.3 s time signal and the Output is a four binary digit label. After sampling, the dataset contains 56100 3-channels signal samples with their corresponding labels, extracted from the 187 stable operating points.

An adequate Model: the CRNN

This work aims to classify multi-channel time-series. The classification of sound-like time-series is a recurrent topic in Deep Learning. A naive approach is to use a Multi-Layer Perceptron (MLP) [START_REF] Rumelhart | Learning internal representations by error propagation[END_REF]. Built with Fully Connected layers, its architecture is simple and easy to set-up, but limitations emerge when considering high-dimensional series where coherent patterns should be aimed at. A better architecture is thus obtained by enhancing the MLP with convolution layers in order to build a Convolutional Neural Network (CNN) [START_REF] Goodfellow | Convolutional Networks[END_REF]. Based on filtering techniques, a CNN learns an optimal set of filters in the form of convolutional kernels. Finally, a potential improvement of the efficiency of common CNN for time-series is obtained by the addition of a recurrent block [START_REF] Nassif | Speech Recognition Using Deep Neural Networks: A Systematic Review[END_REF]. The Long Short Term Memory (LSTM) layer is chosen for this study for its good performance in collaboration with CNN [START_REF] Zhao | Speech emotion recognition using deep 1D & 2D CNN LSTM networks[END_REF][START_REF] Yalta | CNN-based Multichannel End-to-End Speech Recognition for Everyday Home Environments[END_REF]. The CNN generally focuses on learning patterns at different scales. The LSTM layer, when put after the convolutional blocks, learns the time dependence between the patterns [START_REF] Hochreiter | Long Short-Term Memory[END_REF].

The Convolutional Recurrent Neural Network (CRNN) is composed of four stacked convolutional blocks. Each block is made of a 1-Dimensional convolution layer followed by a Batch Normalization [START_REF] Ioffe | Batch normalization: Accelerating deep network training by reducing internal covariate shift[END_REF], a Rectified Linear Unit (ReLU [START_REF] Jarrett | What is the best multi-stage architecture for object recognition?[END_REF]) and a Max-Pooling layer. The four stacked blocks are connected to the LSTM layer. The output of the recurrent layer is finally associated to fully connected layers used to interpret the patterns and their dependence, learnt by the precedent layers and to output the predicted label. The task to perform is more specifically a multi-class-multi-label classification, multi-class as there are four different directions to assert and multi-label as one sample can be part of more than one class when several directions are simultaneously pointing toward instability. Therefore the output activation function is chosen to be four sigmoid functions, one for each direction, thus yielding four scores laying in ]0, 1[. Each score is interpreted as a predicted probability that the concerned direction points toward TAI. The network used in this study is represented in Fig. 6 and detailed specifications are given in Tab. 2.

Definition of the Baseline

Before deploying the Deep Learning signal analysis, a simpler approach has to be implemented. Comparing the performance of this simple baseline strategy to the complete analysis gives the opportunity to emphasize both weaknesses and forces of the learning approach. The Baseline is built as follows. The dataset is composed of 187 stable operating points. Each operating point is associated with a couple (U sw , φ). The objective of the baseline is to associate this couple to the label attributed to the operating point. This association relies on four Logistic Regressions (LR) taken from the Scikit-Learn Python library [START_REF] Pedregosa | Scikit-learn: Machine Learning in Python[END_REF]. Each regression is expected to predict one direction out of the four. The LR consists in optimizing w = (w 1 , w 2 ) and w 0 such that:

ŷi = f (w 0 , w, X i ) = sigmoid(w 0 + (X i ) T w) = sigmoid(w 0 + w 1 U i sw + w 2 φ i ) ( 1
)
where ŷi is the predicted probability that the considered direction of the i th sample is unstable, U i sw = U i sw /U max sw is the normalized bulk velocity of the i th sample, φ i = φ i /φ max is the normalized equivalence ratio of the i th sample, and

X i is the vector (U i sw , φ i ).
The sigmoid function completes the model by ensuring that the score ŷi lies in ]0, 1[ with a strong dichotomy between scores close to 0 and scores close to 1 in order to obtain a clearer idea of the stability of the considered direction of the i th sample. The weights w and w 0 are obtained following the minimization of the L2-regularized logistic cost function:

J = C N s i=1 log exp -2y i -1 (X i ) T w + w 0 + 1 + 1 2 w T w (2)
where C is a user defined constant and y i is the true label. In this expression, y i = 1 if the considered direction of the i th sample is unstable, y i = 0 if it is not. The minimization is performed by the L-BFGS-B algorithm [START_REF] Zhu | L-BFGS-B: Algorithm 778: L-BFGS-B, FORTRAN routines for large scale bound constrained optimization[END_REF]. This strategy only relies on the analysis of the correlation between the two operating parameters and the output labels. For example, if one operating point defined by its bulk velocity and equivalence ratio has the same label for every injection configurations, the Baseline will already perform perfectly on this operating point. The signal analysis is then useless. However if there are operating points on which the signal analysis outperforms the baseline, this would ensure that information is extracted from the signal, and could not be leveraged by the simpler Baseline. Figure 7 summarizes the learning task and the Baseline.

Training the CRNN

Splitting the dataset

In order to verify that the training is successful, the dataset is split into three distinct groups:

1. training set: The label attribution is learnt on the samples from the training set. The 187 operating points are allocated into the three sets so that the distribution of labels fits the natural distribution observed in the entire dataset, thus reducing unbalanced behaviors where only a few different types of labels can be precisely predicted. Figure 8 displays how the distribution into the sets is built. The complete set is composed of 187 stable operating points (accounting for 51600 samples), the training set has 132 stable operating points (39600 samples), the validation set has 27 stable operating points (8100 samples) and the testing set has 28 stable operating points (8400 samples). The evaluation of the Baseline is eventually made following the same dataset splitting such that performances of both methods can be compared after training.

Loss Function

The measure of the performance of the model through the epochs of training is done by the loss function. For this task, it is chosen to work with the L2-Regularized Binary Cross-Entropy (BCE L 2 ) [START_REF] Mitchell | Artificial Neural Networks[END_REF][START_REF] Goodfellow | Deep Feedforward Networks[END_REF][START_REF] Goodfellow | Regularization for Deep Learning[END_REF] . It is computed as follows for the four directions:

BCE L 2 = - 1 N s N s i=1 y i log(ŷ i ) + (1 -y i ) log(1 -ŷi ) +L 2 (θ) (3)
where N s is the total number of samples composing the evaluated set (testing or validation set), ŷi is the predicted probability that the considered direction of the i th sample is unstable, y i is the true label. In this expression, y i = 1 if the considered direction of the i th sample is unstable, y i = 0 if it is stable. The operator L 2 (θ) accounts for the L 2 norm of certain optimizable parameters θ of the model. It avoids that the parameters of the model reach too high values and smooths the training. The four BCE L 2 obtained for each direction are then averaged to produce a single value after each epoch, and this defines the loss. This metric evaluates the confidence of the predictions by strongly penalizing overconfident behaviors. For example, a high testing or validation Binary Cross-Entropy [START_REF] Mitchell | Artificial Neural Networks[END_REF][START_REF] Goodfellow | Deep Feedforward Networks[END_REF] in parallel with a high rate of true predictions can be interpreted as the fact that most of the predictions are right but when a prediction is wrong, it is far from the truth.

Training results

Over the epochs of training, the loss is monitored for both training and validation sets. Figure 9 highlights the fact that an optimum is obtained when the validation loss is minimal after 17 epochs representing a 20 minutes training on a Nvidia Tesla V100 GPU using the Keras framework with TensorFlow as backend [START_REF] Abadi | TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems[END_REF][START_REF] Chollet | Deep Learning with Python[END_REF]. After training, this optimal network, frozen at 17 epochs, is saved for further use, namely predicting the proximity of unstable operating points based on the analysis of samples never seen during the training. Results are then compared to the Baseline.

First performance evaluation: Accuracy measurement

The output of both the model and the baseline is a vector of four values in [0, 1]. Each value is interpreted as the probability that the corresponding direction is unstable. The discriminant threshold of 0.5 is applied to these values such that a probability under (resp. over) the threshold informs that the direction concerned is predicted stable (resp. unstable). The accuracy is then computed as follows:

a = 1 N s N s i=1 n i 4 ( 4
)
where N s is the number of samples in the set considered and n i is the number of directions correctly evaluated for the i th sample. For example, if the i th sample should have been attributed {up:0-left:0-right:1-down:1} but got attributed {up:0-left:0-right:1-down:0}, only three directions are correctly asserted, n i /4 = 0.75. Using the formalism introduced in order to describe the loss function, n i is written as follows:

n i = 4 j=1 1 -1 thr ŷi j -y i j ( 5
)
where ŷi j is the predicted probability that the j th direction of the i th sample is unstable, y i j is the true label, and y i j = 1 if the j th direction of the i th sample is unstable, y i j = 0 if it is stable. The function 1 thr is the indicator function set with the default discriminant threshold of 0.5. Therefore, 1 thr (x) = 1 if x 0.5 and 1 thr (x) = 0 if x < 0.5.

The accuracy of the model and its Baseline on both validation and testing sets is given in Tab 3. First of all it is important to observe that the Baseline is already performing better than randomly guessing the directions pointing toward instability, showing that some basic information is already contained into the operating parameters. Then, the very high accuracy reached by the model is a first step toward the conclusion that the model can efficiently extract information from the signal in order to distinguish the directions of instability. However, as the task is multi-classmulti-label, the sole evaluation of the accuracy is not sufficient to conclude on the performance of the model. Actually it could be possible that three out of four directions are perfectly predicted for every samples and the last one would only be poorly detected. Hence, a Receiver Operating Characteristic (ROC) study is conducted [START_REF] Fawcett | An introduction to ROC analysis[END_REF].

Second performance evaluation: ROC curves

The Receiver Operating Characteristic study is adapted to binary classification. It consists in evaluating both the false positive rate (FPR) and the true positive rate (TPR) given by the model as the discriminant threshold introduced and fixed to 0.5 in the last paragraph is moved. For each value of the threshold, the TPR is plotted versus the FPR, producing a curve bounded by a square of side 1 and linking the point {0, 0} to the point {1, 1}. The rates are obtained on the combined validation and testing datasets, composed of samples not used for training. For example, as the threshold changes, a perfect classifier would not change its results, even if the threshold is moved far from 0.5, i.e. close to 0 or close to 1. The corresponding ROC curve would limit a perfect square. The area under the curve (AUC) is then reaching 1. A random classifier produces a straight line between {0, 0} and {1, 1} with an AUC of 0.5. A good classifier is therefore characterized by a high AUC ranging from 0.8 to 1.0. It is then possible to compare two binary classifiers by observing their respective AUC.

The task performed in this document is not strictly a binary classification. Nevertheless, when each direction is considered separately, the problem becomes a set of four binary classifications, each task being to define the class of the direction: stable (0) or unstable (1). Thus, four ROC curves are obtained, one for each direction, allowing a separated evaluation of the performance of the two methods regarding each direction. Figure 10 plots these curves along with the corresponding AUC. It appears that the directions "UP" (higher equivalence ratio) and "DOWN" (lower equivalence ratio) are well classified by both methods. However, the trained model manifests higher AUC, exceeding 0.9, showing that these directions are more efficiently classified by the model than by the baseline. Results are even clearer regarding the directions "LEFT" (lower bulk velocity) and "RIGHT" (higher bulk velocity). The Baseline laboriously reaches AUCs near 0.6 when the trained model continues to perform around 0.9.

To demonstrate the use of ROC curves to set safety margins, two constraints are imposed. First, the False Negative Rate (FNR) should not exceed 5 %, which corresponds to a True Positive Rate higher than 95 %, then the False Positive Rate (FPR) should be below 20 %. These mimic system constraints such as damaging instabilities that must be avoided, and the wish to thoroughly explore the stability map, respectively. These two different constraints lead to performances represented in Fig. 11. In the first case, the FPR is contained for every direction, regarding the high expectations in terms of safety. In the second case, TPR are elevated enough for a precise exploration. The ROC study thus offers the opportunity for the user to define its tolerances regarding the specifications of the experimental setup, leading to different accuracy trade-offs.

The ROC study thus confirms that the training of the model went successfully, rendering balanced results for the prediction of instability through the analysis of stable state samples of signal. The model adds the possibility to determine the directions concerned by the instability on a bi-parameter plane. User imposed precision and safety tolerances can then be used to enhance the definition of thresholds for a better detection of the instability.

Using the trained CRNN to predict instability

The best trained model obtained in the previous sections is now frozen to predict labels on operating points never seen before. The following paragraphs introduce a possible application for this selected network.

Real-time prediction

The network is trained to attribute a label on a 0.3 s time sample at an acquisition frequency of 10 kHz. In real use, the experiment can be lit for minutes and sensors deliver information continuously. For real-time prediction, the signal from the sensors must be seen by the network regularly under the exact form of a 0.3 s sample acquired at 10 kHz. In order to approach this configuration, the trained network is successively fed by the content of a rolling window passing on a 10.3 s long signal as illustrated in Fig. 12. At each step of the rolling window, the model produces a label predicting the proximity of unstable operating points, based on its training. This operation is finally averaged over the total acquisition time of 10.3 s, yielding an overall prediction compared to the original label determined by the experiment.

An important parameter of the rolling window process is the delay imposed between each translation step. The minimal delay is theoretically imposed by the acquisition frequency and thus reaches 0.1 ms. However it can be changed in order to fit the practical limitations of the acquisition system. In fact, its first order limitation is the inference delay, which is the time needed for the network and the architecture translating its results to produce one label from the analysis of one sample. The longer the inference delay is, the more spaced the predictions are, impeding a real-time use. On an Intel i5-8210Y CPU dedicated to personal computers, the maximum inference time is 19 ms. Thus a minimum of 19 ms is to be imposed as the refreshing time in-between each prediction. This delay is restrained enough to guarantee a fast response of the analyzing system to sudden changes of the operating parameters of the experimental setup, offering the possibility to use the model as a real-time diagnostic when directly wired to the stream of information obtained by the sensors. It should be noted that modern dedicated architectures for network inference, i.e. in embedded systems, can be orders of magnitude faster, thus significantly lowering the limitations in real systems.

Real use scenario

In order to simulate a real-use case, the rolling-window is set with a 20 ms delay between two successive evaluations. Figure 13 shows the rolling window evaluation of the model on the operating point PC1, U sw = 30 m.s -1 , φ = 0.70, taken from the testing dataset. The four outputs of the network, with respect to the four directions assessed, are plotted in four different graphs. The point considered in Fig. 13 is surrounded by two unstable operating points, the first one in the upper direction (U sw = 30 m.s -1 , φ = 0.75), and the second one on the right (U sw = 32 m.s -1 , φ = 0.70). Therefore, the network should predict constant values of 1 for the "UP" and "RIGHT" directions and constant values of 0 for the "DOWN" and "LEFT" directions. In practice, for this operating point, the network is able to give the right prediction nearly constantly. However, due to the intermittent nature of the signal recorded by the sensors for operating points close to instability, it appears that the network is mislead during short periods concerning the higher bulk velocity direction, i.e. "RIGHT" in Fig. 13, leading to the value 0.9. Yet, the tool is still steady enough for the interpretation of the directions pointing toward unstable operating points to remain clear, confirming that fine specific features of the signal have been extracted from the noise and are translated by the CRNN in order to inform the user. The possibility to have an updated prediction every 20 ms gives the opportunity to have a real-time directional information on the proximity of unstable operating points.

Influence of essential choices made during implementation

During the implementation of the final approach, steps have been taken regarding two essential elements of the model, the structure of the input and the architecture of the network. In order to observe the effects of these developments, the accuracy is given along with the AUC from the ROC analysis. A third metric is added in order to measure the confidence of the predictions, the Binary Cross-Entropy (BCE) [START_REF] Mitchell | Artificial Neural Networks[END_REF][START_REF] Goodfellow | Deep Feedforward Networks[END_REF], formerly used with regularization as a loss function during the training phase. Here, only the Cross-Entropy term is kept, such that for each of the four directions:

BCE = - 1 N s N s i=1 y i log(ŷ i ) + (1 -y i ) log(1 -ŷi ) (6)
where N s is the total number of sample composing the evaluated set (testing or validation set), ŷi is the predicted probability that the considered direction of the i th sample is unstable, y i is the true label, and y i = 1 if the considered direction of the i th sample is unstable, y i = 0 if it is not. The four BCE obtained for each direction are then averaged to render the overall Binary Cross-Entropy of the model.

Choice of the set of sensors

A key step when building the tool is to choose the content of the samples placed at the input of the network. The number of channels is essential to define, provided that the more sensors are available the better it usually is in order to capture the dynamic of a system, and that the invasiveness of diagnostics is currently a crucial limitation for real applications. To briefly study the influence of the choice of sensors on the performance of the predictor, three configurations have been tested. In the first configuration, the only sensor studied is the microphone M1. The prediction is made on single-channel samples. In the second configuration both the photo-multiplier PM and the microphone M1 are considered. The prediction is made on two-channel samples. Finally, the third configuration is the three-channels sample configuration, proposed in previous sections. Table 4 shows the difference in precision appearing when the number of channels is changed.

Both the accuracy and the AUC worsen with nearly constant BCE when the number of channels is reduced. However it is interesting to point out that using only the pressure sensor is already enough to reach 87.6 % of validation accuracy and 84.8 % of testing accuracy with acceptable AUC. Furthermore, when coupled with a photo-multiplier, the prediction is improved to reach a precision close to the three-channels configuration. The three-channels configuration however remains the most efficient one. Nonetheless, it appears that the evaluation of the acoustic pressure alone can be sufficient to obtain a fair prognosis of the stability of the combustor.

This succinct comparison is a first step of an ongoing study of the influence of the choice of sensors on the efficiency of the prediction. The next step is to test every combination of sensors and evaluate its influence on the quality of the result, along with the minimization of the invasiveness of the diagnostics involved.

Evolution of the architecture of the model

The CRNN is built by stacking layers of different types. In order to observe the benefits of the addition of each type of layer, three networks are evaluated. The first network is a three layer MLP built with the same layers as the three last Fully Connected layers of the final CRNN (See Tab. 2). The second network is a CNN built with four convolution layers stacked over the MLP. Finally, the third network is the final CRNN architecture. Each network is representative of a step of the evolution of the model through the development process. Evaluated on the same dataset, the Binary Cross-Entropy (BCE) loss and the accuracy for both validation and testing sets along with the AUCs are computed and given in Tab. 5.

The addition of convolution layers on top of the MLP enhances its accuracy but has a limited effect on the confidence of the network, as shown by BCE scores. However, the combination of convolution layers with a recurrent LSTM layer both increases the accuracy and the confidence of the model. Therefore, for this specific task, it appears that an efficient architecture is obtained with a CRNN approach.

The path to the final architecture is therefore marked by two main observations. First, the accuracy of a Multi-Layer Perceptron (MLP) when performing this multi-class-multi-label classification task is low. Enhancing the MLP with four convolution layers and thus building a Convolutional Neural Network (CNN) normally helps to aim for patterns in the signal and, in this case, increases the accuracy and AUC of the model. Then, the CNN displays high validation and testing accuracy when performing the task but tends to do over-confident wrong predictions. Adding the LSTM recurrent layer helps to mitigate this flaw and thus both increases the accuracy and decreases the BCE, leading to better results.

This study does not offer a comprehensive search for the optimal network for the task at hand, and numerous parameters could be explored for each network. The focus here is on the proof of concept. However, further work, notably including more data, could focus on optimizing the neural network, e.g. using Neural Architecture Search [START_REF] Elsken | Neural architecture search: A survey[END_REF].

Critical analysis

Interpretability of a neural network

The promising results obtained in this study suggest that meaningful features are present in the signal at the onset of the instability and that these features seem to translate a directional information on a bi-parameter operating plane. However, the use of Deep Learning raises one crucial criticism. Even if the model efficiently associates features to empirical scores, it does not provide the exact characteristics of these features. This is due to the fact that the model is built by an optimization algorithm which is practically impossible to interpret. This can be a problem when the reliability of the tool has to be evaluated, and could limit the applicability of such tools for real-world combustors, however well they perform during training and testing. Progress in these techniques must be accompanied by the development of tools to interpret them in the future.

Supervised Machine Learning has been described as the ability to produce internal representations from the analysis of inputs and their corresponding outputs [START_REF] Goodfellow | Machine Learning Basics[END_REF]. Interpreting the model is equivalent to finding a way to explain these representations. Several groups have worked on the interpretability of neural networks through the study of representations [START_REF] Gilpin | Explaining explanations: An overview of interpretability of machine learning[END_REF][START_REF] Buhrmester | Analysis of Explainers of Black Box Deep Neural Networks for Computer Vision[END_REF]. Yet, these remain case-by-case approaches adapted to very specific tasks. In this study, there is no simple intuition to explain what allows convergence based on combustion theory. In fact, it is also precisely one of the supporting reasons for the use of Deep Learning. Its power is revealed when input sizes are large and when the features searched in it are unknown a priori by the user.

The context of data-starvation

Using Deep Learning techniques in an experimental context where the dimensions of the dataset can be strongly limited underlines a recurrent problem of such an approach. The learning algorithm demands a certain quantity of training items which is difficult to evaluate beforehand as it strongly depends on the task. In this study, the apparently low number of stable operating points available is mitigated by the possibility to sample each operating point in overlapping windows. Still, it does not protect from having created hidden trivial structures not seen at first glance but found by the model during training.

However, the construction of a simpler baseline helps to deal with dataset bias. On the one hand, it gives an easy-to-setup approach to solve the problem with performances comparable to more complete approach. On the other hand, it helps to see the hidden structures of the dataset. In fact, if the baseline displays results close to or better than the model, it is an obvious sign that the hidden parameter is revealed. For example in this case, if the baseline (focused on the bulk velocity and the equivalence ratio) performed as well as the model, it would emphasize the fact that these operating parameters are detectable in the signal and sufficient to conclude on the proximity of instability. As it is not clearly the case in this study, the intuition that discriminant patterns compose the signal and could be picked up by the CRNN to give a directional information on the stability of the burner still holds.

Nonetheless, increasing the size of the dataset and finding more robust definitions of the proximity of instability are further objectives to reach, this work being a very first step.

Conclusion

The use of a CRNN for the analysis of combustion dynamics shows promising results to change operating conditions without entering an unstable domain. The trained network predicts the proximity of unstable operating points on the bulk velocity-equivalence ratio plane for different configurations by treating short windows of 0.3 s of acoustic pressure, acoustic flux and heat release rate fluctuations recorded during stable states of the system.

The steps to build the learning tool are:

Results show that it is possible to locate unstable regions in the operating map based on records during stable operation and signal treatment from a limited number of sensors, thus avoiding violent instability episodes when operating in unknown conditions. This analysis can be done in real-time, giving a constant help as signal streams from the sensors. Additionally, the CRNN approach appears as an efficient first step towards a data centered tool able to systematically warn the user about impending critical issues.

Several directions must be investigated to improve this initial study. First, it would be necessary to improve the generalization ability of the tool to a wider range of configurations with finer operating grids and output labels containing more information. Evaluating the robustness to transformations of the burner and operating atmospheric conditions is necessary to ensure that the method can learn a general feature about proximity to instability, regardless of the configuration. Therefore, a complete analysis of the resilience of the model to changes in atmospheric conditions and injection configurations is needed. Furthermore, the modularity of the MIRADAS burner allows the transformation of its geometry (swirler, chamber and plenum), a key factor toward an efficient and robust learning system for the prediction of thermo-acoustic instability in different combustion devices. The construction of new metrics to predict, applicable to different combustors and thus less linked to empirical considerations is also required for generalization. In parallel, the influence of changes in the structure of the network and the choice of sensors to analyze should be studied, together with a complete hyperparameter tuning procedure, so that the optimum performances of such an approach can be achieved. Finally, a better interpretation of what the network can unveil from combustion dynamics is needed. It is the key step to both improve the theoretical understanding of precursors of combustion instability and design reliable, and efficient tools for their diagnostic in real-size combustion devices.

2 .

 2 validation set: It is built to evaluate the error during the training on data that are not used for the training and therefore not seen by the model. It allows to monitor how well the model generalizes on unseen data through the epochs of training. Hyperparameters of the training process -how many epochs to train, parameters of the optimization algorithm, details of the network architecture -can be chosen optimally with respect to this generalization task, instead of optimizing only for the training loss and risking overfitting. 3. testing set: The samples of this set are kept completely unseen during the training. They are used to verify that the label attribution learnt on the training set is rightfully performed on new data.

Figure 1 :

 1 Photo-Multiplier

Figure 2 :

 2 Figure 2: Schematic representation of the injection configurations [18].

Figure 6 :Figure 7 :

 67 Figure 6: Schematic representation of the CRNN used in this study [24].

Figure 11 :

 11 Figure 11: Choice of constraints and corresponding performances in terms of False Positive Rate (FPR) and True Positive Rate (TPR) for the four evaluated directions.

Figure 12 :

 12 Figure 12: Representation of the rolling window evaluation made to simulate the prediction on a streaming signal.

Table 1 :

 1 Description of the injection configurations.

			Main		Secondary		N op
			injection		injection
	Ref					X				56
	PC1 pilot CH Total					366
	Main	Main	Main	Pilot	Main	Main	Pilot	Main	Main	Main
	CH4 + Air	CH4 + Air	CH4 + Air	CH4 + H2 + Air
	1. Ref		2. PC1 3. PC2	4. PH1 5. PH2	6. MH1 7. MH2

4 (1 % in Pow.) 55 PC2 premixed pilot CH 4 (2 % in Pow.) 52 PH1 CH 4 /air pilot H 2 (1 % in Pow.) 50 PH2 pilot H 2 (2 % in Pow.) 51 MH1 prem. H 2 (1 % in Pow.) 50 MH2 prem. H 2 (2 % in Pow.) 52

Table 4 :

 4 Influence of the number of channels on the performance of the method.

			1 Ch.	2 Ch.	3 Ch.
	Validation Accuracy 87.6 % 88.4 % 90.6 %
		BCE	0.33	0.28	0.27
	Testing	Accuracy 84.8 % 87.5 % 87.7 %
		BCE	0.35	0.30	0.34
	ROC-AUC UP	0.89	0.91	0.92
		DOWN	0.98	0.98	0.99
		LEFT	0.78	0.84	0.88
		RIGHT	0.92	0.96	0.95

Acknowledgements

This project has received funding from the European Research Council under the European Union's Horizon 2020 research and innovation program Grant Agreement 832248, SCIROCCO.

Source code and data for the deep learning in this study is available at https://gitlab.com/cerfacs/code


François Boyer

Eric Falcon

Wave turbulence on the surface of a ferrofluid submitted to a magnetic field

Keywords: numbers: 47.35.Tv, 47.65.Cb, 47.27.-i

We report the observation of wave turbulence on the surface of a ferrofluid mechanically forced and submitted to a static normal magnetic field. We show that magnetic surface waves arise only above a critical field. The power spectrum of their amplitudes displays a frequency-power law leading to the observation of a magnetic wave turbulence regime which is experimentally shown to involve a 4-wave interaction process. The existence of the regimes of gravity, magnetic and capillary wave turbulence is reported in the phase space parameters as well as a triple point of coexistence of these three regimes. Most of these features are understood using dimensional analysis or the dispersion relation of the ferrohydrodynamics surface waves.

Wave turbulence is an out-of-equilibrium state where waves interact with each other nonlinearly through Nwave resonance process. The archetype of wave turbulence is the random state of ocean surface waves, but it appears in various systems: capillary waves [1,2], plasma waves in solar winds, atmospheric waves, optical waves, and elastic waves on thin plates [START_REF] Zakharov | Kolmogorov Spectra of Turbulence[END_REF]. Recent laboratory experiments of wave turbulence have shown new observations such as intermittency [START_REF] Falcon | [END_REF], fluctuations of the energy flux [5], and finite size effect of the system [2,6]. Some of these phenomenon have recently been considered theoretically [7]. Wave turbulence theory allows to analytically derive stationary solutions for the wave energy spectrum as a power law of frequency or wave number [START_REF] Zakharov | Kolmogorov Spectra of Turbulence[END_REF]. The spectrum exponent and the number N of resonant waves depend on both the wave dispersion relation and the dominant nonlinear interaction. Several theoretical questions are open, notably about the validity domain of the theory [8], and the possible existence of solutions for non dispersive systems [9]. In this context, finding an experimental system where the dispersion relation of the waves could be tuned by the operator should be of primary interest to test the wave turbulence theory.

A ferrofluid is a suspension of nanometric ferromagnetic particles diluted in a liquid displaying striking properties: the Rosensweig instability [10], the labyrinthine instability, magnetic levitation [START_REF] Rosensweig | Ferrohydrodynamics[END_REF]. In contrast with usual liquids, the dispersion relation of surface waves on a ferrofluid displays a minimum which depends on the amplitude of the applied magnetic field [START_REF] Broaweys | [END_REF]13]. Thus, one can easily tune the dispersion relation of surface wave from a dispersive to a non dispersive one with just one single control parameter. To our knowledge, no experimental observation of wave turbulence on a magnetic fluid has been reported. Here, we study the wave turbulence on the surface of a ferrofluid submitted to a normal magnetic field. We observe for the first time a regime of magnetic wave turbulence. We characterize this regime by measuring the power spectrum and distribution of the The experimental setup is shown in Fig. 1. It consists of a cylindrical container, 12 cm in inner diameter and 4 cm in depth, filled with a ferrofluid up to a depth h = 2 cm. The ferrofluid used is a ionic aqueous suspension synthesized with 8.5% by volume of maghemite particles (Fe 2 0 3 ; 7.6 ±0.36 nm in diameter) [START_REF]The ferrofluid synthesis has been performed by the Laboratory LI2C[END_REF]. The properties of this magnetic fluid are: density, ρ = 1324 kg/m 3 , surface tension, γ = 59 × 10 -3 N/m, initial magnetic susceptibility, χ i = 0.69, magnetic saturation M sat = 16.9 × 10 3 A/m, and estimated dynamic viscosity 1.2 ×10 -3 Ns/m 2 . The container is placed between two horizontal coaxial coils, 25 cm (resp. 50 cm) in inner (resp. in outer) diameter, 7 cm far apart. A DC current is supplied to the coils in series by a power supply (50 V/35 A), the coils being cooled with water circulation. The vertical magnetic induction generated is 99% homogeneous in the horizontal plane [START_REF] Broaweys | [END_REF], and is up to 780 G. It is measured by a Hall probe located in the center near the surface of the container. Surface waves are generated on the ferrofluid by the horizontal motion of a rectangular plunging Teflon wave maker driven by an electromagnetic vibration exciter. The wave maker is driven with low-frequency random vibrations (typically from 1 to 5 Hz). The amplitude of the surface waves η(t) at a given location is measured by a capacitive wire gauge (plunging perpendicularly to the fluid at rest) with a 7.1 mm/V sensitivity [2]. η(t) is recorded by means of an acquisition card with a 4 kHz Typeset by REVT E X sampling rate, low-pass filtered at 1 KHz during 300 s, leading to 1.2 × 10 6 points recorded.

In the deep fluid approximation, the dispersion relation of linear inviscid surface waves on a magnetic fluid submitted to a magnetic induction B perpendicular to its surface, reads [START_REF] Rosensweig | Ferrohydrodynamics[END_REF]

ω 2 = gk - f [χ] ρµ 0 B 2 k 2 + γ ρ k 3 (1)
where ω is the wave pulsation, k its wavenumber, g = 9.81 m/s 2 the acceleration of the gravity, µ 0 = 4π × 10 -7 H/m the magnetic permeability of the vacuum, and

f [χ] ≡ χ 2 /[(2+χ)(1+χ)]
. χ is the magnetic susceptibility of the ferrofluid which depends on the applied magnetic field, H, through Langevin's classical theory [START_REF] Abou | [END_REF] χ

(H) = M sat H L 3χ i H M sat (2)
where L(x) ≡ coth (x) -1/x, and thus on the magnetic induction, B, through an implicit equation since

B = µ 0 (1 + χ)H (3)
For B = 0, the dispersion relation of Eq. ( 1) is monotonous, and is dominated by the gravity waves at small k, and by the capillary waves at large k. When B is increased, the quadratic term -B 2 k 2 increases, and the dispersion relation becomes non monotonous: an inflection point appears at B = 0.93B c , then a minimum which leads, at B = B c , to the Rosensweig instability [START_REF] Rosensweig | Ferrohydrodynamics[END_REF]. This stationary instability (a hexagonal pattern of peaks on the ferrofluid surface) occurs when ω 2 (k) becomes negative, that is for a critical induction

B 2 c = 2µ 0 √ ρgγ/f [χ(Hc)]
where χ(H c ) is determined using Eqs. ( 2) and (3) [10]. This leads to the theoretical value of the critical induction B c = 292.3 G for the Rosensweig instability of our ferrofluid. When B is slowly increased, a direct visualization of the surface leads to an experimental value of B c = 294 ± 2 G which is close to the above expected value.

The power spectrum of the wave amplitude on the surface of the ferrofluid is shown in Fig. 2 for different applied magnetic induction B. For B = 0 (see the inset of Fig. 2), it displays similar results than those found with a usual fluid [2]: two power laws corresponding to the gravity and capillary turbulence wave regimes. The capillary regime is found to scale as f -2.9±0.1 in good agreement with the prediction of weak turbulence theory in f -17/6 [16], and the gravity regime is found in f -4.6 . The exponent of the gravity cascade being known to depend on the forcing parameter [2,6] (in contrast with the theoretical prediction ∼ f -4 [17]), it is thus only fitted to measure the crossover frequency f gc between gravity and capillary regimes. As previously reported with usual fluids [2], f gc is also found here to decrease (from 26.6, 21 to 17.2 ± 0.3 Hz for a random forcing of 1 -6 Hz, 1 -5 Hz, and 1 -4 Hz). The expected value is given by f gc = 1 π 2g/l c ≃ 15.2 Hz where l c = γ/(ρg) [2]. Consequently, in order to study the crossover frequency dependence with the magnetic induction B, one has to rescale it by its value at B = 0 such as f (B) ≡ f gc (0) × f (B)/f (0) (see below).

For B = 0, the power spectrum of wave amplitudes shows two striking new results. As shown in Fig. 2, the crossover frequency is strongly decreased down to roughly 5 Hz, and a power law in f -3.3 appears in roughly all the accessible frequency range. Let us try first to understand this latter observation. The power spectrum of η(t) can be derived by dimensional analysis for the gravity and capillary wave turbulence regimes, respectively as [9]

S grav η (f ) ∼ ǫ 1/3 gf -4 and S cap η (f ) ∼ ǫ 1/2 γ ρ 1/6 f -17/6 ,
where ǫ is the energy flux [dimension (L/T ) 3 ]. Since gravity wave turbulence is a 4-wave interaction process, and capillary waves a 3-wave one, it gives the dependence of the energy flux exponent in ǫ 1/(N -1) for a N -wave process [9]. One can also derive dimensionally the power spectrum for the magnetic wave turbulence regime as

S mag η (f ) ∼ ǫ α B 2 ρµ 0 2-3α 2 f -3 , (4)
In contrast with the above dispersive systems, this -3 frequency exponent does not depend on the energy flux exponent α, that is on the number N of resonant waves. The frequency exponent predictions for the magnetic and capillary regimes, respectively -3 and -17/6 ≃ -2.8 can not be distinguished experimentally within our experimental accuracy. This could explain that only one single slope is observed on Fig. 2 which thus corresponds to a "magneto-capillary" wave turbulence regime. Figure 3 shows the frequency exponent of the magnetocapillary spectrum when B is increased. For small B, capillary waves are dominant and the exponent is found roughly constant ≃ -2.9 in good agreement with the capillary prediction. When B/B c ≥ 0.65, magnetic waves becomes dominant (see below), and the exponent begins to slightly increase with B up to -3.1 in rough agreement with the -3 prediction of Eq. ( 4). For B/B c ≥ 1, the Rosensweig instability occurs, and the spectrum exponent strongly changes with B. This could be attributed to the growing of the hexagonal pattern. The inset of Fig. 3 shows the power spectrum amplitude (averaged between 8 and 18 Hz) as a function of B. The spectrum amplitude is found roughly constant when B is increased till B/B c ≃ 0.65, onset of the magnetic waves. Above this value, the spectrum amplitude increases roughly linearly with B up to B/B c = 1 where the instability occurs. Using Eq. ( 4), one can thus deduce from this linear dependence that 2 -3α = 1, and thus α = 1/3. This experimentally shows that the magnetic wave turbulence observed involves a 4-wave interaction process.

The crossover frequency between the gravity and magneto-capillary regimes decreases when B is increased (see Figure 2). Figure 4 then shows the evolution of the rescaled crossover frequency f (B) (see above) as a function of B for 3 different frequency bandwidths of the random forcing. When B is increased, f is found to decrease with the same law whatever the forcing frequency. Beyond B = B c , f is roughly of the same order than the upper frequency of the forcing, and consequently can not be measured anymore. However, when B > B c and for the 1 -6 Hz forcing, the power spectrum displays three power laws (not shown here) showing a magneto-capillary crossover. This new crossover frequency is reported in Fig. 4 with (•)-symbols without rescaling. For lower fre- 5), ( 6) and (7). The triple point (ft = 10.8 Hz , Bt/Bc = 0.65) is from Eq. ( 8).

quency bandwidths of vibration, the slope breaking is too small to allow an accurate measurement of the transition between the magnetic and capillary regime.

The evolutions of these crossover frequencies with B are described as follows. Whatever B, Eq. ( 1) is dominated, at small k, by the linear term (gravity waves) and, at high k, by the cubic term (capillary waves). One can assume that the quadratic term (magnetic waves) dominates when it is greater than the linear and cubic terms. This arises when f [χ]B 2 > µ 0 √ ρgγ, that is, using Eqs. ( 2) and (3) and the ferrofluid properties, when B > 0.65B c . Thus, when B < 0.65B c , no wavelength exists for which the magnetic term dominates in Eq. (1). When B > 0.65B c , magnetic waves exist on a range of wavelength between the gravity and capillary ones. This explains why both the spectrum exponent and its amplitude shown in Fig. 3 change for B ≃ 0.65B c . Surprisingly, this critical magnetic induction has never been reported previously. The crossover frequencies between the gravity, the magnetic and the capillary regimes are derived by balancing the dispersion relation terms each to each. For the gravity-capillary transition, one balances the first and the third terms of the second-hand of Eq. ( 1), that is gk gc = (γ/ρ)k 3 gc , thus for k gc = ρg/γ which substituted into Eq. ( 1) gives

ω 2 gc = 2 g 3 ρ γ - gf [χ]B 2 µ 0 γ , for f [χ]B 2 < µ 0 √ ρgγ . (5)
Similarly, by balancing the first and second terms, the gravity-magnetic crossover frequency reads

ω 2 gm = γ ρ µ 0 ρg f [χ] 3 B -6 , for f [χ]B 2 > µ 0 √ ρgγ . (6)
Finally, by balancing the second and the third terms, the

ω 2 mc = gf [χ] µ 0 γ B 2 , for f [χ]B 2 > µ 0 √ ρgγ . (7)
Using Eqs. ( 2) and ( 3) and the ferrofluid properties, these crossover frequency curves are plotted in Fig. 4 as a function of B/B c and show a good agreement with the experimental data. This plot also shows the range of existence of magnetic waves for B/B c > 0.65 and of a triple point corresponding to the coexistence of the three domains (gravity, magnetic and capillary ones). It can be derived by balancing the three terms of Eq. ( 1), which leads to

f t = 1 2π g 3 ρ γ 1/4 and B 2 t = µ 0 √ ρgγ/f [χ(H t )] , (8)
corresponding to f t = 10.8 Hz and B t /B c = 0.65 in good agreement with the data of Fig. 4. Finally, the probability density functions (PDFs) of the wave amplitudes are shown in Fig. 5 for different value of B/B c , at high enough amplitude of forcing. For B = 0, the PDF is asymmetrical due to the strong steepness of the waves as for usual fluids [2]. This means that the deep through are rare, whereas high crests are much more probable, thus showing the nonlinear nature of the wave interactions in a wave turbulence regime. The PDF asymmetry is enhanced when B is increased. Note also that the most probable value of the PDF is more and more negative as B increases, although its mean value η is zero as expected. The bottom inset of Fig. 5 shows that all these PDFs normalized to its standard deviation value, σ η , roughly collapse on a single non Gaussian distribution. This means that the distribution depends only on σ η . The top inset of Fig. 5 then shows the evolution of σ η as a function of B/B c . For B/B c ≤ 0.63, the rms value amplitude of the waves does not depend on the magnetic induction. Note that this value is very close to the above predicted onset of magnetic waves B t /B c = 0.65. When this onset is exceed, σ η increases roughly linearly with B up to the occurrence of the Rosensweig instability at B c .

FIG. 1 :

 1 FIG. 1: Experimental setup

FIG. 2 :

 2 FIG. 2: Power spectrum of η(t) for two values of B. Inset: B = 0: Gravity and capillary wave turbulence regimes. Dashed lines have slopes -4.6 and -2.9. Crossover: fgc ≃ 20 Hz. Main: B = 0.9Bc: Magneto-capillary wave turbulence. Dashed line has slope -3.3. Crossover: fgc ≃ 5 Hz. Forcing parameters: 1 -5 Hz. Bc = 294 G.

FIG. 3 :

 3 FIG.3: Exponent of the magneto-capillary spectrum as a function of the dimensionless magnetic induction B/Bc. Inset: Amplitude of the power spectrum (averaged between 8 and 18 Hz) as a function of B/Bc. Forcing parameter: 1 ≤ f ≤ 4 Hz. Bc = 294 G.

FIG. 4 :

 4 FIG. 4: Rescaled crossover frequencies, f (B) ≡ fgc(0) × f (B)/f (0) as a function of B for different frequency bandwidths of the random forcing: (×) 1 to 4 Hz, (⋄) 1 to 5 Hz, and (+ or •) 1 to 6 Hz. Theoretical curves fgc, fgm and fmc are respectively from Eqs. (5), (6) and(7). The triple point (ft = 10.8 Hz , Bt/Bc = 0.65) is from Eq. (8).

FIG. 5 :

 5 FIG.5: Probability density functions of wave amplitude, η, for different values of the dimensionless magnetic induction, from B/Bc = 0, 0.3, 0.54, 0.77, 0.92 to 1.2 (see the arrows). Forcing parameter: 1 ≤ f ≤ 4 Hz. Bottom inset: Same PDFs displayed using the reduced variable η/ση. Gaussian fit (dashed line). Top inset: Standard deviation of wave amplitude ση as a function of B/Bc

We thank D. Talbot for the synthesis of the ferrofluid, J.-C. Bacri and A. Cebers for fruitfull discussion, A. Lantheaume, and C. Laroche for technical assistance. This work has been supported by ANR Turbonde BLAN07- 3-197846.