ML for fast assimilation of wall-pressure measurements from hypersonic flow over a cone

Data assimilation (DA) integrates experimental measurements into computational models to enable high-fidelity predictions of dynamical systems. However, the cost associated with solving this inverse problem, from measurements to the state, can be prohibitive for complex systems such as transitional hypersonic flows. We introduce an accurate and efficient deep-learning approach that alleviates this computational burden, and that enables approximately three orders of magnitude computational acceleration relative to variational techniques. Our method pivots on the deployment of a deep operator network (DeepONet) as an accurate, parsimonious and efficient meta-model of the compressible Navier–Stokes equations. The approach involves two main steps, each addressing specific challenges. Firstly, we reduce the computational load by minimizing the number of costly direct numerical simulations to construct a comprehensive dataset for effective supervised learning. This is achieved by optimally sampling the space of possible solutions. Secondly, we expedite the computation of high-dimensional assimilated solutions by deploying the DeepONet. This entails efficiently navigating the DeepONet’s approximation of the cost landscape using a gradient-free technique. We demonstrate the successful application of this method for data assimilation of wind-tunnel measurements of a Mach 6, transitional, boundary-layer flow over a 7-degree half-angle cone.

the disturbance waves, and downstream measurements in the nonlinear regime are more challenging to interpret.In addition, some flow quantities are impossible to measure directly.
Unlike experiments, computational models enable non-intrusive access to all flow quantities of interest [12][13][14][15] .The shortcoming of computations, even direct numerical simulations (DNS) that resolve all the flow scales, is that they invoke modeling assumptions, most importantly idealized inflow/boundary conditions that determine the downstream flow, including laminar-to-turbulence transition 16 .Such idealizations are a main source of uncertainty since changing the oncoming waves can quantitatively, and also qualitatively, change the transition process [16][17][18][19][20] .
The field of data assimilation (DA) fuses experiments and simulations.Available measurements are adopted as targets that computations aim to reproduce quantitatively.By adjusting uncertain parameters until their optimal values are discovered, simulations become significantly more realistic.Such computations then provide non-intrusive access to all flow quantities of interest.
In fluid dynamics, DA has been applied to a wide range of flow regimes using various strategies, including filtering and nudging [21][22][23][24][25] , adjoint-variational techniques [26][27][28][29][30] , ensemble-variational methods [31][32][33][34] and neural networks [35][36][37][38][39][40][41][42] .In the context of high-speed, transitional boundary layers, the most relevant work for our purposes is by Buchta et al. 34 .They employed an ensemble-variational (EnVar) technique to deduce the upstream perturbations from downstream wall-pressure measurements on a cone.Their recent investigation was the first-of-itskind to accurately interpret sensor data from within the nonlinear, transitional flow regime.A critical aspect to bear in mind is that for each new set of measurements, the EnVar process must be repeated.Pre-trained neural networks, capable of swift evaluation hold the promise of expediting the solution of these inverse problems.
In this work, we introduce a strategy that exploits deep learning to expedite the interpretation of experimental observations.The measurements are wall-pressure spectra acquired from PCB sensors on a 7-degree cone, at freestream Mach number M = 6.14 , as in Buchta et al. 34 .We estimate the upstream instability waves, quantitatively reproduce the wall-pressure spectra, and discover the full spatio-temporal flow field that led to the measured data.We compare our results with the recent and only known solution, by Buchta et al. 34 , and show that within the time needed for EnVar to solve one data assimilation problem our approach computes about 2000 solutions, with commensurate accuracy.
The flow configuration of interest and the essence of our proposed strategy are presented schematically in Fig. 1.We are provided with experimental measurements from Mach M = 6.14 flow over a 7 • half-angle cone.The measurements are the signals from four PCB pressure sensors ( s 1 , s 2 , s 3 , s 4 ), flush-mounted on the surface of a cone.As depicted in the schematic, we are not provided any information regarding the spatio-temporal details of the flow, and rely solely on the measurements from these four probes.Using Welch's method, we analyze these experimental signals and evaluate their spectral content.The spectral information is then adopted in a data assimilation process to discover the unknown oncoming disturbances.Once discovered, using the upstream conditions as input enables us to simulate and to accurately predict the entire flow which is consistent with the measurements, and thus to have access to the full spatio-temporal fields.Further information regarding the experimental setup can be found in the "Methods" section.
Our data assimilation procedure, named FastDA (Fast assimilation through DeepONet acceleration), adopts a pre-trained DeepONet, in lieu of the computationally demanding compressible Navier-Stokes equations.The DeepONet serves as a computationally efficient meta-model that establishes a connection between the oncoming flow disturbances input and the wall-pressure spectra output, as depicted in the assimilation loop in Fig. 1.Prior to deployment, the DeepONet is efficiently trained using supervised learning, with labeled data obtained from direct numerical simulations of the compressible Navier-Stokes equations, as illustrated in the training loop in the figure.Once trained, the same DeepONet is adopted for many sets of measurements, and we predict the oncoming disturbances associated with each set.
The major advantage of our proposed approach, relative to the previous ensemble-variational method 34 , is the significant reduction in computational expense.Including the initial training cost, we are able to perform approximately 2000 data-assimilation problems in the same wall-clock time that was required to solve just one inverse problem using EnVar.This cost reduction is primarily attributed to the application of the DeepONet as a meta-model for the compressible Navier-Stokes equations.Once trained, the DeepONet is three-orders-ofmagnitude more efficient in terms of wall-clock time for each data assimilation solution.
The utilization of DeepONet requires that we address two notable challenges.The first is to optimize the selection of training data, and specifically to minimize the computational cost for data generation by expensive direct simulations of the Navier-Stokes equations.To address this issue, we efficiently parameterize the space of possible oncoming disturbances, model its statistics using a physics-based approach that exploits the available measurements, and use Latin hypercube sampling as a strategy to cover relevant ranges of parameters.The second challenge has to do with solving the inverse problem for a nonlinear and chaotic system which amplifies errors, in particular when the DeepONet introduces model errors.How do we navigate the tortuous optimization landscape as we search for the optimal oncoming disturbances that reproduce the measurements?For this task, we utilize an enhanced version of downhill simplex optimization, which exhibits improved performance for high-dimensional optimization problems.

Problem formulation and results
Wall-measurements of transitional flows in high-speed flight tests and experiments are valuable.However, they are far from describing the rich and complex dynamics of the flow since the data are often collected by a limited number of localized sensors, and capture only one quantity rather than the entire flow state.The present study uses experimental measurements of wall pressure, which were described in detail by Kennedy et al. 10 .The main question that we address is the following: 'given the wall-pressure spectra at a few probes, can we accurately and efficiently reconstruct the full flow field?'

Measurements and unknown inputs
The time series of wall pressure recordings from the experiments 10 are processed to compute the spectra (Fig. 1) and to form the measurement vector, where p(s N s , f N f ) is the wall-pressure spectral amplitude at sensor s i and frequency f j , N s is the total number of sensors, and N f is the total number of discrete frequencies of interest.The measurements are initially assumed to be statistically stationary, and the entire 80-ms time series is used to compute the spectra using Welch's method.This assumption enables us to compare to previous work, and is subsequently relaxed to demonstrate the efficiency and impact of our method.In all the data assimilation problems presented m ∈ R 76×1 contains spectra with f ∈ [50, 75, ..., 500] kHz, N f = 19 , from four sensors, N s = 4 .Details on the computation of the spectra can be found in the "Methods" section.
In order to reconstruct the flow corresponding to these measurements, it is necessary to estimate the spectral makeup of the oncoming disturbances that amplify within the boundary layer, disturb the flow away from its laminar state, and ultimately lead to the recorded measurements.These disturbances are boundary-layer instability waves, which may be linear at their onset 16 , but quickly become nonlinear, interacting with one another and with the base state itself.The amplitudes of the individual instability waves thus have a profound impact on the flow evolution and, hence, on the sensor signals.Our objective is to identify the amplitudes of these waves, which together form our control vector c ≥ 0 that we aim to optimize such that our computational predictions reproduce the measurements.Only then can we discover the experimental flow field.
We assume that the upstream flow, at the start of the domain of interest, is a linear superposition of a laminar state and fifty two instability waves with frequencies and integer azimuthal wavenumbers (f, k); details are provided in the "Methods" section.Thirteen frequencies are considered f ∈ [50, 75, ..., 350] kHz, and four integer azimuthal wavenumbers k ∈ [0, 20, 40, 60] .The control vector is therefore a fifty-two dimensional, c ∈ R 52×1 , unknown input and searching for its optimal value can be challenging since the system is nonlinear and chaotic.

Fast assimilation
Data assimilation is the procedure that attempts to minimize the difference between available experimental measurements m and their computational prediction m S .The latter depend on the amplitudes of the oncoming disturbances c by the seemingly, and deceptively, simple expression m S = M S (c) , where M S the combination of two steps.First, starting from c the nonlinear, compressible, Navier-Stokes equations q = N (c) predict the flow state q .Direct numerical simulation of this relation is considered the highest fidelity computational model.
Second, measurements are extracted at the locations of the wall-pressure probes and the pressure spectra are evaluated.This step is encapsulated in the observation operator H , where m S = H(q) = H(N (c)) = M S (c) .
The overall computational cost of these two steps is almost entirely due to the first one.When the computational model is direct numerical simulation, each evaluation N (c) of the flow over the cone is a resource-intensive task, which requires 2880 CPU hours (5 wall-clock hours using 576 CPU cores) in the high performance computing cluster with Intel Xeon Gold Cascade Lake 6248R CPUs at Johns Hopkins University.This cost multiplies quickly because data assimilation involves an iterative minimum-seeking procedure where the estimate of c is successively updated until the computational model outputs sufficiently approach the measurements.The minimization is challenging in part due to the nonlinear nature of the compressible Navier-Stokes operator within M S (c) , which requires an iterative optimization procedure.Iterative methods must grapple with the substantial computational demands of numerically solving the compressible Navier-Stokes equations at each iteration.This challenge is particularly pronounced when dealing with experiment-and flight-relevant Reynolds and Mach numbers, as is the case here with the measurements m .In addition, transitional flows are chaotic, and hence lead to an oscillatory cost landscape that is difficult to navigate.
The first, and to our knowledge only, instance in the literature, which attempts the assimilation of physical measurements m from a hypersonic boundary-layer flow, is the recent work by Buchta et al. 34 .They were able to identify the optimal control vector c ⋆ S shown in Fig. 2, using an ensemble variational method (EnVar).EnVar addresses the nonlinear dependence of the cost function on the control parameters by adopting a local quadratic approximation, estimating the local gradient and Hessian using an ensemble of DNS evaluations, and iteratively repeating the gradient-descent procedure.The presence of chaotic dynamics, which leads to an oscillatory landscape of the cost function, is addressed with a number of countermeasures, e.g.attempting multiple initial guesses, inflation of the ensemble covariance matrix, randomizing the search directions.Each iteration of EnVar is computationally costly since it involves an ensemble of direct numerical simulations (DNS).Only the final estimated flow is retained and all the simulations that were performed prior to convergence are discarded.Should new experimental data be collected with a slightly different flow condition, the entire process must be repeated, thus incurring the same high computational cost.These properties render the method primarily suitable for application to a few high-valued experimental configurations rather than an efficient approach that can be widely adopted.
Here, we deploy a fast data assimilation algorithm (FastDA), which takes advantage of the capacity of a deep operator network (DeepONet) to learn nonlinear functionals.The DeepONet model M D (c) efficiently replaces the Navier-Stokes operator M S (c) and significantly reduces computational time to O(0.01) seconds per evalu- ation.FastDA comprises two steps.The first is the DeepONet training, involving minimal dataset generation, supervised learning, and storage of the trained network in memory.The second is the assimilation step, where the optimal control vector is computed through Bayesian estimation using the DeepONet and a gradient-free minimum-seeking technique.
The choice of a gradient-free method mitigates the effect of prediction errors introduced by a meta-model, which pose a challenge to gradient-based approaches due to the high sensitivity of derivatives to fluctuations in the function.Here, our meta-model is a DeepONet due to the desire to adopt an operator network that maps an input function, namely the inflow spectra c , to an output function, namely the wall-pressure spectra m .DeepOnets are specifically designed for this purpose 43,44 .These operator networks can be supplemented with governing equations in the training loss to make up for missing data in training 45 , in which case they are

Minimal dataset and training
The network is shown in Fig. 1, and is pre-trained before deployment as illustrated by the 'DeepONet training' loop, when the arrows are connected to 1 .Training data are generated using expensive DNS.Our goal is to minimize the number of samples to keep computational expenses low, while ensuring that data are both diverse and representative of probable realizations.To achieve this goal, we sample the subspace comprised of the most uncertain directions of the input.Specifically, we map c ∈ R 52×1 to a normal distribution, and reduce the dimension to the leading eigen-directions of the covariance matrix that collectively account for 95% of the total variance.The result is a ten-dimensional subspace, akin to a rank-ten dimensionality reduction by principal component analysis (PCA) 46 .This subspace is parameterized as c = φ(w 0 + W w) , where W ∈ R 52×10 is the matrix of the ten eigen-directions and w ∈ R 10×1 is sampled using LHS from a zero-mean normal distribution with diagonal covariance w (see "Methods").
A total of 68 samples were sufficient to limit the DeepONet's prediction error, ε D , on validation data (not used in training) to 5% (mean-square-error, MSE).Therefore, we perform 68 DNS to gather the necessary measurements M S (c) .To place this cost in perspective, EnVar required 110 evaluations of N for a single data assimilation problem.Moreover, to facilitate efficiency and accuracy of DeepONet training, we included a mapping from the 52-dimensional c to the 10-dimensional w in the input layer of the DeepONet.This mapping can be motivated by the recent work by Kontolati et al. 47 on the superior predictive accuracy of DeepONets when the inputs are mapped to reduced, or latent, spaces.

Bayesian estimation
Once pre-trained, the DeepONet is deployed to predict the measurements M D (c) , which facilitates swiftly search- ing the high-dimensional input space for the optimal amplitudes c that reproduce the experimental measurements m (the ' Assimilation' loop in Fig. 1, when the arrows are connected to 2 ).
We consider uncertain prior knowledge regarding c represented by a likelihood function, P(c) , which is assumed to follow a log-normal distribution logN(c 0 , c ) .This choice implicitly enforces the constraint that the vector of amplitudes c remains non-negative.Additionally, we introduce the transformation c = φ(w) to represent the log-normal distribution using its corresponding normal distribution N(w 0 , w ) and likelihood P(w) .Detailed information about this process is available in the "Methods" section.
Given experimental measurements m , we introduce the conditional likelihood P(M D (w)|w) with normal distribution N(m, m ) .With these likelihoods in place, we leverage Bayesian estimation to compute the optimal c ⋆ = φ(w ⋆ ) that best reconstructs m as the maximizer of the likelihood P(M D (w)|w)P(w) , or equivalently in the context of normal distributions as the minimizer of the negative log-likelihood, The second term in J D (w) can be interpreted as a regularization, and arises naturally in the Bayesian interpreta- tion of the cost function as the negative log-likelihood of Bayes formula for a Gaussian distribution.
The conventional utilization of the Navier-Stokes operator M S (c) in place of the DeepONet model M D (c) results in a data assimilation problem involving the minimization of a different scalar function, J S (w) .This distinction arises because the DeepONet M D (c) acts as a meta-model for M S (c) , rendering J D (w) an approxi- mation, or a noisy reproduction, of J S (w).

Gradient-free minimum seeking
Minimizing J D (w) introduced an additional challenge relative to J S (w) , since computing the gradient of the former amplifies errors in the approximation of the Navier-Stoke operator by the DeepONet.For this reason, we opt for gradient-free techniques which, while less efficient than gradient-based methods, offer greater robustness against false fluctuations in the cost landscape.We employ an advanced version of the Nelder-Mead simplex algorithm 48 , which is also effective at avoiding local minima and has improved convergence in high-dimensional search spaces.This algorithm can thus effectively navigate a cost function landscape that is tortuous, which is typical of chaotic dynamics as in our transitional boundary layer.Convergence is achieved with roughly 3000 evaluations of M D , equivalent to under 2 min on a standard one CPU desktop equipped with an Intel Xeon W-1250 CPU, and yields the inflow amplitudes c ⋆ D = φ(w ⋆ D ) based on The optimal control vector c ⋆ D may deviate from the inflow amplitudes c ⋆ S computed through the minimization of J S (w) due to the discussed differences between J D (w) and J S (w) , arising from the discrepancy between M S and M D .For our strategy to be successful, the deviation between the original solution involving J S (w) and our assimilation with FastDA needs to be accurate up to the validation error of the DeepONet.We now proceed to assess this aspect, and leave further details on dataset creation, neural network training and architecture, and the minimization process of J D (w) for the "Methods" section.

Results: accuracy and cost
Utilizing the same set of experimental measurements m as in the recent study by Buchta et al. 34 , we solve the data assimilation problem by determining the inflow amplitudes through the minimization process described in Eq. ( 3).The corresponding inflow c ⋆ D = φ(w ⋆ D ) is visualized in Fig. 2. Notably, our results agree well with the recent predictions using the EnVar approach, which is based on direct numerical simulations of the Navier-Stokes equations as the forward model.This general agreement is evident when comparing the amplitude distribution across frequencies and azimuthal wave numbers (panels a and b of Fig. 2).
Similar to the EnVar prediction, the majority of the energy is distributed on three-dimensional waves rather than two-dimensional ones.This observation is important because it dispels an entrenched view that transition over cones is predominantly due to planar waves, which are the most exponentially unstable disturbances according to linear stability theory.However, this view is at odds with the measurements, which are most accurately reproduced by three-dimensional waves.One may remark that the importance of three-dimensionality should be expected since transition to turbulence requires it, but the extent to which early three-dimensional waves are relevant to the pre-transitional flow measurements was certainly underestimated prior to application of data assimilation to this problem.
The estimated inflow which is given by c ⋆ D is subsequently adopted in the simulation of the complete threedimensional hypersonic boundary-layer flow, using direct computation of N (c ⋆ D ) .An instantaneous visualization of the field is shown in Fig. 2d, and for comparison the prediction from the EnVar approach is shown in panel c.The three-dimensional iso-surfaces show pressure fluctuations, and agreement between the two figures is evident.The flow exhibits the anticipated streamwise and spanwise fluctuations.In effect, this type of result is the most valuable outcome of solving the inverse problem, since we can now evaluate any flow quantity of interest at full resolution, when the measurements were limited to a few pressure sensors along the wall.
The fidelity of our predictions in reproducing the measurements is important, and is examined in Fig. 3.The figure compares the available experimental data m , their reproduction using DeepONet M D (c ⋆ ) , and their reproduction by Navier-Stokes simulation of the DeepONet-predicted inflow M S (c ⋆ ) .For reference, the figure also shows two limiting cases: the prediction if a linear model is assumed for the forward dynamics and the prediction from EnVar.Notably, the DeepONet effectively captures the primary and secondary peaks in the spectra, which correspond to the primary instability waves and the nonlinear generation of the harmonic.While the primary wave has its origin in linear stability theory, predicting its correct amplitude downstream requires an accurate nonlinear forward model.This point is evident by sensors 3 and 4, where the linear model deviates appreciably from the experiments.In contrast, the DeepONet-dependent results M D (c ⋆ ) and M S (c ⋆ ) accurately reproduce the measurements, within the validation error of the network and as accurate as EnVar.The prediction of the secondary peak due to the nonlinear generation of harmonic frequencies is also within the same level of accuracy as EnVar (note that the figure is in logarithmic scale).To quantify the discrepancies between these curves, we compute the normalized error, The results reveal small differences in the accuracy of EnVar, δ(M S (c ⋆ S )) = 9.98% , compared to the DeepONet outcomes, δ(M D (c ⋆ D )) = 10.67% and δ(M S (c ⋆ D )) = 10.83%.In effect, the reproduction of the measurements is similar in accuracy with either the EnVar or FastDA.This agreement and the accuracy of the predicted state confirm that deploying DeepONet as an approximation of the compressible Navier-Stokes equation is an effective and, importantly, efficient approach to accelerate data assimilation.Additionally, the limitations in achieving a more accurate assimilated flow does not appear to hinge on improving the DeepONet model which in our FastDA yields δ(M S (c ⋆ D )) = 10.83% .In the limit of a perfect DeepONet prediction, we recover the Navier-Stokes operator which in EnVar led to δ(M S (c ⋆ S )) = 9.98% .As such, a more accurate meta-model, although it may reduce the prediction errors, would not lead to substantial improvement in the accuracy of solving the assimilation problem.Achieving a more accurate data-assimilation solution may require a different strategy.Given the dimensionality of the problem, the most effective approach is likely to increase the number of sensors, both in the streamwise and azimuthal directions.Such increase would in fact benefit any assimilation method.Exploring this avenue is, however, beyond the scope of our effort which focuses primarily on expediting the data assimilation process, while preserving accuracy against the most recent advances in the field.Our results demonstrate that we have achieved the accuracy goal, and we now comment on the efficiency of our approach.
In terms of computational cost, our FastDA takes a very small fraction of the wall-clock time requirements of EnVar to solve the same data-assimilation problem.The computational cost of both approaches is reported in Fig. 4, in wall-clock time to emphasize the time to solution independent of the number of computational cores used in the algorithm.Using this metric, FastDA far outperforms EnVar, even if we include the time for the generation of the training data and for training.If we consider the CPU time, the difference is more striking.The wall-clock time of the EnVar algorithm corresponds to 316, 800 CPU hours (110 DNS, 2, CPU hours each), since members of each ensemble are simulated simultaneously and each of them employs many processors.For FastDA, the offline one-time generation of the training data with DNS required the majority of the time, 195, 840 CPU hours (68 DNS, 2, 880 CPU hours each), and an additional two GPU hours were required for training the network.The live data-assimilation itself, which is the search for the optimal control vector that minimizes J D (w) , was performed in 1.73 min on one CPU.In addition, the pre-trained DeepONet can be re- used to solve new data-assimilation problems associated with new experimental measurements m , each within the same 2-min timescale.Specifically, we can estimate the oncoming disturbance spectra c for approximately 2, 000 additional experimental measurements m on a single CPU in the same wall-clock time required for EnVar to solve one data assimilation problem on 576 CPUs (Fig. 4).We now provide examples of the additional data assimilations that we performed.www.nature.com/scientificreports/

FastDA on a cohort of measurements
The time advantage achieved with our FastDA enables efficient solution of a large number of data-assimilation problems, which is impractical using conventional methods.For example, due to the cost of EnVar, Buchta et al 34 used the entire 80 ms time series from the PCB sensors to compute the wall-pressure spectra.Their implicit assumption is that the experimental data is statistically stationary during this duration.Dividing the time series into shorter horizons and assimilating the spectra associated with different sub-intervals was not conceivable given the associated computational cost.Our FastDA approach enables us to consider the many separate cases that result from dividing the 80 ms time series into hundreds of locally stationary segments.We extracted 392 such intervals, applied Welch's method within each of these segments, computed the associated 392 wall-pressure spectra m , and assimilated each of them using FastDA.Our approach thus unlocks the possibility to perform a comprehensive study of inflow variability during the experiment.The 392 wall-pressure spectra m are illustrated in Fig. 5 together with the boundaries of the training data.The orders-of-magnitude variability of the spectra m shown in the figure demonstrates that the inflow c may actually not be stationary, which further underscores the need for a fast data assimilation approach.Moreover, it is noteworthy that, in the frequency range f = [50, 350] kHz, the spectra of the first sensor are perfectly contained within the training dataset.This demonstrates the efficacy of the design of our training dataset with sufficient variance.
The discrepancy δ(M D (c ⋆ D,i )) measured according to Eq. (4) for each of the i = {1, . . ., 392} newly solved data-assimilation problems is presented in Fig. 5b.The predicted control vectors using FastDA reproduce the measurements with accuracy on the order of 10% and with mild variability, similar to the earlier results.The similarity indicates that our method is both effective and robust in performance for measurements within the physics-based design of the training space.The wall-clock time to solution for the complete set of 392 problems is approximately 11 h on a desktop with a NVIDIA RTX A4000 GPU, 32GB of RAM, and an Intel Xeon W-1250 CPU, with solutions computed in series each requiring less than 2 min.
Two solutions from the set of data-assimilation problems are presented in detail in Fig. 6.The inflows c in panels (a) and (b) present a considerable difference in the amplitudes of all frequency and wave-number pairs, especially the most energetic pairs.The range of these energetic waves is the same in both cases.However, since transitional flows are extremely sensitive to variations of the inflow conditions, differences are amplified downstream and the pressure fluctuations at the sensors deviate by orders of magnitude (panels (c) and (d)).A comparison of the respective pressure fields is shown in panels (e) and (f), which highlights the qualitative differences in the two states of the flow, during the same experiment.
These appreciable changes in the flow state are important in the study of hypersonic flows, and are not possible to discern if the entire measurement horizon is assumed to be statistically stationary in a single data assimilation.Our FastDA enables us to assimilate valuable measurements in much finer detail, previously not possible due to computational cost of conventional methods.The present results thus demonstrate the potential role of deep learning in data assimilation for nonlinear, transitional, hypersonic flows.

Conclusions
We have shown how deep learning can be leveraged to drastically accelerate the solution of a data-assimilation problem based on a realistic scenario.The configuration is a Mach 6, transitional, boundary-layer flow over a cone with half-angle 7 • .Only one previous study has successfully performed assimilation of such experimental measurements 34 , using an ensemble-variational method that is computationally costly.Since the majority of the computational cost is associated with performing high-fidelity simulations of the compressible Navier-Stokes equations, we put forward the idea to replace the original high-fidelity model with a computationally efficient deep operator network (DeepONet) previously trained on far fewer high-fidelity simulations.
The proposed approach required us to tackle two key challenges.The first is the design of the training data, which are expensive to compute using high-fidelity simulations, must span a 52 dimensional input space, and ensure that the corresponding output space spans potential measurements of interest.We overcame this challenge by reducing the input space to the most uncertain search directions, using Latin hypercube sampling to experimental measurements; each bar spans a 0.3581% interval of discrepancy error δ (Eq. 4); each column represents the amount of different data assimilation problems within an interval of error values; discrepancy error computed as in Eq. ( 4).
span the inputs, and by adopting a physics-based design of the distribution of samples to span the output space of interest.The second challenge is to solve the data-assimilation minimization problem when the cost-function landscape is affected by the unknown errors of the DeepONet.We facilitated the training of the network by mapping the input layer to the reduced search space, effectively mimicking the so-called L-DeepONet 47 , and demonstrated the capability of an accelerated gradient-free minimization method to find the solution c of the data-assimilation process.
The trained low-cost DeepONet is then deployed to solve the data-assimilation problem.Starting from the experimental wall-pressure spectra at the sensors, our FastDA accurately predicts the amplitudes of the oncoming disturbances, in comparison to the only available data-assimilation result in the literature.The computational cost is significantly reduced, from tens of hours to minutes of wall-clock time (from 65 h to 1.73 min).The time savings allow us to solve many new data-assimilation problems, efficiently.As such, we examined the experimental measurement in fine detail, by dividing the measurement horizon into segments and assimilating each of them independently.The results demonstrated the significant variability in the flow state during the experiment, as well as the robustness and efficiency of our FastDA approach.
Finally, we illustrated the determination of full 3D flow fields from the discovered inflow condition for some select cases by visualizing iso-pressure surfaces.These fields are consistent with the surface pressure data.
Future efforts can integrate uncertainty quantification into the deep-learning framework to provide estimates of prediction uncertainty.This step could be particularly helpful when the training data are contaminated with noise, for example if they are sourced from experiments or lower fidelity computations.Additionally, optimizing sensor placement can enhance the accuracy of data assimilation, especially in situations with limited experimental measurements.Leveraging the DeepONet model can accelerate this process by efficiently mapping inflow spectra to potential sensor placements downstream.

Experimental setup
The experimental measurements used in this work were detailed in the work by Kennedy et al. 10 .The experiment was performed at the Air Force Research Laboratory, in a Mach 6 Ludwieg tube facility [see Kimmel et al. 49 for details].When a fast-acting valve is opened, compressed and heated gas from the driver section of the tube accelerates through the converging-diverging nozzle.The free-stream conditions in the test section are U ∞ = 904 m s −1 , T ∞ = 54 K, and ρ ∞ = 0.0274 kg m −3 , and the associated per-meter Reynolds number is Re ∞ /L = 7.11 × 10 6 m −1 .Each tunnel firing yields two intervals of measurements, each 100 ms, that are assumed to be stationary.The data we used are 80 ms within the second interval.The cone cross-section is circular, its half-angle is 7 • , its length is 414 mm, and its leading-edge has radius r n = 0.508 mm which is considered sharp.At the start of the experiment, the cone is at room temperature, which is essentially maintained during the short-duration test.The experiment was performed at nominally zero angle of attack, and wall-pressure measurements were recorded at 5 MHz using six PCB model 132A pressure sensors.Four of the sensors were placed along a streamwise ray from the leading edge, at positions (s 1 , s 2 , s 3 , s 4 ) = (215, 241, 266, 291, 316) mm.The fourth sensor was flanked in the azimuthal direction by two sensors, each displaced by ±8.5 • .

Wall-pressure spectra of experimental measurements
An 80 ms recording of the wall-pressure was used from each sensor data to compute the wall-pressure spectra m that are the measurements for the data assimilation.The evaluation of the spectra was performed using Welch method 50 , with non-overlapping, Hann-windowed bins of 0.08 ms.These parameters were unchanged throughout the study.
Earlier published papers used the entire 80 ms to compute the spectra 10,34 , effectively assuming statistical stationarity for the entire duration.Here, we have additionally computed the wall-pressure spectra from shorter intervals of the entire 80 ms time-trace, and assimilated these measurements.Specifically, we computed the spectra from intervals [t 0 , t 0 + T l ] , with T l = {4, 6, 8, . . ., 80} ms.The wall-pressure spectra computed from the shortest interval T l = 4 ms contains 50 bins, and that from the longest interval T l = 80 ms contains 1000 bins.For a given T l , we can obtain multiple spectra by shifting t 0 which we have varied with an increment 0.25T l .Using this approach, the total number of wall-pressure spectra that we can extract from the time series is N exp = 392 , each of which can be assimilated to predict the flow within the associated time interval within the experiment.
For each spectra that we adopt as measurements, we quantify the uncertainty of the Welch estimate using the 99% confidence interval of a χ 2 distribution 50 .This confidence interval varies with the number of bins.When only 50 bins are available, the interval is log 10 |p| ± 0.15 , and reduces to log 10 |p| ± 0.03 when there are 1000 bins.These values should be viewed against the magnitude of log 10 |p| , which is order one.

Compressible Navier-Stokes
The flow considered in this work satisfies the compressible Navier-Stokes equations.We adopt as reference scales the free-stream fluid properties (density ρ * ∞ , shear viscosity µ * ∞ , specific heat at constant pressure C * p,∞ ), the free-stream speed of sound c * ∞ , and the length L * = 1 m.Using these reference quantities, the non-dimensional Navier-Stokes equations are, The velocity vector in the i-direction is u i , the pressure is p, the total energy is E = p/(γ − 1) + 0.5ρu i u i , the temperature is T = (γ p/ρ)/(γ − 1) , and γ is the ratio of specific heats.The viscous stress tensor is defined as, where µ b is the bulk viscosity, and the dependence of the shear viscosity on temperature is given by µ = [(γ − 1)T] n .The Reynolds and Prandtl numbers are where k * is the conductivity.In compact operator notation, the governing equations are expressed as q = N (c) , where c are the simulation parameters (e.g.inflow disturbance amplitudes) and q = [ρ, ρu, E] ⊤ is the flow state.Observa- tions are extracted from the state according to m S = H(q) = H(N (c)) = M S (c) , where m S are the wall-pressure spectra from probes located on the wall as in the experimental setup.The high-fidelity simulations that generate the training data solve Eqs.(5-7) downstream of the cone leading-edge shock (see e.  www.nature.com/scientificreports/instability waves n,m c n,m q n,m , where c n,m are the amplitudes and q n,m are the discrete, slow-mode linear instability waves at the frequency and integer azimuthal wavenumber pair ( ω n , k m ).We consider 52 pairs, with frequencies f ∈ [50, 75, . . ., 350] kHz and integer azimuthal wavenumbers k ∈ [0, 20, 40, 60] .The data assimilation searches for the vector c = [. . ., c n,m , . ..] ⊤ of modal amplitudes at the inlet, which accurately reproduces available measurements.

The distribution of the input space
The distribution of c ∈ R 52×1 is logN(c 0 , c ) , with corresponding normal distribution N(w 0 , w ) .Here, c ∈ R 52×52 , w 0 ∈ R 52×1 , and w ∈ R 52×52 .The following element-wise relationships for the means and covari- ances hold, with i = 1, . . ., 52 and j = 1, . . ., 52 .Each element of the vector c can be expressed as which defines the one-to-one mapping c = φ(w).
Buchta et al. 34 generated their initial estimate of c 0 = c ⋆ N using a physics-based approach, which amounts to using the linearized Navier-Stokes operator in the cost function J S (w) and performing a one-shot, least-squares solution.For the nonlinear optimization, they reduced the dimensionality of the search space by focusing only on the subspace where the input has the majority of uncertainty, which corresponds to the space spanned by the ten leading eigenfunctions of their ensemble covariance c .
We search the same sub-space: We select c 0 = c ⋆ N and the original c as the mean and covariance of the distribution in the input, control-vector space, which we map using c = φ(w) , or more specifically Eqs. ( 9)- (10), to w 0 and w .We perform an eigen-decomposition w = W w W ⊤ , and write w = w 0 + W w where W ∈ R 52×10 are the first ten leading eigenvectors of w and w ∈ R 10×1 normally distributed N(0, w ) with zero mean and covariance w ∈ R 10×10 .
Retaining these ten eigen-directions is equivalent to a rank-ten dimensionality reduction through principal component analysis (PCA), or in the continuous form Karhunen-Loève (KLD) or proper orthogonal (POD) decomposition 46 .The eigenvalues associated with these ten leading directions satisfy, which is over 95% accuracy in approximating the covariance matrix.Since the eigenvalues represent the variances along the eigen-directions, this level of accuracy implies that samples of the control vector from this sub-space explore the major uncertainty.
We inflate the rank-ten covariance, w , by a factor of one hundred.This factor was determined using a physics-based approach.The variance of the input distribution is designed such that the training data would likely span an ensemble of available experimental measurements in the output space.This step would normally be computationally costly because statistics of the measurements are nonlinear functions of those of the inputs, due to the nonlinearity of the Navier-Stokes equations.To mitigate this computational cost, we adopt the linearized operator N , and only consider the first sensor location.As previously illustrated in Fig. 3, the linear predictions are accurate at this most upstream sensor, s 1 .We subsequently verified that the choice is adequate, and that the training data satisfactorily samples the output space of the measurements, from all the sensors and in the full nonlinear problem.Finally, we sample N(w 0 , w ) , which has an associated log-normal distribution logN(c 0 , c ).

Latin hypercube sampling
In order to generate the training data, we sample the input distribution using Latin hypercube sampling (LHS).Since LHS is built to sample points in a hypercube, which correspond to sampling multivariate uniform distributions, we apply LHS to sample a uniform distribution U(0, I) with mean 0 ∈ R 10×1 , covariance I ∈ R 10×10 , and each element of the random vector defined on the interval [0, 1].The samples are then mapped to be samples z ∈ R 10×1 of the normal distribution N(0, I) by using the inverse of the cumulative distribution function (CDF).Note, the CDF of independent identically distributed (normal) random variables is the product of the CDF of each random variable, and the CDF is bound in the interval [0, 1].Therefore, the inverse of the CDF can be applied to each element of the random vector separately.The samples of the N(0, I) are then mapped to be samples of N(w 0 , w ) with the linear transformation w = w 0 + W 1/2 w z , which ultimately lead to samples from the desired log-normal logN(c 0 , c ).
To determine the optimal number of samples, we adopt the discrepancy criterion 52 , which is a uniformity criterion that assesses the space filling by a number of samples in a hypercube.The discrepancy can also be interpreted as a measure of the distance between the analytic cumulative distribution function and its samplebased counterpart.Here, we sample 68 different c = φ( w) , corresponding to a discrepancy of approximately (9) c 0,i = e w 0,i + 1 2 w,ii ,

Learning the parameter space
Supervised learning is employed to train the DeepONet.The network models a nonlinear functional that maps the amplitudes of the oncoming disturbance amplitudes c to the wall-pressure spectra M D (c)(x) , where x = (x i , f j ) includes the sensor position and the frequency of interest.The network output is selected to directly predict the wall-spectra at four locations, rather than predicting the entire three-dimensional flow field, with the purpose of solving the data assimilation problem using the available wall-mounted sensor data.The DeepONet architecture is comprised of branch and trunk networks (see Fig. 7), each with 10 layers and 60 neurons per layer.Every neuron but the output layer of the branch network include an ELU activation function.The input to the branch is the vector of oncoming disturbance amplitudes c , and the input to the trunk is the downstream location and frequency x = (x i , f j ) where we wish to predict the measurement.The output to the DeepONet is a dot-product of the outputs of the branch and trunk networks.The final outcome is therefore a weighted superposition of a learned basis, where the weights are the output of the branch (and hence the output layer of the branch does not include activation function) and the basis is the output of the trunk.
Since the equations that govern the nonlinear operator are known, namely the compressible Navier-Stokes equations, a set of labelled data is generated using direct numerical simulations.This dataset is used for supervised learning with a loss function: The first term represents the error between DeepONet predictions and training data, with N D is the cardinality of the dataset.The second term, L WD , is the typical weight decay regularization term 53 .The optimization is performed with the Adam optimizer 54 .The learning rate is 5 × 10 −4 for the first 10,000 epochs and is changed to 5 × 10 −5 for the remainder 5000 epochs.
Figure 7 presents the structure of the dataset, the DeepONet architecture, and the loss function during supervised learning.The data is split into a training set with ≈ 90 % of the input vectors (61 of the 68 c ) and a validation set with the remaining ≈ 10 % (7 different c ).Since the output at the sensor position is M D (c) ∈ R 76×1 , the train- ing and validation sets comprise 4636 and 532 points, respectively.The training points are used for supervised learning, and the validation points are used to compute the generalization error on unseen data during training.The input layer of the DeepONet includes a mapping from the 52-dimensional c to the 10-dimensional w .This step enhances convergence during training, consistent with recent work by Floryan and Graham 55 and Kontolati et al. 47 , where the dimension reduction in the data space was shown to facilitate learning.For example, Kontolati et al. 47  of M D requires O(0.01) seconds.Our DeepONet reaches a mean-square validation error, ε D , of ≈ 5% (Fig. 7), which is also a measure of the average error of the DeepONet for the prediction of unseen in-distribution c.

Searching the parameter space
Replacing M D (c) with M S (c) into Eq. ( 2 gives the cost function J S (w) .DeepONet M D is a meta-model for direct numerical simulations M S , and therefore M D (c) = M S (c) + ε where ε is the error.The associated cost functions for data assimilation are therefore related according to, In effect, J D (w) is a noisy representation of J S (w) .Since derivatives of J D (w) amplify erroneous oscillations in this approximation, use of gradient-based methods to minimize J D (w) is fraught with difficulty.A gradient-free approach much still be robust to noise in the DeepONet predictions, which can lead to an oscillatory cost landscape.For these reasons, we adopt an enhanced Nelder-Mead simplex algorithm 48 .The method is based on the evaluation of J D (w) on the vertices of a simplex, is gradient-free, and is robust to noise.Since the computational cost of the original Nelder-Mead algorithm scales with the dimension of the search space, the version deployed here includes adaptive parameters to maintain computational efficiency, also in a ten-dimensional search space.Therefore, it is possible to navigate the noisy J D (w) to converge to a solution of the minimization problem, The solutions to the data-assimilation problem are accurate up to the validation error of the DeepONet.On average, solution requires O(3000) evaluations of M D , for a total time to solution of O(2) minutes.This average was computed by solving 392 different data-assimilation problems which are discussed in the Results section, each with independent measurements m.

WFigure 1 .
Figure 1.Schematic of the proposed method, FastDA (Fast assimilation through DeepONet acceleration), consisting of two loops.DeepONet training loop (1): train the neural network only once, with available DNS data.Fast assimilation loop (2): search the inflow c which results in the experimental measurements m with a DeepONet as meta-model of the compressible Navier-Stokes equations.

Figure 2 .
Figure 2. Comparison of inflows, the solution to the data assimilation problem, and the corresponding flow fields discovered with direct numerical simulations of the compressible Navier-Stokes equations.Left: FastDA (proposed method); right: EnVar (reference).Top: Inflow spectra c .Bottom: Pressure fluctuations iso-surfaces at −0.1, −0.05, 0.05, 0.1 (blue to red) at matching times.

Figure 3 .
Figure 3. Pressure spectra at the four sensors s 1 − s 4 ; measurements m corresponding to different solutions c ⋆ of the data assimilation problem.Circles: experimental data (99% confidence interval is log 10 |p| ± 0.03 ).Grey line: output m corresponding to c ⋆ N , explicit solution based on the linearized Navier-Stokes operator N .Blue line with dots: output M S (c ⋆ S ) , EnVar solution (reference based solely on direct numerical simulations).Red solid line: output M D (c ⋆ D ) , FastDA solution (proposed method).Red transparent band: ± three times the validation error of the DeepONet, M D (c ⋆ D ) ± 3ε D .Black dashed line: output M S (c ⋆ D ) , FastDA solution simulated with the Navier-Stokes operator N .

Figure 4 .
Figure 4. Wall-clock time to solution for the same DA problem.EnVar: 65 wall-clock hours (316, 800 CPU hours) on a HPC cluster with Intel Xeon Gold Cascade Lake 6248R CPUs (Data Assimilation).FastDA: 5 wall-clock hours (195, 840 CPU hours) on same HPC cluster (generation of training data), 2 GPU hours on one NVIDIA RTX A4000 GPU (DeepONet training), 1.73 CPU minutes on one Intel Xeon W-1250 CPU (Data Assimilation).

Figure 5 .
Figure 5. Left (a): cohort of 392 experimental measurements (black solid lines); boundaries of dataset (red dashed lines).Right (b): performance of fast data assimilation with FastDA on the cohort of 392experimental measurements; each bar spans a 0.3581% interval of discrepancy error δ (Eq. 4); each column represents the amount of different data assimilation problems within an interval of error values; discrepancy error computed as in Eq. (4).
g. Fig 1), and are designed to accurately model the experimental configuration.The inflow Reynolds number based on the post-shock, boundary-layer, edge conditions ( U e = 880 m/s , T e = 65 K , and ρ ∞ = 0.0458 kg/m 3 ) is Re o ≡ ρ e U e L o /µ e = 1,345 using the length scale L o = √ µ e ξ o /ρ e U e .The flow at the inlet to the computational domain is a superposition of an equilibrium (steady) base state q B and (5) ) i , 10% .This number can be put in context by noting that EnVar required a total of 110 evaluations of N to solve a single data-assimilation problem.

2 BranchFigure 7 .
Figure 7. Left: training data obtained via DNS, with 68 inflows chosen with Latin Hypercube Sampling on the prior distribution of c .Center: DeepONet architecture: branch and trunk sub-networks with 10 layers each and 60 neurons per layer.Right: data loss L for the training data (blue) and the validation data (red).