Antonios Giannopoulos

J.-L Aider

Data-driven order reduction and velocity field reconstruction using Neural Networks: the case of a turbulent boundary layer

We present a data-driven methodology to achieve identification of coherent structures dynamics and system order reduction of an experimental turbulent boundary layer (TBL) flow. The flow is characterized using time-resolved Optical Flow Particle Image Velocimetry, leading to dense velocity fields that can be used both to monitor the overall dynamics of the flow and to define as many local visual sensors as needed. A Proper Orthogonal Decomposition (POD) is first applied to define a reduced-order system. A non-linear mapping between the local upstream sensors (inputssensors) and the full-field dynamics (POD coefficients) as outputs is sought using an optimal Focused Time-Delay (FTD) Artificial Neural Network (ANN). The choices of sensors, ANN architecture and training parameters are shown to play a critical role. It is verified that a shallow ANN, with the proper sensor memory size, can lead to a satisfying full-field dynamics identification, coherent structure reconstruction, and system order reduction of this turbulent flow.

I. INTRODUCTION

Current experimental methods such as Particle Image Velocimetry (PIV) or numerical methods such as Large Eddy Simulations (LES), used to characterize 3D non-stationary flows, produce a large amount of data leading to high-order systems. It is generally necessary to use an order reduction to derive some information that could be handled either to model or control the targeted flow. Data-driven methods are nowadays becoming more and more efficient and reliable even for fluid mechanics research [START_REF] Deng | Time-resolved turbulent velocity field reconstruction using a long short-term memory (lstm)based artificial intelligence framework[END_REF][START_REF] Guzmán-Iñigo | Data-based, reducedorder, dynamic estimator for reconstruction of nonlinear flows exhibiting limit-cycle oscillations[END_REF] . Among successful applications, one can cite statistical learning [START_REF] Inj | A dynamic observer to capture and control perturbation energy in noise amplifiers[END_REF] or machine learning [START_REF] Gautier | Closedloop separation control using machine learning[END_REF][START_REF] Wang | Model identification of reduced order fluid dynamics systems using deep learning[END_REF] algorithms. More recently Artificial Neural Networks (ANNs) have been used extensively for fluid mechanics, because of their architecture flexibility, needed for reduced-order modelling 29 , as well as for field reconstruction of complex non-linear flows [START_REF] Liu | Deep learning methods for super-resolution reconstruction of turbulent flows[END_REF] . A short review on different applications of Deep Learning to Fluid Dynamics can be found in [START_REF] Kutz | Deep learning in fluid dynamics[END_REF] . Nevertheless, applying identification or modelling algorithms to experimental fluid mechanics is still challenging because of various noises leading to non physical errors as well as a limited number of sensors and the complexity of the experimental setup. Large datasets with both good spatial and temporal resolutions are needed for successful applications of these data-hungry methods to experiments.

In an inverted flag flow experiment, Deng et al. [START_REF] Deng | Time-resolved turbulent velocity field reconstruction using a long short-term memory (lstm)based artificial intelligence framework[END_REF] applied an ANN identification method to reconstruct time-resolved velocity fields from a handful of velocity sensors. It was also proven recently that ANNs can be used to predict the dynamics and reconstruct the time-resolved fields of an experimental Backward-Facing Step (BFS) flow [START_REF] Giannopoulos | Prediction of the dynamics of a backward-facing step flow using focused time-delay neural networks and particle image velocimetry data-sets[END_REF] . We apply in this study a similar methodology to an experimental TBL flow, with the objective to identify correctly the dynamics, reconstruct the a) Electronic mail: antonios.giannopoulos@espci.fr b) Electronic mail: jean-luc.aider@espci.psl.eu FIG. 1. Control scheme for minimizing perturbation kinetic energy using a neural-network identifier velocity field and predict the total Perturbation Kinetic Energy (PKE) of the full velocity field. The ANN algorithm will be used first to identify the global dynamics from upstream local sensors and second to reduce drastically the order of the system to a handful of easily measurable quantities. Such an algorithm would be an ideal part of a model-free control loops [START_REF] Gautier | Closedloop separation control using machine learning[END_REF] aiming at reducing the PKE field, as shown in the Fig. 1.

A. Artificial Neural Networks

Artificial Neural Networks (ANNs) can provide a nonlinear mapping between one set of inputs-signals (that are easily measurable) and some corresponding output signals (that may be hard to measure in a final application, but still valuable for monitoring or control purposes). Great progresses have been made lately due to the ability to handle very large datasets, the advancements of Graphics Processing Units (GPU) parallel programming, and the increasing availability of multiple optimized toolboxes. A short review of applications of deep ANN to fluid mechanics can be found in Nathan Kutz [START_REF] Kutz | Deep learning in fluid dynamics[END_REF] . ANNs can be divided into Feed-Back (or recurrent) and Feed-Forward. They can also be discriminated according to their depth, either shallow or deep, depending on the number of hidden layers (one or more). Finally, they can be divided into static or dynamic, if the output of the current step depends on the previous steps as well, giving it a notion of memory.

In the case of a feed-forward ANNs, the output of any layer only modifies the next layer, not the previous ones. On the other hand, in a recurrent ANN the system output is calculated from its previous past time-steps along with the input at the current time-step (eq. 1), hence introducing a notion of output memory in the network.

y t = f w,b (y t-1 , x t ) (1)
One category of recurrent ANNs are the Long-Short Term Memory (LSTM) networks. They are gradient-based recurrent NNs used for supervised learning both in classification and prediction problems [START_REF] Vlachas | Datadriven forecasting of high-dimensional chaotic systems with long-short term memory networks[END_REF] . They were first developed by Hochreiter and Schmidhuber 18 to solve the vanishing or exploding gradient problem of the back-propagated error. In an LSTM architecture, the network is left to learn alone the size of the memory of each neuron during the training process. They do so using a sophisticated gate-neuron that determines if the input is important enough or if it should be forgotten and when it should output its value.

B. Focused time-delay Artificial Neural Networks

The neural network used in the present System Identification (SI) study is a fully-connected FTD ANN, which was first introduced for speech recognition [START_REF] Waibel | Phoneme recognition using time-delay neural networks[END_REF] . It is basically a standard feed-forward architecture along with a tapped constant timedelay (of time-step size k) in the input (see Fig. 2). The term "focused" comes from the fact that the notion of memory is introduced only in the input (sensors), and not in the output like in LSTM ANNs. They are used to model long-range temporal dependencies by keeping a number of k past measurements of the input at each time step x t . This leads to the following expression for the output of the system :

y t = f w,b (x t , x t-1 , ..., x t-k ) (2)
where w and b are weight and bias parameters. They have been used for rainfall prediction [START_REF] Charaniya | Focused time delay neural network model for rainfall prediction using indian ocean dipole index[END_REF] as well as for effective acoustic modelling [START_REF] Huang | Deeper time delay neural networks for effective acoustic modelling[END_REF] , in deeper architectures.

For the above network with a non-linear activation function f 1 in the hidden layer and a linear activation function f 2 in the linear layer, the equation giving the k th neuron output of a single hidden layer network connected to the j th neuron of the previous layer is: where N 1 is the number of neurons in the first (hidden) nonlinear layer and N 2 the number of neurons in the second (linear) layer. Regarding the choice of the non-linear activation function f 1 , the tan-sigmoid or hyperbolic tangent function is used:

y k = f 1 N 2 ∑ j=0 w (2) k j f 2 ( N 1 ∑ i=0 (w (1) k j x i + b i )) + b j (3)
tanh(x) = e 2x -1 e 2x + 1 (4)
The correct training process of the ANN (meaning finding the optimum weights and biases connecting the neurons of different layers) consists in dividing the data-set between a training data-set and a validation data-set. For the training data-set we chose one set of weights (randomly initialized in the first iteration) connecting the layers and we use the same weights for each time step to calculate the output of the model. The error of the real vs the model-generated output signal is computed and the weights and biases are updated according to different back-propagation schemes (in our case the Scaled Conjugate Gradient method). The above process is called one "epoch". We continue the process for as many epochs as needed until a satisfactory fit error is achieved. The second data-set is used to test the performance of the network on new data and verify the achieved error, hence avoiding overfitting.

In the following, we first present the experimental setup as well the main characteristics of our Turbulent Boundary Layer (TBL) flow. The choice of the sensors of the ANN is then discussed and parametric studies of different training parameters are presented. The efficiency of the ANN for the dynamics identification and velocity field reconstruction is illustrated before turning to the discussion and conclusion.

II. EXPERIMENTAL SETUP A. Hydrodynamic channel

Experiments were carried out in an hydrodynamic channel in which the flow is driven by gravity. The flow is stabilized by divergent and convergent sections separated by honeycombs leading to a turbulence intensity of 0.6 %. A NACA 0020 profile is used to smoothly start the boundary layer. The test section is 80 cm long with a rectangular cross-section 15 cm wide and 8.5 cm high as shown in Fig. 3.

B. Time-resolved Particle Image Velocimetry measurements

The flow is seeded with 20 µm neutrally buoyant polyamide particles, illuminated by a laser sheet created by a 10 W continuous laser (MX-6185, Coherent, USA) operating at 532 nm. A thin layer of fluorescent paint (FP Rhodamine 6G, Dantec) was applied to the illuminated surface, absorbing the laser wavelength to avoid reflections and to allow cor-rect near-wall measurements. The camera used was a PCO DIMAX-cs with a sampling frequency of 150 Hz. An narrowband optical filter centered around 532 nm was mounted on the camera to avoid any optical noise.

The velocity field is calculated from the acquisition of successive snapshots of the vertical laser sheet in the middle of the test section using a home-made optical-flow algorithm. The first version of the code has been developed at ONERA 5 and later modified, optimized and adapted to the constraints of real-time measurements [START_REF] Gautier | Real-time planar flow velocity measurements using an optical flow algorithm implemented on gpu[END_REF] . The advantage of this algorithm compared to a standard FFT-PIV algorithm is its high computational speed when implemented on GPUs with CUDA functions [START_REF] Cpan | Evaluating the accuracy performance of lucas-kanade algorithm in the circumstance of piv application[END_REF] . The code has been used many times both for timeresolved PIV measurements with a high spatial resolution [START_REF] Varon | Chaotic dynamics of large-scale structures in a turbulent wake[END_REF] as well as for closed-loop flow control experiments [START_REF] Gautier | Flow control using optical sensors[END_REF][11][START_REF] Gautier | Frequency-lock reactive control of a separated flow enabled by visual sensors[END_REF][START_REF] Varon | Adaptive control of the dynamics of a fully turbulent bimodal wake using real-time piv[END_REF] .

The vertical (x, y) PIV plane was chosen for the present identification because a larger experimental dataset was available. Indeed, only the boundary layer region contains relevant information for the identification algorithm so that only the lower part of the snapshots are used, which allows a significant reduction of computational time. Also, the distribution of streaks in the spanwise direction was found homogeneous for the specific forcing. The extension of the method in the horizontal plane and for different boundary layer forcing scenarios could also be investigated. vector per pixel.

III. FLOW CHARACTERISATION

In the present study we focus on a specific transition scenario of the zero pressure gradient flat plate boundary layer flow. This flow can experience various mechanisms for its transition to turbulence. In one scenario, when the incoming turbulence intensity is lower than 1% initial environment disturbances are amplified via a receptivity process 20 and are carried downstream by the flow with increasing amplitude. A breakdown of the laminar flow occurs due to a selective amplification of a narrow band of frequencies which creates the so-called Tollmien-Schlichting waves. The amplitude of the waves increases up to the point where unstable 3D peaks and valleys occur until they breakdown into hairpin eddies and turbulent spots. These wave-like structures are periodic in space and time. The critical Reynolds number based on the local displacement thickness is defined as:

Re δ * 0 = δ * 0 (x)U ∞ ν (5)
where U ∞ is the free-stream velocity, ν is the kinematic viscosity, x the streamwise coordinate (y and z being respectively the vertical and spanwise axis) and

δ * 0 = ∞ 0 (1 - u U ∞ )dy (6)
is the local displacement thickness at the inlet of the PIV window, at x 0 = 37 cm from the leading edge. The critical Reynolds number is Re crit δ * = 520. A different transition scenario may occur if the turbulence intensity of the incoming flow is order of 1% 26 or higher: the bypass transitions. Groups of streaks are usually observed, which are streamwise structures periodically distributed along the spanwise direction 1,24 , for Reynolds numbers as low as Re δ * (x) = 300. It is called bypass, because it can bypass the TS instability mechanism altogether.

The dynamics identification and coherent structure reconstruction of a subcritical bypass transition is investigated in the present experiment. The transition is triggered using a finite amplitude pulsating jet (described in the following section). Groups of wavy streaks are observed as shown on Fig. 4. The shape factor H of the forced boundary layer is H = 1.31, which is typical of a turbulent flow. The coherent structure advection velocity u adv has been estimated around 0.75 U ∞ (based on a space-time diagram of a horizontal profile of the streamwise velocity perturbation inside the boundary layer). In the following all the spatial variables are nondimensionalized using the displacement thickness at the beginning of the PIV window δ * 0 (at x 0 = 37 cm), for the unperturbed flow and noted (x * , y * , z * ). The velocity components are non-dimensionalized with the freestream velocity U ∞ and noted (u * , v * , w * ). Consequently, we refer to the dimensionless time t * = tU ∞ /δ * 0 .

A. Transition forcing mechanism

The unperturbed flow is a laminar boundary layer. The Reynolds number measured at the inlet of the PIV window is Re δ * = 409 for a freestream velocity U ∞ = 0.188 m.s -1 . No natural instability was observed without forcing. The displacement thickness at x 0 = 37 cm is δ * 0 = 2.5 mm leading to a shape factor H 0 = 2.34, typical of a laminar boundary layer.

In order to evaluate the potential of ANN on a non-linear flow, the transition to turbulence of the boundary layer is forced using a finite amplitude pulsating jet perturbation. The perturbation is induced by a spanwise slotted inclined (45 • ) jet as shown on Fig. 3. The injection system is optimized to obtain an homogeneous velocity along the spanwise direction. It is located at x = 5 cm downstream the leading edge, as shown in Fig. 3. The jet outlet has a rectangular cross-section which is 0.4 cm long (streamwise direction) and 9 cm wide (spanwise direction). The jet is supplied by a pressurized water tank, monitored by a SMC ITV 1010 electropneumatic regulator to control the jet frequency f f orcing and amplitude U f orcing . The frequency of the perturbation was random (0 < f f orcing < 5Hz), for a constant maximum amplitude. The experiments were repeated multiple times show- ing the same Proper Orthogonal Decomposition (POD) mode structures.

B. Proper Orthogonal Decomposition

Decomposing a dynamical system in modes of decreasing importance can help reducing the order of the variables of the system. N consecutive instantaneous velocity fields {U (n) = (u * , v * )} n=1...N were computed from consecutive flow snapshots with an acquisition frequency f ac = 150 Hz. By calculating the mean field [u, v] we were able to obtain the velocity fluctuations u (t) = u * (t)u and v (t) = v * (t)v, which contained all the dynamics and were used to create the reduced-order system. The fluctuation matrices organised in columns for each time-step were used to form the so-called "snapshot matrix" to be decomposed. The reduced-order system is obtained using POD, which has been used extensively in fluid mechanics [START_REF] Bore | Extended proper orthogonal decomposition: a tool to analyse correlated events in turbulent flows[END_REF][START_REF] Mendez | Pod-based background removal for particle image velocimetry[END_REF] . It allows us to build a ranked and orthonormal basis containing N modes [START_REF] Lumley | The structure of inhomogeneous turbulent flows[END_REF]30 . The first M modes {Φ m } m=1...M with M ≤ N containing a sufficient percentage of the total energy is then chosen to compute the approximated velocity field Ũ (n):

Ũ (n) = M ∑ m=1 Φ m , U (n) Φ m = M ∑ m=1 a m (n)Φ m (7)
where the scalar product •, • is the energy-based inner product. The system output to be identified is obtained through the reduced state vector containing the M POD coefficients a m (n):

Y (n) = [a 1 (n) a 2 (n) ... a M (n)] T (8)
The full-field dynamics are now contained in their POD coefficients a m (t). The balance between the order and accuracy of the POD reduced-order system is crucial, because for a large number of POD modes the SI methods are much more likely to fail.

C. Overview of the method

In the present SI method we first apply a POD to the timeresolved PIV snapshot matrix, to perform a first system order compression. Then we proceed to train a FTD ANN for the non-linear mapping between local upstream sensors and the dynamic POD coefficients. All the dynamics information is then included in a 2 easily measurable quantities. More specifically, the available time-resolved PIV dataset is divided into the learning or training part and the validation part. During the training step, the parameters of the model (weight, bias) are updated accordingly using an appropriate learning algorithm which takes into account the error between the expected output and the model-generated output via a backpropagation learning algorithm, in our case the scaled conjugate gradient. After each weight and bias update, an epoch has passed. As many epochs as needed are used to obtain a satisfying fit error. Then, during the validation step, the model generated outputs are compared with the new data, this time expecting equivalent fit results as during the training step. Larger validation fit errors yield over-fitting and make the network useless. The correctly trained network will allow the accurate prediction of the PKE at any moment using just a few local PIV sensors. The different steps of the overall identification process are summarized on the diagram shown in Fig. 5.

IV. VALIDATION CRITERION

To evaluate the efficiency of the identification, one has to define a relevant quantitative criterion to compare the POD coefficient time-series results obtained with the different ANN architectures to the ones obtained experimentally. In the present study, we compute the mean-squared error (MSE) at each time-step n for each POD coefficient a m (n):

MSE m = 1 N N ∑ n=1 (a exp,m (n) -a NN,m (n)) 2 (9)
Then the averaged MSE for all the coefficients (M = 17) time-series gives the final evaluation error for the specific ANN architecture:

MSE = 1 M M ∑ m=1 MSE m (10)
V. RESULTS

First, we obtained the full dataset of N = 6400 consecutive snapshots, acquired with a sampling frequency of 150 Hz. A POD was applied using the snapshot method. 17 POD modes were needed to obtain satisfying results with at least 80 % the total energy. The energy of the individual POD modes as well as their cumulative energy are shown in Fig. 7. The spatial structure of the third and 13 th POD modes can be found in Fig. 8 while their time evolution is shown in Fig. 9. We clearly see the similarity of the 3 rd POD mode with the largescale coherent structures visualized in the instantaneous field of Fig. 4a).

VI. SENSORS DEFINITION A. Choice of the sensor

The visual sensors will be the inputs in the identification process. The choice of the sensor(s) is then a critical step. There are many possible type of sensors as well as many possible size and locations.

The first step is to choose the physical nature of the sensor which should be based on the two components of the instantaneous 2D velocity field measured at each time-step. It can be simply one of, or the two components of the velocity field, or the velocity magnitude measured in a given window defined in the PIV field. Because of their simplicity, they are the first sensors we will test.

As we are interested in the dynamics of coherent structures advected into the boundary layer, one can also choose to compute the vorticity field. As the vortical structures are embedded into the boundary layer, it is better to use more refined de- FIG. 10. Time-series and power spectrum of velocity (s 1 ) and swirling strength (s 2 ) sensor.

tection criteria, like the Q criterion or the λ Ci criterion, which are well adapted to the identification of vortical structures inside a shear layer.

In the following, we will use the swirling strength λ Ci criterion, which was first introduced by [START_REF] Chong | A general classification of 3dimensional flow fields[END_REF] who analyzed the velocity gradient tensor D = -→ ∇ -→ u and proposed that the vortex core could be defined as a region where ∇u has complex conjugate eigenvalues. It was subsequently improved by [START_REF] Zhou | Mechanisms for generating coherent packets of hairpin vortices[END_REF] and by [START_REF] Chakraborty | On the relationships between local vortex identification schemes[END_REF] . It was also successfully applied by [START_REF] Cambonie | Transition scenario of the round jet in crossflow topology at low velocity ratios[END_REF] to vortices created by a Jet in Cross-Flow measured by Volumetric Velocimetry or by [START_REF] Gautier | Frequency-lock reactive control of a separated flow enabled by visual sensors[END_REF] in a closed-loop flow control experiment using a similar visual sensor. For 2D data, λ Ci can be computed quickly and efficiently using eq. ( 11) when such a quantity is real (else λ Ci = 0):

λ Ci = 1 2 4 det(∇u) -tr(∇u) 2 (11)

B. Number of sensors

On one hand, reducing the number of sensors leads to a reduction of the number of inputs (and outputs) which generally would make the training of the system simpler. On the other hand we may loose valuable dynamics information, helpful for the identification process. A compromise has to be found. A combination of two inputs has been tested: first, the wallnormal velocity measured at one point (averaged over the five neighbouring pixels) and second, the swirling strength vortex identification criterion λ Ci (as a sum of all the pixels in a 15 pixels-wide window). The combination of the 2 inputs was found to give better results compared to the velocity sensors alone or the swirling strength alone or the local vorticity and velocity combination. The main MSE results can be found in table.

The swirling strength window allows a good detection of every vortices passing into the boundary layer. The spectra of the two sensors (Fig. 10) show that the velocity sensor detects the large structures passing in the neighbourhood with a frequency close to 1 Hz, while the swirling strength criterion gives us more information about the high frequency dynamics up to 4.8 Hz. Some frequencies can be found in both spectra, like a clear 2.32 Hz peak. The efficiency of the combination of the two sensors can be explained through the cross-correlation of the sensor signals with the time-series signals of the POD coefficients, as illustrated in Fig. 11. One can clearly see that the 3 rd POD coefficient is weakly correlated with the swirling strength sensor, while it is strongly correlated to the velocity sensor. On the opposite the 13 th POD coefficient, which contains higher frequencies (see Fig. 9), is strongly correlated with the swirling strength sensor. This could be one of the reasons why the training using these two sensors is much more accurate than using a single sensor.

C. Position and size of the visual sensor(s)

The width of the swirling strength window also plays an important role to obtain good SI results. If it is too large it creates an unnecessary smooth event, while if too small it can be too noisy, especially for gradient variables computed from experimental data. A good compromise was found for 0.7 y * . The height of the was chosen large (5 y * ), so that it contains entirely the advected structures in the wall-normal direction, as observed in the POD modes.

The choice of the position of the sensors is also crucial. In general, the sensors should be located as upstream as possible in order to measure and predict early the downstream dynamics. The final position of the chosen sensors is shown in Fig. 6.

VII. PARAMETER SENSITIVITY STUDIES

The architecture and parameters of the ANN can be changed to avoid overfitting and improve the training process, mainly by choosing the appropriate number of neurons in the hidden layer. The depth of the network (meaning number of hidden layers) can also be increased if needed. Generally a deeper network (using the same number of total neurons) will allow the identification of more complex features and information. The drawback is that the successful training can become very tricky, especially since it increases significantly the computational time for well as the complexity of the network. This makes training parameter sensitivity analysis almost impossible.

In the FTD ANNs used in the present study the input timedelay (meaning the size of the sensor memory) can also be chosen accordingly to improve training and validation fit. It was found critical for the success of the training. The training process had a very low computational time: it was of the order of O(1) minute using a Intel Xeon E5-2630 CPU running at 2.2 GHz. This allowed for a full parametric study to find the optimal time delay k for the inputs and the minimum number of neurons for the hidden layer.

For all the ANN calculations (training, validation) the MATLAB Deep Learning Toolbox was used. A basic scheme of the final FTD ANN can be shown in Fig. 12. An anti-causal zero phase low-pass moving average (over four time-steps) filter has been applied to each pixel time-series. The maximum frequency of the forced BL dynamics was 5 Hz, so seven times lower than the low-pass filter cutoff frequency.

We mention that training with LSTM recurrent architectures was found inefficient in our case, with a strong increase of the computational time. A simple shallow Time-Delay ANN scheme with 75 neurons gave the best results. This is in contradiction with the result of Ogunmolu et al. [START_REF] Ogunmolu | Nonlinear systems identification using deep dynamic neural networks[END_REF] who have found out that deep architectures are more efficient than single hidden layer models for noisy data sets. In our case a single hidden layer model was sufficient.

The ANN used contained a tan-sigmoid transfer function. As a standard time-delay neural network it also contained an output layer with a linear transfer function to the output. The training time using a scaled conjugate gradient backpropagation algorithm was around 1 minute for a shallow network using a Intel Xeon E5-2630 CPU running at 2.2 GHz.It allowed a full parametric study to find the optimal time-delay for the inputs and the appropriate number of neurons in the hidden layer.

A. Training parameters

The full data-set is divided into training (85% of the snapshots), validation (10% of the snapshots) and over-fitting check data-set (5% of the snapshots) as summarized in Table IV. The third data-set is used as an early stopping criterion to test if a rising error is observed during the training process. As a common practice the data are shuffled randomly before the beginning of the training process. It avoids bad models if the data are initially classified. It also makes fitting faster because neighbouring points are not similar. Block consecutive data-set division failed to give good fit results. We also mention that a data standardization was applied before the training process which also helps when different inputs or outputs have different scales. Lastly, we mention that the weights and biases of the network are initialized randomly for the first epoch.

The optimum number of neurons for the hidden layer was found to be 75 and the best sensor time-delay corresponded to k = 300.8 × t snapshots (or around 4 s).

The evolution of MSE with the epoch number is shown in Fig. 13 In Fig. 14 we can see hat the MSE for each POD coeffi-

B. Influence of the sensor memory

The influence of the time-delay in the inputs, is illustrated on the plot of the evolution of the MSE as a function of the time-delay shown in Fig. 16. We observe that the sensor memory is critical for the correct dynamics identification: a minimum of k = 150 × t time units is required to obtain for the first time a correct identification with an MSE lower than 50%. A sensor memory of 300.8 × t breaks the error limit of 5%. The interpretation of the physical meaning of an optimal neural network is often useful, but not frequently studied. In our case we observe that an abrupt error drop occurs in the region of the sensor memory size of t * = 150 × t ; using the structure advection velocity we obtain a characteristic length scale of 112.2 x * , which is 1.22 times the size of our PIV domain (92 x * ). This shows that the optimal memory size is actually a function of the sensor position, the structure advection velocity and the size of the domain of which we are identifying the dynamics.

C. Influence of the number of neurons in the hidden layer

Keeping an optimal constant time delay of 600 time-steps or 300.8 t*, we can study the sensitivity of the validation errors as a function of the number of neurons the hidden layer, as shown in Fig. 17. A minimum of 70 neurons is needed to achieve an error of the order of 10 %. We also should keep in mind that increasing the number of neurons in the hidden layer increases drastically the number of parameters and the complexity of this fully connected architecture, even if the network is shallow. In Fig. 18 we show that decreasing the number of snapshots used for the training of the network leads to an increase of the number of neurons in the hidden layer needed for an equivalent validation error (keeping a constant time-delay of 301 time units t * ). Above a given number of snapshots, increasing the number of neurons does not help: it introduces overfitting hence increasing the validation error.

VIII. OPTIMAL PARAMETER RESULTS

Using the optimal time delay, the right number of neurons in the hidden layer and the proper combination of double upstream sensors, we can check the efficiency of the network in the POD coefficient signal prediction and the time-resolved instantaneous fields reconstruction. The dynamic coefficient signals are accurately predicted. In Fig. 15 it is clear that the fit error is similar for both training (left) and validation (right) datasets, which ensures we avoided overfitting. In Fig. 19 we see that the main large coherent structures are well reconstructed.

Finally, the kinetic energy of these perturbations (PKE) computed for each pixel in the PIV window is defined as: The sum of this quantity for each pixel is an interesting control quantity for the future model-free control scheme, since it is an indication of the total advected turbulent fluctuations intensity. So it is calculated and a compared, between ANN generated results and experiment. For the validation data-set, the optimal ANN architecture gives us a MSE of 10.46 %. The model and experimental time-series are in very good agreement as shown in Fig. 20.

PKE(t) = 1 2 u (t), v (t)

IX. CONCLUSION

A successful application of an ANN SI framework was presented for the case of an experimental Turbulent Boundary Layer flow. FTD ANNs were proven flexible enough to achieve sufficient prediction precision, even for a highorder multiple-input multiple-output non-linear system. A non-linear mapping between upstream sensors and the global dynamics of the flow (in the form of POD coefficients) was presented. The cross correlation between sensors and coefficient signals was proven a useful tool for the choice of the optimal sensors. Different sensors correlate well with different coefficient signals. It was revealed why the collaboration of a swirling strength sensor and a local velocity sensor resulted in an accurate network training. The dynamics were identified, the turbulent coherent structures (wavy streaks) were reconstructed, and the system order was reduced to two easily measurable signals. A physical interpretation was given for the optimal sensor memory size in the final ANN architecture.

The difficultly in the estimation of the coherent structures in this TBL flow arises due to streak wavyness and weak threedimensionality. Still the ANN was proven efficient if the sensor memory size and number of neurons are properly chosen. The achieved fit accuracy is around 5%. The present study is the first step in identifying an efficient model-free control law like in Gautier et al. [START_REF] Gautier | Closedloop separation control using machine learning[END_REF] , Li et al. [START_REF] Li | Drag reduction of a car model by linear genetic programming control[END_REF] , that comes from deep knowledge and precise prediction of the dynamics of each point in the velocity field. It will target the minimization of the PKE, since this quantity can be tracked in real-time. We also mention that the Reynolds number could also be integrated as an input, while the output could include POD coefficients for the different Reynolds numbers, hence going towards a system of variable velocity in the style of San et al. 29 .

The presented methodology could also be applied to the same turbulent flow but based on the horizontal (spanwisestreamwise) PIV plane (Fig. 4). Spanwise-distributed sensors and deeper networks are possibly needed for this complicated task but it clearly is now feasible. Different transition scenarios in the boundary layer flow can also be investigated.

X. ACKNOWLEDGEMENTS

This research was supported by CNRS and Photon Lines SAS.

The data that support the findings of this study are available from the corresponding author upon reasonable request.

FIG. 2 .FIG. 3 .

 23 FIG.2. An example of a FTD ANN to perform a non-linear mapping from N inputs to M outputs.

FIG. 4 .

 4 FIG. 4. Instantaneous streamwise velocity fluctuation field of the TBL experiment. Vertical field for z * = 0 (top) and horizontal field for plane y = 2 d * 0 (bottom). Streamwise velocity streaks can be clearly observed.

FIG. 5 .

 5 FIG. 5. Block diagram showing the different steps used for the experimental perturbation energy identification through POD and training/validation steps.

FIG. 6 .

 6 FIG. 6. Location of the double visual sensors shown over the time-averaged unperturbed streamwise velocity field u.

FIG. 7 .

 7 FIG.7. Energy of each POD modes together with the cumulative total energy of the POD modes up to 80 % of the total energy.

5 FIG. 8 . 3

 583 FIG. 8. 3 rd (a) and 13 th (b) POD mode structures.

2 FIG. 9 .

 29 FIG. 9. Evolution of 3 rd (a) and 13 th (b) dynamic POD coefficient.

  FIG. 11. Cross-correlation of the velocity (s1) and swirling strength (s2) sensors with 3rd (a,b) and 13th (c,d) POD coefficient.

FIG. 12 .

 12 FIG. 12. Structure of the FTD ANN used for this study, with 2 inputssensors (v and λ Ci ), a k time-delay, 1 hidden layer with 75 neurons and M outputs (in our case M = 17 POD modes).

FIG. 13 .FIG. 14 .

 1314 FIG. 13. Training and validation mean-squared error for each epoch of the training process.

4 FIG. 15 .

 415 FIG. 15. Comparison of the 1st to 6th POD coefficients time-series obtained from experimental data with the ones obtained with the ANN for training(a) and validation (b) datasets. The network zero signal due to time-delay is evident in the training set.

FIG. 16 .

 16 FIG. 16. Network validation performance as a function of the timedelay in the sensors. Increasing the size of the sensor memory decreases the error in a linear way, after approximately 150 time units.

8 FIG. 17 .FIG. 18 .

 81718 FIG. 17. Network validation performance as a function of the number of neurons in the hidden layer.

FIG. 19 . 6 FIG. 20 .

 19620 FIG. 19. Instantaneous horizontal fluctuation velocity field comparison between experiment (a) and ANN-reconstructed field (b). Movie online. Time-steps are not consecutive due to the dataset shuffling process during the training-validation division.

TABLE I .

 I The parameters of the PIV experiments used in the following are summarized in table I. One should keep in mind that the optical flow algorithm leads to dense vector fields with one PIV experiment parameters.

	PIV length(cm)	23.39
	Snapshot size	162x1978
	Spatial resolution 0.116 mm/pixel
	# Snapshots	6400
	Acquisition time(s)	41.3

TABLE II .

 II Comparison of the validation data-set fit error obtained with a single velocity input, a single swirling strength input and their combination. Using the combination of the inputs minimises the MSE and the number of neurons.

	n 1 MSE Sensor choice
	75 0.05	[ λ ci , v ]
	110 0.11	v
	90 0.32	λ ci

TABLE III .

 III Neural network training parameters configuration

TABLE IV .

 IV for both training and validation. The MSE decreases sharply for the first epochs and more smoothly for increasing epochs. The optimum architecture and training parameters are summarized in the table IV. Final choice for the neural network parameters.

	Network layer structure	2-75-17-17
	Activation function	Hyperbolic Tangent
	Loss function	MSE
	Training method	Scaled Conjugate Gradient
	Time-delay(s)	4

Beneitez M, Duguet Y, Schlatter P, Henningson DS (2019) Edge tracking in spatially developing boundary layer flows. Journal of Fluid Mechanics 881:164-181, doi:10.1017/jfm.2019.763


Antonios Giannopoulos

Jean-Luc Aider
email: jean-luc.aider@espci.fr

Prediction of the dynamics of a backward-facing step flow using Focused Time-Delay Neural Networks

Keywords: Backward-facing step flow, Neural Networks, Particle Image Velocimetry, Machine learning, System Identification

The full field dynamics of a separated, noise-amplifier flow, the Backward-Facing

Step at Re h = 1385, have been identified by probe-like, upstream measurements using an artificial Neural Network. Local visual sensors, coming from time-resolved Particle Image Velocimetry, were used as inputs and the dynamic Proper Orthogonal Decomposition coefficients were defined as goals-outputs for this non-linear mapping. The coefficients time-series were predicted and the instantaneous velocity fields were reconstructed with satisfying accuracy. The choices of inputs-sensors, training data-set size, hidden layer neurons and training hyperparameters are discussed for this experimental fluid system.

Introduction

Noise-amplifier flows are fluid systems which are globally stable, but which selectively amplify the upstream perturbations coming from random environment noise by convective instability mechanisms [START_REF] Huerre | Local and global instabilities in spatially developing flows[END_REF][START_REF] Chomaz | Global instabilities in spatially developing flows: nonnormality and nonlinearity[END_REF]. Typical examples are the flat plate Boundary Layer (BL) and Backward-Facing Step (BFS) flows [START_REF] Beaudoin | Three-dimensional stationary flow over a backwards-facing step[END_REF][START_REF] Dergham | Stochastic dynamics and model reduction of amplifier flows: the backward facing step flow[END_REF]. Noise-amplifier flows play an important role in many industrial flows, like separated flows around airfoils [START_REF] Darabi | Active management of naturally separated flow over a solid surface. part 2. the separation process[END_REF] and the complex 3D wakes of ground vehicles [START_REF] Aider | Experimental and numerical investigation of the flow around a simplified vehicle model[END_REF][START_REF] Beaudoin | Cavitation as a complementary tool for automotive aerodynamics[END_REF]. Understanding and controlling these flows will be crucial for the drag minimisation and reduction of green-house gas emissions, which have to be decreased in the European Union of at least 40 % by 2030, compared to the levels of 1990 [START_REF] Eu-Council | European council, european council[END_REF]. A large portion of these emissions are due to ground vehicles, and a large portion of the ground vehicle emissions is due to aerodynamic drag [START_REF] Aider | Experimental and numerical investigation of the flow around a simplified vehicle model[END_REF][START_REF] Grandemange | Characterisation of the flow past real road vehicles with blunt afterbodies[END_REF]. Controlling the shear layers to reduce the wake of a bluff-body or a ground vehicle has been proved to be an efficient strategy to reduce the aerodynamic drag of ground vehicles [START_REF] Aider | Drag and lift reduction of a 3D bluff-body using active vortex generators[END_REF][START_REF] Aider | Active flow control using pulsed micro-jets on a full-scale production car[END_REF][START_REF] Eulalie | Active flow control analysis at the rear of an suv[END_REF][START_REF] Li | Feedback control of bimodal wake dynamics[END_REF][START_REF] Grandemange | Characterization of the flow past real road vehicles with blunt afterbodies[END_REF].

In modern experimental and testing / measuring techniques, data-driven methods are becoming of great interest, since they don't require a priori knowledge of a model and since the data-sets available are becoming larger and larger. This is especially true when studying thoroughly non-stationary flows for different Reynolds numbers with various sensors such as Particle Image Velocimetry (PIV) or multiple Pitot / multiple hot-wire probes. The increasingly powerful algorithms and computers let us handle and process quickly a large amount of data. Methods like statistical /regression and machine learning algorithms, supervised or unsupervised, are becoming efficient and reliable for both academics or industrial applications. Neural networks particularly, are nowadays in the center of attention in this machine learning revolution we experience.

For complex flows, the number of degrees of freedom obtained from a 2D-2C (2-Components in a 2D velocity field) optical-flow Particle Image Velocimetry (PIV) measurement of a few millions pixels image is millions. Such a large system is impossible to control; a reduced-order model has to be identified. A dynamic observer can identify such a model based only on input-output measurements from measurable system quantities. As proposed firstly by [START_REF] Guzmán | A dynamic observer to capture and control perturbation energy in noise amplifiers[END_REF] and verified experimentally for PIV data by [START_REF] Varon | Experimental application of a dynamic observer to capture and predict the dynamics of a flat-plate boundary layer[END_REF][START_REF] Varon | Closed-loop control of separated flows using real-time piv[END_REF] it is possible to predict the full dynamics of a transitional flat-plate BL flow, in the form of Proper Orthogonal Decomposition (POD) coefficients, from local upstream measurements. The first step in their method was to create a successful reduced-order system using POD. The second step was to identify an optimal state-space model using a statistical learning process (the so-called N4SID algorithm), in order to predict at any moment all the POD coefficients (output) by measuring one or two local variables upstream (sensors or inputs) in the flow. A similar approach was also presented in a paper from Beneddine et al. [START_REF] Beneddine | Conditions for validity of mean flow stability analysis[END_REF], where the full frequency spectrum was obtained from local frequency information of the flow. Time-resolved field reconstruction was also successfully obtained for time-resolved PIV data of a round jet flow using a point sensor and the mean field from [START_REF] Beneddine | Unsteady flow dynamics reconstruction from mean flow and point sensors: an experimental study[END_REF].

In the present study, we explore means of performing this step of local-toglobal dynamics system identification (SI) using a Neural Network (NN) architecture. In this machine-learning data-driven identification process, we use a given data-set to learn the relation between local sensors and the global fluctuation dynamics of the system (in the form of POD coefficients) and a new data-set to validate the learning-training step. The identified NN coefficients will then allow the reconstruction of the full flow field, which would help design controllers targeting the kinetic energy of the full perturbation field. We will show the importance of the various choices (from sensors, to NN parameters) in a successful SI of an experimental time-resolved separated flow.

The paper is organised as follows. First, we present the various artificial NNs that can be found in the literature, and a brief review on their use in fluid mechanics problems. Then we present the experimental setup used to study the BFS flow. The choice of the NN architecture is then be discussed, together with the choice of the sensors characteristics and the training data size. The efficiency of the chosen shallow NN architecture for such a SI is then demonstrated through the time-resolved velocity field reconstruction before turning to the conclusions.

Artificial Neural Networks

Definitions

Artificial Neural Networks (NN) are a subset of supervised Machine Learning. They can provide a non-linear mapping between one set of inputs and a set of corresponding outputs. Great progresses have been made lately because of the availability of large data-sets, the increasing availability of multiple optimised toolboxes and also the progress in Graphics Processing Units (GPU) parallel programming, which improves the computational speed. This is the reason why NNs are becoming more and more popular nowadays.

The fundamental part of a NN is the neuron or "perceptron". In general, one defines a weight w i associated to each i th neuron of the previous layer. To obtain the output of a perceptron from a set on N neurons, one computes the sum of the N inputs multiplied by their corresponding weight and adds a given bias b i . An activation function f is then used to compute the output. A classic activation function is a step-function, but more refined functions are usually needed. To improve the efficiency of classic NNs, it is possible to add one or many "hidden layers" between the inputs and the outputs. Each hidden layer represents a non-linear function of linear combinations of its inputs, using a weight and a bias for each input. The inputs are then coming from the previous hidden layer. It is theoretically proven that any continuous function can be approximated with a single hidden layer [START_REF] Cybenko | Approximation by superpositions of a sigmoidal function[END_REF][START_REF] Funahashi | On the approximate realization of continuous mappings by neural networks[END_REF].

The simplest shallow, fully-connected NN architecture consists of the input layer (with n neurons), a single hidden layer (with an arbitrary number of n 1 neurons) and a linear layer (with n 2 = m neurons) connected to the output, as shown in Fig. 1. A standard rule for the linear layer is to have the same number of neurons as the output layer. Regarding the choice of the number of neurons in the hidden layer, it can be as high as needed to increase accuracy, but without over-fitting. For a NN with a non-linear activation function f 1 in the hidden layer and a linear activation function f 2 in the linear layer, the equation giving the k th neuron output of a single hidden layer network connected to the j th neuron of the previous layer is:

y k = f 1   n2 j=0 w (2) kj f 2 ( n1 i=0 (w (1) kj x i + b i ) + b j   (1)
where n 1 is the number of neurons in the first (hidden) non-linear layer and n 2 the number of neurons in the second (linear) layer. Regarding the choice . . . of the non-linear activation function f 1 , usually the popular tan-sigmoid or hyperbolic tangent function is used:

f 1 (x) = e 2x -1 e 2x + 1 (2)
The correct training process of the network (i.e. finding the optimal weights and biases connecting the neurons of different layers) consists of dividing the data-set into a training data-set and a validation data-set. Starting with the training data-set, the first set of weights connecting the layers is randomly initialised for this first iteration. The error of the real vs the model-generated output is computed and the weights and biases are updated according to different back-propagation schemes, in our case the Scaled Conjugate Gradient (SCG) method. In this case, the step size is adjusted at each iteration in order to minimise the performance function. The above process is called an "epoch". We continue the process for as many epochs as needed until a satisfactory fit error is achieved (Fig. 2). The second data-set, is the "validation" data-set and is used to evaluate the performance of the network on new data and verify the achieved error, hence avoiding over-fitting.

Neural Network types

NNs can be divided into Feedback (or recurrent) and Feed-Forward. They can be discriminated according to their depth, either shallow or deep, depending on the number of hidden layers (one or more). Finally, they can also be divided into static or dynamic, if the output of the current step depends on the previous steps as well, giving it a notion of memory. In the case of a feed-forward NNs the output of any layer only modifies the next layer. System identification using feed-forward NNs have been a common practice since the 90s [START_REF] Narendra | Identification and control of dynamical systems using neural networks[END_REF][START_REF] Chu | Neural networks for system identification[END_REF]. Recently, deeper feed-forward NNs have also been implemented in complex SI cases [START_REF] Olekan | Nonlinear systems identification using deep dynamic neural networks[END_REF].

The network used in the present SI study is a fully-connected Focused Time-Delay Neural Network (FTDNN) and is basically a standard feed-forward architecture along with a tapped time-delay (of time-step size k) in the input. The term "focused" comes from the fact that the notion of memory is introduced only in the input and not in the output. They were first introduced for speech recognition [START_REF] Waibel | Phoneme recognition using time-delay neural networks[END_REF]. They are used to model long-range temporal dependencies by keeping a number of past measurements of the input at each k previous time steps x t-k leading to the following expression for the output of the system [START_REF] Charaniya | Focused time delay neural network model for rainfall prediction using indian ocean dipole index[END_REF]:

y t = f w,b (x t , x t-1 , ..., x t-k ) (3)
where w and b are weight and bias parameters. On the other hand, in a recurrent NN, the system output is calculated from its previous past time-steps along with the input at the current time-step like in equation 4:

y t = f w,b (y t-1 , x t ) (4)
Recurrent NNs also introduce a notion of memory in the output of the system. One specific family of recurrent NN is the Non-linear Auto-Regressive eXogenous (NARX) models. They are autoregressive because the outputs of the current time step depends on the output of a number of previous steps. They are exogenous because the output depends also on a number of inputs. The NARX models were first introduced by [START_REF] Leontaritis | Input-output parametric models for nonlinear systems part i: deterministic non-linear systems[END_REF] and used with NN with significant success by [START_REF] Chen | Non-linear system identification using neural networks[END_REF] for multiple non-linear SI cases. It is a natural extension of the Autoregressive Exogenous model (ARX), which has been extensively used in linear black-box SI. Another category of recurrent NNs are the Long-Short Term Memory (LSTM) networks. They are gradient-based recurrent NNs used for supervised learning both in classification and prediction problems [START_REF] Vlachas | Datadriven forecasting of high-dimensional chaotic systems with long-short term memory networks[END_REF]. They were first developed by [START_REF] Hochreiter | Long short-term memory[END_REF] to solve the vanishing or exploding gradient problem of the back-propagated error. The difference with the FTDNNs, which also introduce an arbitrary long-term time dependence of the input to its previous moments, is that the network is left to learn alone the size of the memory of each neuron during the training process. They do so using a sophisticated gate-neuron, that determines if the input is important enough or if it should be forgotten and when it should output its value. They are more computationally expensive to train than the FTDNNs.

Neural Networks in Fluid Mechanics

In fluid systems, feed-forward artificial NNs have been used for data-driven reduced-order modelling [START_REF] Müller | Application of machine learning algorithms to flow modeling and optimization[END_REF][START_REF] San | Neural network closures for nonlinear model order reduction[END_REF][START_REF] Wang | Model identification of reduced order fluid dynamics systems using deep learning[END_REF][START_REF] Mohan | A deep learning based approcah to reduced order modeling for turbulent flow control using lstm neural networks[END_REF][START_REF] Pan | Long-time predictive modeling of nonlinear dynamical systems using neural networks[END_REF] with many results showing better field reconstruction than traditional POD methods [START_REF] San | An artificial neural network framework for reduced order modeling of transient flows[END_REF][START_REF] Erichson | Shallow learning for fluid flow reconstruction with limited sensors and limited data[END_REF]. They were also used by [START_REF] Mi | Vertical two-phase flow identification using advanced instrumentation and neural networks[END_REF] for experimental flow regime identification in multiphase flows. We mention that it has also been proved by [START_REF] Baldi | Neural networks and principal component analysis: Learning from examples without local minima[END_REF], that a linear NN can be equivalent to a POD basis structure. Convolutional NNs have also been used for the efficient real-time 2D and 3D inviscid simulations [START_REF] Tompson | Accelerating eulerian fluid simulation with convolutional networks[END_REF] or along with pressure measurements for the velocity field prediction around a cylinder [START_REF] Jin | Prediction model of velocity field around circular cylinder over various reynolds numbers by fusion convolutional neural networks based on pressure on the cylinder[END_REF].

Furthermore, if there are hints or deeper understanding of the underlying physics, simple shallow NNs can provide very good results in SI [START_REF] Guzmán | A dynamic observer to capture and control perturbation energy in noise amplifiers[END_REF][START_REF] Huang | Control and system identification of a separated flow[END_REF] as well as for control laws creation [START_REF] Lee | Application of neural networks to turbulence control for drag reduction[END_REF][START_REF] Herbert | Laminar flow control with neural net[END_REF][START_REF] Rabault | Artificial neural networks trained through deep reinforcement learning discover control strategies for active flow control[END_REF]. Their potential has been demonstrated early for modelling surface pressure and aerodynamic coefficients of 3dimensional unsteady cases of aircrafts flows [START_REF] Faller | Neural networks: Applications and opportunities in aeronautics[END_REF]. Finally, deep NNs are increasingly more important in the fluid mechanics community, especially for the modelling of complex turbulent flows. Srinivasan et al. [START_REF] Srinivasan | Predictions of turbulent shear flows using deep neural networks[END_REF] compared deep feed-forward and recurrent LSTM networks for turbulent shear flows prediction. Recently Deng et al. [START_REF] Deng | Time-resolved turbulent velocity field reconstruction using a long short-term memory (lstm)-based artificial intelligence framework[END_REF], used LSTM networks to reconstruct the POD coefficient time series using sub-sampled distributed velocity sensors in an inverted flag flow PIV experiment. A short review of applications of deep NN to fluid mechanics can be found in [START_REF] Kutz | Deep learning in fluid dynamics[END_REF].

Experimental setup

Hydrodynamic channel

Experiments have been carried out in a hydrodynamic channel in which the flow is driven by gravity, with a maximum free-stream velocity U ∞ = 22 cm.s -1 . The flow is stabilised by divergent and convergent sections separated by honeycombs leading to a turbulence intensity of 0.8 %. A NACA 0020 profile is used to smoothly start the boundary layer. The test section is 80 cm long with a rectangular cross-section w = 15 cm wide and H = 10 cm high (Fig. 3). The step height is h = 1.5 cm. The maximum Reynolds numbers based on the step height is Re h,max = U ∞ h/ν ≈ 3000. The vertical expansion ratio is A y = H/(h+H) = 0.82 and the spanwise aspect ratio is A z = w/(h+H) = 1.76.

Lucas-Kanade Optical Flow PIV measurements

The flow is seeded with 20 µm neutrally buoyant polyamide particles, which are illuminated by a laser sheet created by a 2 W continuous laser (MX-6185, Coherent, USA) operating at 532 nm. A thin layer of fluorescent paint (FP Rhodamine 6G, Dantec) was applied to the illuminated surface to absorb the laser wavelength and avoid reflections and to allow correct near-wall measurements. The Camera used was a 4 Mp PCO DIMAX-cs with an acquisition frequency f ac = 150 Hz. An narrow-band optical filter was mounted on the camera to visualize only the laser light reflected by the particles. The length of the PIV window 11.2 h and its height is 3.7 h.

The time-resolved velocity fields are calculated from the acquisition of successive snapshots in the vertical symmetry plane at z = w /2, using a homemade Lucas-Kanade Optical Flow (LKOF) algorithm. The first version of the code has been developed at ONERA [START_REF] Champagnat | Fast and accurate PIV computation using highly parallel iterative correlation maximization[END_REF] and later modified, optimised and adapted to the constraints of real-time measurements by Gautier & Aider [START_REF] Gautier | Control of the separated flow downstream of a backward-facing step using visual feedback[END_REF]. Among the advantages of the LKOF algorithm compared to a standard FFT-PIV algorithm is the calculation of a dense velocity field with one vector per pixel. It also allows for high computational speed when implemented on a GPU with CUDA functions [START_REF] Pan | Evaluating the accuracy performance of lucas-kanade algorithm in the circumstance of piv application[END_REF]. The high spatial resolution is important for nearwall measurements while the high computation speed is important for real-time measurements that can be used as inputs in closed-loop flow control experiments [START_REF] Gautier | Real-time planar flow velocity measurements using an optical flow algorithm implemented on gpu[END_REF]. The code has been used many times both for time-resolved PIV measurements with a high spatial resolution [START_REF] Varon | Chaotic dynamics of large-scale structures in a turbulent wake[END_REF], as well as for closed-loop flow control experiments [START_REF] Gautier | Control of the separated flow downstream of a backward-facing step using visual feedback[END_REF][START_REF] Gautier | Closedloop separation control using machine learning[END_REF]. The PIV calculations in the present study were performed on a NVIDIA TESLA K80 GPU.

Backward-Facing Step flow

Characterisation of the BFS flow

The objective of this study is to evaluate the potential of a NN SI method on experimental data of a shear-layer flow. We focus on a BFS flow which is a typical case of noise-amplifiers. Upstream perturbations are amplified in the shear layer leading to significant downstream disturbances. Separation is imposed by a sharp edge creating a strong shear layer where Kelvin-Helmholtz instability leads to vortex shedding (Fig. 3). Another important feature is the creation of a large separation bubble, which is usually associated to pressure drag [START_REF] Dahan | Feedback control for form-drag reduction on a bluff body with a blunt trailing edge[END_REF]. Its reduction is then a common objective to most flow control experiments targeting drag reduction. It is also considered as a benchmark case for the study of separated flows. For this reason it has been extensively studied both numerically and experimentally [START_REF] Beaudoin | Three-dimensional stationary flow over a backwards-facing step[END_REF][START_REF] Hung | Direct numerical simulation of turbulent flow over a backward-facing step[END_REF][START_REF] Armaly | Experimental and theoretical investigation of backward-facing step flow[END_REF]. The Reynolds number in the present experiment is Re h = 1385, corresponding to a free-stream velocity u ∞ = 11 cm.s -1 . The incoming boundary layer, downstream from the leading edge, is laminar and follows a Blasius profile. The boundary layer thickness just upstream the step edge is δ 0 = 7mm = 0.5 h corresponding to a shape factor H = 2.3, typical of a laminar boundary layer. The vortex shedding frequency is f shed = 1.13Hz, which corresponds to a Strouhal number St h = f shed h U∞ = 0.15. As we are interested in the growth and dynamics of coherent structures, one can also choose to monitor the vorticity field. Since vortical structures are embedded into the shear layer, it is better to use more refined criteria, like the Q criterium or the λ Ci criterium, which are well adapted to the identification of vortical structures inside a shear layer. The λ Ci criterium was first introduced by Chong et al. [START_REF] Chong | A general classification of 3-dimensional flow fields[END_REF], who analysed the velocity gradient tensor D = -→ ∇ -→ u and proposed that the vortex core could be defined as a region where ∇u has complex conjugate eigenvalues. It was subsequently improved by Zhou et al. [START_REF] Zhou | Mechanisms for generating coherent packets of hairpin vortices[END_REF] and by Chakraborty et al. [START_REF] Chakraborty | On the relationships between local vortex identification schemes[END_REF]. It was also successfully applied by [START_REF] Cambonie | Transition scenario of the round jet in crossflow topology at low velocity ratios[END_REF], to visualise the 3D vortices created by a Jet in Cross-Flow measured by Volumetric Velocimetry, and by Gautier et al. [START_REF] Gautier | Frequency-lock reactive control of a separated flow enabled by visual sensors[END_REF] in a closed-loop flow control experiment using a similar visual sensor. For 2D data, λ Ci can be computed quickly and efficiently using equation ( 5), when such a quantity is real (else λ Ci = 0):

λ Ci = 1 2 4 det(∇u) -tr(∇u) 2 (5)

Proper Orthogonal Decomposition

Decomposing a dynamical system in modes of decreasing importance can help reducing the order of the variables of the system. N = 4197 consecutive velocity fields {U (n) = (u *

x , u * y )} n=1...N were computed from consecutive flow snapshots ,with a sampling frequency of f ac = 150 Hz. The size of each snapshot-velocity field is 346 × 1010 pixels, with a spatial resolution of 0.166 mm/pixel . By calculating the mean field [U x , U y ] we were able obtain the velocity fluctuations u x (t) = u * x (t) -U x and u y (t) = u * y (t) -U y , which contained all the dynamics and were used to create the reduced-order system. The fluctuation matrices organised in columns for each time-step were used to form the so-called "snapshot matrix" to be decomposed. The reduced-order system is obtained using POD, which has been used extensively in fluid mechanics [START_REF] Borée | Extended proper orthogonal decomposition: a tool to analyse correlated events in turbulent flows[END_REF][START_REF] Mendez | Pod-based background removal for particle image velocimetry[END_REF]. It allows us to build a ranked and orthonormal basis containing N modes [START_REF] Lumley | The structure of inhomogeneous turbulent flows[END_REF][START_REF] Sirovich | Turbulence and the dynamics of coherent structures[END_REF]. The first K modes {Φ k } k=1...K with K ≤ N containing a sufficient percentage of the total energy is then chosen to compute the approximated velocity field Ũ (n):

Ũ (n) = K k=1 Φ k , U (n) Φ k = K k=1 a k (n)Φ k (6)
where the scalar product •, • is the energy-based inner product. The system output to be identified is obtained through the reduced state vector containing the K POD coefficients a k (n):

Y (n) = [a 1 (n) a 2 (n) ... a K (n)] T (7)
The full-field dynamics are now contained in their POD coefficients a k (t). The balance between the order and accuracy of the POD reduced-order system is crucial, because it was seen that for a large number of POD modes the SI methods are much more likely to fail.

The energy of the individual POD modes as well as their cumulative energy are shown in Fig. 5. One can see that 50 % of the energy is contained in the three first modes. The system size containing at least 80 % of the total energy was used; that corresponding to K = 10 modes Φ k and 10 POD coefficients a k . The evolution of the 1 st , 3 rd and 5 th POD coefficients are shown in Fig. 6. The characteristic frequency of the 1 st POD mode corresponds to the Kelvin-Helmholtz frequency, i.e. the shedding of the vortices in the shear layer.

System Identification steps

In the present study we explore the potential of Artificial NNs for local-toglobal dynamics SI applied to a shear layer flow. A full scheme of the identification process is summarised in Fig. 7. First, in the full data-set of time-resolved PIV experiment is decomposed to identify the dynamics in the form of POD coefficients. Then, in this data-driven identification process, we just rely on the input (sensors) -output offline measurements for a period of time from the operation of the plant. Once the machine (in our case an artificial NN) has monitored a sufficient number of realisations (training data-set), it will identify the relationship between the given input-sensor and the goal-output. If the method is successful, then monitoring the sensor will allow us to predict successfully the full global dynamics (in the form of POD coefficients) of the system in a new, unknown to the machine, data-set (validation), with no further need for a time consuming field decomposition analysis. Then the time-resolved field can be reconstructed from the identified NN-generated POD coefficients, using only the local sensor measurements.

In the following sections we discuss the influence of the main parameters on the SI process. We will especially show the importance of the choice of the physical parameters (nature, number and location of the sensors) in the ability of the SI learning process to find the proper NN weights and biases for each neuron. We will also show how the choice of the NN parameters (number of neurons in the hidden layer, time-delay) plays a crucial role in the efficiency of the SI process.

Validation Criterion

To evaluate the efficiency of the identification, one has to define a quantitative criterium to compare the POD coefficient time-series results obtained with the different NN architectures to the ones obtained experimentally. In the following, we compute the mean-squared error (MSE) at each time-step n for each POD coefficient a k (t):

M SE i = 1 N N n=1 (a exp,i (n) -a N N,i (n)) 2 (8)
Then the averaged MSE for all the coefficients (m = 10) time-series gives the final evaluation error for the specific NN architecture:

M SE = 1 m m k=1 a k (9)
7. Influence of the sensors definition

Choice of the sensor(s)

First, it is necessary to choose the physical quantities measured by the sensor(s). As the visual sensors are the inputs in the identification process, their choice is a critical step. The sensor(s) are necessarily based on the measurements of the two components of the instantaneous 2D velocity fields at each time-step. Our challenge is to identify the birth and growth of vortical structures in the reciptivity process of the noise amplifier; in our case inside the shear layer right after separation. For this reason, possible sensors can be one or the two components of the fluctuation velocity field, the velocity magnitude or other vortex identification criteria, like the the λ ci criterion.

Sensor position

The location of the sensor is critical. It should allow the detection of the perturbations during the initial phase of the reciptivity-amplification process of the shear layer [START_REF] Varon | Closed-loop control of separated flows using real-time piv[END_REF][START_REF] Guzmán | A dynamic observer to capture and control perturbation energy in noise amplifiers[END_REF]. So in such a noise-amplifier flow the sensor was placed as upstream in the velocity field as possible, right after the separation in the highshear region, as shown in Fig. 4. The proximity to the wall also is important for the method to be realistically applicable if other measuring devices, like hot-wires, were to be tested in the future. The exact sensor location was finally chosen to be kept slightly away from the wall (x = 0.25 h), in order to avoid possible noisy near-wall measurements. We mention that the flow field may be difficult to measure experimentally because of the large velocity gradients in the sensor region. Nonetheless, the good spatial resolution of the LKOF PIV measurements allows for the computation of gradient-based quantities like the vorticity or the swirling strength criterion.

Number and size of sensors

Reducing the number of sensors, as well as reducing the number of outputs, generally should make the training of the network simpler (lower order multiple input-multiple output regression). On the other hand, using less sensors may lead to a loss of valuable information about the flow, so a compromise has to be found. Single sensor configurations have been tested, i.e. either only vertical fluctuation velocity u y or λ ci , as well as combinations of two or three sensors. The velocity component sensor s 1 is computed as an average in a five neighbouring pixels window. The swirling strength vortex identification sensor s 2 is defined as the sum of λ Ci computed over all the pixels in a dx = 15 pixels-wide window:

s 2 = y2 y1 x2 x1 λ Ci dxdy ( 10
)
The height of the swirling strength window d y has to be close to the thickness of the shear, so that the vortex activity can be properly monitored. A good compromise has been obtained with dy = 0.5 h. The width of the swirling strength window also influences the quality of identification results. If it is too large, it creates an unnecessary smooth result, while if too small it can be too noisy, especially for gradient variables computed from experimental data. Finally, a window width of d x = 0.15 h was used. Table 1: Comparison of the validation data-set fit error obtained with a single velocity sensor, a single swirling strength sensor and their combination. Using both sensors minimises the MSE and the number of neurons.

Choice of the sensor

The SI error obtained using the different sensors separately or combined can be found in Table 1. One can see that the choice of the sensor(s) is indeed critical. Using a single swirling strength sensor leads to a large MSE (37.87%) together with a large number of neurons (105). Using a single velocity sensor reduces the MSE (8.5%) but still needs a large number of neurons (90). Finally, the best results are obtained with the combination of the two sensors: the MSE is minimum (6.5%) and the number of neurons is strongly reduced [START_REF] Champagnat | Fast and accurate PIV computation using highly parallel iterative correlation maximization[END_REF].

This result can be explained. Indeed, both sensors being placed inside the shear layer just after separation, they contain a lot of information. This region is rich in frequencies coming from the the initial receptivity and amplification process together with the shear layer instability. This is also a region where measurement noise coming from the PIV measurement will be maximum. Using both sensors helps the identification from the NN to be successful in separating the physical from the unwanted frequencies, leading to a strong improvement of the MSE fit (6.5%).

Results and discussion

Optimal NN identification procedure

The goal of the NN identification method is to predict at each time step t the POD coefficients a k (t) of the full field using local upstream sensors s j (t). The POD coefficients time-series have been calculated based on the PIV velocity fields. The sensors s j (t) were also monitored at the same time steps. The pairs [a k (n), s j (n)] (n = 1 : 4197) is our identification data-set.

A FTDNN algorithm has been chosen. Our objective was then to identify the optimal sensors s j , the time delay k and the number of neurons in the hidden layer n 1 . A basic scheme of the network is shown in Fig. 1. An anti-causal zero phase low-pass moving average (over four time-steps) filter has been applied to each pixel time-series, to reduce slightly the measurement noise. We mention that the cut-off frequency of the filter is more than 10 times the vortex shedding frequency. For all the NN calculations (training and validation) the MATLAB Deep Learning Toolbox was used [START_REF] Beale | Deep learning toolbox user's guide[END_REF].

A NN should be as efficient as possible (according to the chosen criteria) and at the same time as simple as possible to easily check robustness and reduce computational time. Simplicity means minimising the number of layers and the number of neurons in each layer. The FTDNN networks tested used a tansigmoid transfer function and had a single hidden layer. In this case the NN is considered as "shallow". The computational time for the training process, using a scaled conjugate gradient back-propagation algorithm, was of the order of O(1) minute using a Intel Xeon E5-2630 CPU running at 2.2 GHz. The low computational time allowed a full parametric test to find the optimal time delay k for the input and number of neurons n 1 for the hidden layer.

For each NN architecture, the full data-set [a k (t), s j been divided into training (85% of snapshots), validation (10% of snapshots) and over-fitting check data-set (5% of snapshots). The third data-set is used to monitor if high over-fitting is observed during training. In this case the process is stopped immediately. As a common practice, the data-set is shuffled randomly before the beginning of the training process. It avoids bad models if the data are initially classified. It also makes fitting faster because neighbouring points are not similar. Block consecutive data-set division failed to give good fit results. We also mention, that a data standardisation was applied before the training process to help handling different inputs or outputs with different scales. Lastly, the weights and biases of the network are initialised randomly for the first epoch. No early stopping criterion has been used. The number of neurons of the hidden layer n 1 was always changed iteratively, in order to find the optimal number leading to a minimization of the meansquared error (MSE).

The same parametric study was carried out for the time-delay k, leading to a double loop parametric study from which the optimal combination of [k, n 1 ] was found. The optimal number of neurons for the hidden layer was found to be n 1 = 50 and the best sensor time-delay corresponded to k = 400 time steps (2.66 s). A more detailed discussion about the time-delay parameter can be found below.

The results of the training performance as a function of the increasing epochs of the optimal architecture are shown in Fig. 9. 900 epochs were proved sufficient for a validation MSE error lower than 10 % for all POD coefficients. The optimal architecture and training parameters are summarised in the Table 2. This simple shallow FTDNN architecture (see Fig. 8) avoids the complexity of deeper architectures which increases the computational time for the training process and the need for a large amounts of training data [START_REF] Canziani | An analysis of deep neural network models for practical applications[END_REF], which are often difficult to obtain experimentally. This result verifies that shallow architectures can still provide satisfying results in many fluid mechanics applications [START_REF] Srinivasan | Predictions of turbulent shear flows using deep neural networks[END_REF][START_REF] Erichson | Shallow learning for fluid flow reconstruction with limited sensors and limited data[END_REF][START_REF] Tracey | A Machine Learning Strategy to Assist Turbulence Model Development[END_REF].

Influence of the Time-Delay

The time-delay k in a FTDNN gives a notion of constant size memory in the sensors-inputs. It is crucial for the identification of downstream dynamics, while our sensors contain only information about the upstream dynamics (in the Kelvin-Helmholtz instability region). A sensitivity parametric study regarding the optimal time-delay in the inputs can be found by looking at the plot of the MSE of the validation data-set, as seen in Fig. 10. Keeping a time delay of k = 400 steps or 2.66 s, we can achieve a validation fit error less than 7%. k = 400 time-steps correspond in roughly 3 oscillations of the 1st mode ( which contains the highest amount of energy and corresponds to the Kelvin-Helmoltz vortex shedding). This notion of "memory" of these three events is crucial for the success of the algorithm. It explains why FTDNN with less memory or a standard feed-forward NN with zero time delay (k = 0) failed to identify the dynamics (M SE > 90%). An intermediate transitional region with k = 170 up to k = 300 might actually give good identification results, but it strongly depends on the random weight-bias initialisation and random shuffling of data. Convergence is achieved only for k > 300 time-steps.

Influence of the size of the training data-set

As expected for any data-driven method, the more informations are given to the machine the better. It will figure out more easily the optimal weights and biases for the mapping we ask it to perform. In Fig. 11 we can see that decreasing the number of snapshots used for the training of the network leads to an increase of the number of neurons in the hidden layer for an equivalent validation error (keeping a constant time-delay of k = 400 time steps). Above a given number of snapshots, increasing the number of neurons does not help, it just introduces over-fitting hence increasing the validation error. The BFS flow PIV experiment for Re h = 1385 is dominated by the vortex shedding frequency (more than 40% of energy in the first two POD modes). Comparing the number of events (or oscillations) with the accuracy of the identification shows that to obtain a MSE lower than 10% we need to record at least 20 events.

Reconstruction of the instantaneous velocity fields

In Fig. 12 we present the results obtained for the ten first POD coefficients with the optimal double sensor and optimal time-delay configuration. We can see that the MSE for each POD coefficient is equivalent for all the coefficients at both training and validation, which proves we avoided over-fitting. In Fig. 13, we also present the time-series of three POD coefficients calculated from the experiment compared with the time-series obtained in the training step (a) and in the validation step (b) using the optimal NN architecture. Both training and validation data are in very good agreement with the experimental data and training and validation fit errors are very similar.

Projecting the NN-predicted POD coefficients a i to the POD modes Φ i allows the reconstruction of the fluctuation velocity fields of all time steps of the validation data-set. It is then possible to compare the predicted velocity field with the real experimental fields. The instantaneous experimental velocity field is compared to the one reconstructed with the optimal double visual sensor information as shown in Fig. 14 and corresponding video. In general, all the main features from the vortex shedding are apparent with correct amplitudes throughout the velocity field. Some smaller structures in the shear layer have reduced amplitude, due to the 80% energy cut during the time-resolved field decomposition, as well as the NN identification process error itself.

Conclusions

A successful application of a NN System Identification method to a timeresolved PIV experiment of a typical noise-amplifier flow has been presented. We were able to predict with satisfying precision the global dynamics of the flow (in the form of POD coefficients), using only two upstream visual sensors coming from local velocity measurements. A shallow FTDNN architecture was sufficient to recover the overall dynamics of the flow. There was no need for sophisticated LSTM gates or more than one hidden layers, which would increase the training complexity and the computational time.

It demonstrates the feasibility to reduce the order of such fluid systems from O(10 6 ) (for typical LKOF PIV measurements) to only a handful of useful variables, which is crucial for control purposes. The final NN architecture allows us to predict the dynamics of the flow using local, upstream visual probes, without the need for time-consuming full-field decomposition analysis or intrusive measuring devices like hot-wires or Pitot tubes. It was demonstrated that the combination of the swirling strength and local velocity sensors lead to a satisfying training and validation fit in the predicted POD coefficient timeseries, even though these sensors were located in a high-gradient and difficult to measure region of the flow (early in the receptivity-amplification process of the noise-amplifier). The accuracy of the method was also illustrated in the instantaneous velocity field reconstruction. The NN training process was found to be very fast on a standard desktop computer (O(1) minutes). The double sensor approach is simple and fast to calculate and is ideal for a model-free closed-loop control scheme, like in [START_REF] Gautier | Closedloop separation control using machine learning[END_REF][START_REF] Li | Drag reduction of a car model by linear genetic programming control[END_REF], with the objective to reduce the turbulent kinetic energy of the flow.

We also mention that, in a future implementation, the Re h number could also be integrated as an input, while the output could include POD coefficients from the different Re h numbers, making it even more useful for flow control purposes [START_REF] San | An artificial neural network framework for reduced order modeling of transient flows[END_REF]. Distributed sensors, upstream and/or downstream, could also help achieve this complicated task. A larger number of POD modes could be included for more precise representation of the initial experimental data-set as well. Finally, the same methods can be also applied to other noise-amplifier flows, like the transitional flat plate boundary layer flow under the Tolmien-Schlichting instability.

Figure 1 :

 1 Figure 1: An example of a shallow NN non-linear mapping to monitor m outputs Y using the n sensors X in a FTDNN architecture with k steps of delay and n 1 neurons in the hidden layer.

Figure 2 :

 2 Figure 2: Block diagram of the NN System Identification training step.

Figure 3 :

 3 Figure 3: Sketch of the BFS geometry and main flow features (shear layer and recirculation bubble). The PIV window is shown in grey and the visual sensor location shown as a red dot.

Figure 4 :Figure 5 :

 45 Figure 4: Averaged streamwise velocity field (a), example of an instantaneous streamwise velocity field (b) and of an instantaneous swirling strength field (c) together with the position of the tested visual sensors. The flow is going from left to right. Movie online.

Figure 6 :

 6 Figure 6: Time-series of 1 st , 3 rd and 5 th dynamic POD coefficients.The evolution of the 1 st coefficient has a peak frequency close to the K-H vortex shedding frequency.

Figure 7 :

 7 Figure 7: Sketch summarising the main steps of the identification procedure: from the time resolved PIV experiment to a convincing reconstruction of the velocity fields using local upstream sensors.

Table 2 :

 2 Final choice for the optimal NN architecture and its training hyperparameters.

Figure 8 :

 8 Figure 8: Structure of the optimal FTDNN non-linear mapping, with 2 inputs-sensors (velocity and λ Ci ), k = 400 steps of time-delay in the input, 1 hidden layer with n 1 = 50 neurons and m = 10 POD coefficients as outputs.

Figure 9 :

 9 Figure 9: Training and validation MSE for each epoch of the training process.

Figure 10 :

 10 Figure 10: Evolution of the MSE as a function of the time-delay k in the sensor. Convergence is achieved for k higher than 300 time-steps.

Figure 11 :

 11 Figure 11: Network validation data-set performance as a function of the number of snapshots used for the training and the number of hidden layer neurons (n 1 ) used.

Figure 12 :

 12 Figure 12: Training and validation MSE for each POD coefficient time series.

17

 17

Figure 13 :

 13 Figure 13: Comparison of the 1st ,3 rd and 5 th POD coefficient time-series obtained from experimental data with the ones obtained with NN for training (a, c, e) and validation (b, d, f) data-sets.

Figure 14 :

 14 Figure 14: Comparison between experimental instantaneous streamwise fluctuation velocity field (a) and NN-reconstructed fields (b) using the double visual sensor and optimal NN architecture with n 1 = 70 neurons in the hidden layer and N = 2728 training snapshots. Movie online.

Acknowledgements

The authors would also like to thank M. Aris Kanellopoulos from Georgia Institute of Technology for his valuable comments and discussions.