Obtaining predictive low-order models is a central challenge in fluid dynamics. Data-driven frameworks have been widely used to obtain low-order models of aerodynamic systems; yet, resulting models tend to yield predictions that grow unbounded with time. Recently introduced stability-promoting methods can facilitate the identification of bounded models, but tend to require extensive brute-force tuning even in the context of simple academic systems. Here, we show how recent theoretical advances in the long-term boundedness of dynamical systems can be integrated into data-driven modelling frameworks to ensure that resulting models will yield bounded predictions of incompressible flows. Specifically, we propose to solve a particular set of convex semidefinite programming problems to (i) certify whether a system admits a globally attracting bounded set for the chosen modelling parameters, and (ii) compute a model with the optimal (tightest) bound on this globally attracting set. We demonstrate the approach via integration within the sparse identification of nonlinear dynamics (SINDy) modelling framework. Application on two low-order benchmark problems establishes the merits of the approach. We then apply our approach to obtain a low-order (six-mode) model of unsteady separation over a NACA-65(1)-412 aerofoil at ${ extit{Re}}=20\,000$ - a flow that has been notoriously difficult to model using data-driven methods. The resulting model is found to accurately predict the dynamics of unsteady separation, with model predictions remaining bounded indefinitely. We anticipate this work will benefit future efforts in modelling strongly nonlinear flows, especially in settings where physically viable long-term forecasts are paramount.

A spectral element solver is developed for the high-fidelity simulation of the unsteady flow over an aerospike nozzle. The Navier-Stokes solver is a kinetic-energy-preserving, discontinuous Galerkin spectral element method (DGSEM) combined with a hybridization of an entropy viscosity (EV) and a finite-volume subcell element (FVSE) shock-capturing scheme. The diffusive FVSE method is locally called only at locations where the EV method cannot sufficiently smooth the sharp solution gradients that suddenly appear in the supersonic, vortex-dominated jet generated by the aerospike nozzle. Two-dimensional tests of a perfectly expanded and an underexpanded nozzle flow demonstrate that the method is high-order accurate and captures unsteady flow phenomena at supersonic and hypersonic conditions. A resolved three-dimensional simulation at a Reynolds number of 95,000 shows that the solver implicitly models turbulent dissipation at the subgrid scales. To the authors' knowledge, these simulations represent the first DGSEM computations of the resolved, unsteady flow over an aerospike nozzle.

When a chordwise Reynolds number (Re) falls below about 105 the performance of wings and aerodynamic sections become sensitive to viscous phenomena, including boundary layer separation and possible reattachment. Here, detailed measurements of the flow inside the boundary layer on the suction surface are shown for an aspect ratio 3 wing with wall boundaries. The separation lines and recirculation zones are shown on the wing and on the wall junction as Re and angle of incidence (alpha) are varied. There is good agreement on the lowest Re case which has also been computed in direct numerical simulation. Though the flow at midspan may sometimes be described as two-dimensional, at alpha 6 degrees it is unrepresentative of the remainder of the wing, and the influence of the wall is seen in strong spanwise flows aft of the separation line. The geometry of the NACA 65(1)-412 section, used here, promotes a substantial chord length for the development of the recirculating regions behind separation making it apt for their study. However, the phenomena themselves are likely to be found over a wide range of airfoils/wings with moderate thickness at moderate alpha.

The flow over a cambered NACA 65(1)–412 airfoil at $Re\,=\,2\times 10^4$ is described based on a high-order direct numerical simulation. Simulations are run over a range of angles of attack, $\alpha$ , where a number of instabilities in the unsteady, three-dimensional flow field are identified. The balance and competing effects of these instabilities are responsible for significant and abrupt (with respect to $\alpha$ ) changes in flow regime, with measurable consequences in time-averaged, integrated force coefficients, and in the far-wake footprint. At low $\alpha$ , the flow is strongly influenced by vortex roll-up from the pressure side at the trailing edge. The interaction of this large-scale structure with shear and three-dimensional modal instabilities in the separated shear layer and associated wake region on the suction side, explains the transitions and bifurcations of the the flow states as $\alpha$ increases. The transition from a separation at low $\alpha$ to reattachment and establishment of a laminar separation bubble at the trailing edge at critical $\alpha$ is driven by instabilities within the separated shear layer that are absent at lower angles. Instabilities of different wavelengths are then shown to pave the path to turbulence in the near wake.

A data assimilation approach is coined that enables the discovery of forcing functions in Lagrangian, point-particle models from limited measurements of trajectory coordinates. Central to the proposed formulation of this inverse problem is the expression of the forcing function in terms of modal basis functions that are dependent on the relative velocity difference between a known carrier flow and the particle solution weighted with coefficients that are known within confidence intervals. The probability density function of the random forcing coefficients is inferred using a combination of the forward, particle model and its adjoint dynamics, which calculates the gradient of the cost function defined as the distance between the measured and predicted particle locations. To ensure convergence of the gradient-based optimization, multiple measurements may be required. If the measurements are noisy, samples of the forcing model within an assumed Gaussian distribution of the confidence interval of the measurement are computed using a Hamiltonian Monte Carlo method. The method is verified to correctly infer the forcing function of particles traced in the Arnold–Beltrami–Childress flow and a homogeneous isotropic turbulence. The confidence interval of the inferred forcing function with respect to a flow condition is improved if the particle is exposed more frequently to the flow condition. The forcing coefficients adapt the model to flow conditions that are outside of the limited range for which the point-particle models are typically known only empirically or within confidence intervals.

Eulerian-Lagrangian models describe fluid flow and particle dynamics in the Eulerian and Lagrangian frameworks, respectively. In chaotic systems the particle dynamics are stochastic because the suspended particles are subjected to random forces. We use a polynomial chaos expansion (PCE), rather than a postulated constitutive law, to capture structural and parametric uncertainties in the particles' forcing. The stochastic particle dynamics is described by a joint probability density function (PDF) of a particle's position and velocity and random coefficients in the PCE. We deploy the method of distributions to derive a deterministic (Liouville-type) partial-differential equation for this PDF. We reformulate this PDF equation in a Lagrangian form, obtaining PDF flow maps and tracing events and their probability in the phase space. That is accomplished via a new high-order spectral scheme, which traces, marginalizes and computes moments of the high-dimensional joint PDF on conformally mapped hypercubes and comports with high-order carrier-phase solvers. Our approach has lower computational cost than either high-order Eulerian solvers or Monte Carlo methods, is unaffected by the Courant-Friedrichs-Lewy (CFL) stability condition, does not suffer from Gibbs oscillations and does not require (order-reducing) filtering and regularization techniques. These features are demonstrated on several test cases.

A predictive, point-cloud tracer is presented that determines with a quantified uncertainty the Lagrangian motion of a group of point-particles within a finite region. The tracer assumes a random forcing within confidence intervals to account for the empiricism of data-driven force models and stochasticity related to the chaotic nature of the subcloud scale dynamics. It builds on the closed Subgrid Particle-Averaged Reynolds Stress-Equivalent (SPARSE) formulation presented in Domínguez-Vázquez et al. (2023) that assumes a deterministic forcing. SPARSE–R describes the first two moments of particle clouds with moment equations in closed-form, with a theoretical third-order convergence rate with respect to the standard deviations of the cloud variables. The cloud model alleviates computational cost and enhances the convergence rate as compared to Monte Carlo (MC) based point-particle methods. The randomness in the forcing model leads to virtual stresses that correlate random forcing and field fluctuations. These stresses strain and rotate the random cloud as compared to a deterministically forced cloud and thus determine to what extent the random forcing propagates into the confidence intervals of the dispersed solution. In symmetric flows the magnitude of the virtual stress is zero. Tests in analytical carrier fields and in a decaying homogeneous isotropic turbulence flow computed with a discontinuous Galerkin (DG) compressible DNS solver are performed to verify and validate the SPARSE–R method for randomly forced particles.

A method to infer and synthetically extrapolate roughness fields from electron microscope scans of additively manufactured surfaces using an adaptation of Rogallo's synthetic turbulence method [R. S. Rogallo, NASA Technical Memorandum 81315, 1981] based on Fourier modes is presented. The resulting synthetic roughness fields are smooth and are compatible with grid generators in computational fluid dynamics or other numerical simulations. Unlike machine learning methods, which can require over twenty scans of surface roughness for training, the Fourier mode based method can extrapolate homogeneous synthetic roughness fields using a single physical roughness scan to any desired size and range. Five types of synthetic roughness fields are generated using an electron microscope roughness image from literature. A comparison of their spectral energy and two-point correlation spectra show that the synthetic fields closely approximate the roughness structures and spectral energy of the scan.

Langevin (stochastic differential) equations are routinely used to describe particle-laden flows. They predict Gaussian probability density functions (PDFs) of a particle's trajectory and velocity, even though experimentally observed dynamics might be highly non-Gaussian. Our Liouville approach overcomes this dichotomy by replacing the Wiener process in the Langevin models with a (small) set of random variables, whose distributions are tuned to match the observed statistics. This strategy gives rise to an exact (deterministic, first-order, hyperbolic) Liouville equation that describes the evolution of a joint PDF in the augmented phase-space spanned by the random variables and the particle position and velocity. Analytical PDF solutions for canonical models of particle-laden flows serve to establish a relationship between the Langevin and Liouville approaches. Finally, our framework is used to derive a new analytical PDF model for fluidized homogeneous heating systems.

We develop a frame-invariant theory of material spike formation during flow separation over a no-slip boundary in three-dimensional flows with arbitrary time dependence. Based on the exact evolution of the largest principal curvature on near-wall material surfaces, our theory identifies fixed and moving separation. Our approach is effective over short time intervals and admits an instantaneous limit. As a byproduct, we derive explicit formulas for the evolution of the Weingarten map and the principal curvatures of any surface advected by general three-dimensional flows. The material backbone we identify acts first as a precursor and later as the centrepiece of Lagrangian flow separation. We discover previously undetected spiking points and curves where the separation backbones connect to the boundary and provide wall-based analytical formulas for their locations. We illustrate our results on several steady and unsteady flows.

An active control of the flow in a turbine is investigated with the objective of an optimal throughflow mass reduction. An array of jets placed at the suction side of a nozzle vane is intended to enhance the effective camber of the vane, thereby reducing the throat area and corresponding throughflow mass. The effect of the jet array’s location is studied parametrically. The jet array’s configuration is also investigated by spatially perturbing the array with harmonic modes. Steady-state nozzle flow simulations at engine conditions are conducted with ANSYS CFX. The injector jet flows are modeled through source points. The jet flow normal to the vane surface leads to a recirculation region that displaces the inviscid flow outside of the viscous, near-wall region. For a jet location upstream of 37% of the chord length, the recirculation region reattaches downstream on the nozzle vane. For a more downstream array location, the flow separates globally. The throat area and throughflow mass reduction is larger for locations near the throat, but relative to the losses measured through the nozzle loss coefficient, the upstream jet array location is more effective and preferred. A jet array placed at the pressure side has an effect similar to a Gurney flap and modifies the Kutta condition. A jet array that is spatially perturbed with a single mode such that the array’s midspan location is more upstream as compared to the location near the end-walls, homogenizes the flow effect of the jet array, i.e., the recirculation region is more uniform.

A closed and predictive particle cloud tracer method is presented. The tracer builds upon the Subgrid Particle Averaged Reynolds Stress Equivalent (SPARSE) formulation first introduced in [Davis et al., Proceedings of the Royal Society A, 473(2199), 2017] for the tracing of particle clouds. It was later extended to a Cloud-In-Cell (CIC) formulation in [Taverniers et al., Journal of Computational Physics, 390, 2019] using a Gaussian distribution of a cloud's influence over a mesh-based, velocity field solution. SPARSE corrects the cloud's trace to second order by combining a Taylor series expansion of the drag coefficient and Nusselt number correction factors around the mean relative velocity of a cloud of particles with a Reynolds decomposition of the particle equations to obtain a governing system for the first two statistical moments of the cloud's position, velocity and temperature. Here, we close the thus far unclosed SPARSE formulation by determining the velocity field in the vicinity of the mean cloud location using a truncated Taylor series velocity representation and by combining that with averaging. The resulting tracer is predictive. It enables the tracing of a cloud of particles through a single point and so reduces the required degrees of freedom in the accurate tracing of groups of particles. We demonstrate the accuracy and convergence of the method in several one-, two- and three-dimensional test cases.

View Video Presentation: https://doi.org/10.2514/6.2023-1469.vid Viscous effects on a truncated aerospike nozzle flow are systematically investigated by comparing flow fields calculated with a method of characteristics (MoC) for linearized, inviscid flow model, an inviscid Euler flow model, and a Reynolds Averaged Navier-Stokes (RANS) model with slip and no-slip walls. Euler and RANS model simulations are conducted using ANSYS Fluent. The MoC is verified against the Euler model results. Two computational configurations are considered. The first geometry includes a combustion chamber, internal expansion region, and a truncated aerospike geometry while the second configuration is reduced to only the truncated aerospike geometry. Boundary conditions for the latter configuration are specified according to averaged flow profiles extracted from RANS simulations of the former. Expansion waves originating at the sharp geometric throat in the internal expansion region are shown to affect the jet expansion along the aerospike as well as the wake dynamics. Inviscid models compare well with MoC results for the configuration without the nozzle. Viscous model simulations show that the development of a shear layer and boundary layers delay wake closure. The primary vorticity generation mechanism in an inviscid model is baroclinicity in the wake region while in the viscous wake the primary mechanism is shear vorticity generation concentrated in shear layers and boundary layer upstream of the wake region. An increase of the nozzle exit area as compared to the base height is shown to induce wake closure.

Eulerian-Lagrangian models of particle-laden (multiphase) flows describe fluid flow and particle dynamics in the Eulerian and Lagrangian frameworks respectively. Regardless of whether the flow is turbulent or laminar, the particle dynamics is stochastic because the suspended particles are subjected to random forces. We use a polynomial chaos expansion (PCE), rather than a postulated constitutive law, to capture structural and parametric uncertainties in the particles' forcing. The stochastic particle dynamics is described by a joint probability density function (PDF) of a particle's position and velocity and random coefficients in the PCE. We deploy the method of distributions (MoD) to derive a deterministic (Liouville-type) partial-differential equation for this PDF. We reformulate this PDF equation in a Lagrangian form, obtaining PDF flow maps and tracing events and their probability in the phase space. That is accomplished via a new high-order spectral scheme, which traces, marginalizes and computes moments of the high-dimensional joint PDF and comports with high-order carrier-phase solvers. Our approach has lower computational cost than either high-order Eulerian solvers or Monte Carlo methods, is not subjected to a CFL condition, does not suffer from Gibbs oscillations and does not require (order-reducing) filtering and regularization techniques. These features are demonstrated on several test cases.

A discrete analysis of the phase and dissipation errors of an explicit, semi-Lagrangian spectral element method is performed. The semi-Lagrangian method advects the Lagrange interpolant according the Lagrangian form of the transport equations and uses a least-square fit to correct the update for interface constraints of neighbouring elements. By assuming a monomial representation instead of the Lagrange form, a discrete version of the algorithm on a single element is derived. The resulting algebraic system lends itself to both a Modified Equation analysis and an eigenvalue analysis. The Modified Equation analysis, which Taylor expands the stencil at a single space location and time instance, shows that the semi-Lagrangian method is consistent with the PDE form of the transport equation in the limit that the element size goes to zero. The leading order truncation term of the Modified Equation is of the order of the degree of the interpolant which is consistent with numerical tests reported in the literature. The dispersion relations show that the method is negligibly dispersive, as is common for semi-Lagrangian methods. An eigenvalue analysis shows that the semi-Lagrangian method with a nodal Chebyshev interpolant is stable for a Courant-Friedrichs-Lewy condition based on the minimum collocation node spacing within an element that is greater than unity.

A comprehensive and systematic, computational investigation is presented on the effect of wall roughness on the supersonic flow over a circular cylinder with a Reynolds number of 500. Flow simulations are conducted using ANSYS Fluent. Wall roughness is modeled by a perturbation of the cylinder geometry with harmonic modes of varying amplitude and frequency. Validated smooth cylinder flow simulations for a range of Mach and Reynolds numbers with slip and no-slip wall serve as a reference. Roughness is shown to increase the effective diameter of the cylinder and the drag by displacing the outer flow along the peaks of the roughness elements. For lower frequencies, this effect is less pronounced than for higher roughness frequencies. While for smooth cylinders the vorticity is mostly generated by viscous shear forces, for rough cylinder the baroclinic vorticity generation is shown to be dominant and shown to determine the topology of the recirculating region.

The forcing of particles in turbulent environments influences dynamical properties pertinent to many fundamental applications involving particle-flow interactions. Current study explores the determination of forcing for one-way coupled passive particles, under the assumption that the ambient velocity fields are known. When measurements regarding particle locations are available but sparse, direct evaluation of the forcing is intractable. Nevertheless, the forcing for finite-size particles can be determined using adjoint-based data assimilation. This inverse problem is formulated with the framework of optimization, where the cost function is defined as the difference between the measured and predicted particle locations. The gradient of the cost function, with respect to the forcing can be calculated from the adjoint dynamics. When measurements are subject to Gaussian noise, samples within the probability distribution of the forcing can be drawn using Hamiltonian Monte Carlo. The algorithm is tested in the Arnold-Beltrami-Childress flow as well as the homogeneous isotropic turbulence. Results demonstrate that the forcing can only be determined accurately for particle Reynolds number between 1 and 5, where the majority of Reynolds number history along the particle trajectory falls in.

A comprehensive and detailed overview of the flow topology over a cambered NACA 65(1)-412 airfoil at Re = 20,000 is presented for angles of attack ranging from 0 to 10 using high-order direct numerical simulations. It is shown that instabilities bifurcate the flow and cause it to change at a critical angle of attack from laminar separation without reattachment over a laminar separation bubble at the trailing edge to a bubble at the leading edge. The transition of the flow regimes is governed by the Karman vortex shedding of the pressure side boundary layer at the trailing edge, Kelvin-Helmholtz instabilities within the separated shear layer on the suction side, as well as three-dimensional instabilities of elliptic flow within the vortex cores and hyperbolic flow in the shear layer between subsequent Karman vortices. As the suction side shear layer transitions and reattaches, the interaction of the two and three-dimensional instabilities results in three-dimensional tubular structures and large-scale turbulent puffs. The formation and shifting of the laminar separation bubble defines the far-wake topology several chord-lengths behind the airfoil and is accompanied by a sudden increase of the lift force and decrease in the drag that underscores the sensitive nature of low-Reynolds number airfoil aerodynamics. Lift and drag polars are presented for direct numerical simulations, wind tunnel experiments, and simplified numerical procedures where incorrect prediction of the force coefficients is caused by the failure to correctly model the low-pressure region at the trailing edge that is caused by the time-dependent generation of the Karman vortices.