Inspired by biological swimming and flying with distributed sensing, we propose a data-driven approach for load estimation that relies on complex networks. We exploit sparse, real-time pressure inputs, combined with pre-trained transition networks, to estimate aerodynamic loads in unsteady and highly separated flows. The transition networks contain the aerodynamic states of the system as nodes along with the underlying dynamics as links. A weighted average-based (WAB) strategy is proposed and tested on realistic experimental data on the flow around an accelerating elliptical plate at various angles of attack. Aerodynamic loads are then estimated for angles-of-attack cases not included in the training dataset so as to simulate the estimation process. An optimization process is also included to account for the system's temporal dynamics. Performance and limitations of the WAB approach are discussed, showing that transition networks can represent a versatile and effective data-driven tool for real-time signal estimation using sparse and noisy signals (such as surface pressure) in realistic flows.

The instantaneous loads in biological swimming and flying are highly sensitive to environmental perturbations, such as the wakes of other animals, or gusts in the atmosphere, respectively. Despite challenging boundary conditions, animals control the flow over their propulsors (i.e., flippers or wings) with ease and even utilize unsteady flows to their advantage.1,2 Biological sensory systems monitor the flow in real-time by gathering feedback at multiple locations on the propulsors. By combining the sensor input with their experience, animals instantaneously estimate and control their present aerodynamic state (e.g., the aerodynamic loads).3,4 These insights have inspired a series of studies adapting the multi-sensor principle for the control of autonomous swimming and flying vehicles.5–7 In the absence of simple aerodynamic models for three-dimensional (3D) and highly separated flows, data-driven methods have been proposed that utilize sparse pressure data to characterize the instantaneous aerodynamic state on an arbitrary body (e.g., a wing) under a variety of conditions [Figs. 1(a) and 1(b)]. Examples include attached flow, separated two-dimensional (2D) flows, weakly separated, and highly separated three-dimensional flows.6–9 

FIG. 1.

Workflow for the estimation process: (a) complex (bio-inspired) dynamical system with sparse sensors; (b) training datased made up of a collection of time series for different configurations, αi [each collection comprises N time series, ηn(t), from sparse sensors and M time series, ξm, corresponding to the variables to estimate]; (c) phase-space definition of dimension equal to N + M, and construction of the transition networks; (d) collection of a testing dataset of input time series, η̃n; and (e) estimation of the M unknown signals, ξm, using the real-time input in panel (d) and a pre-trained transition networks of panel (c).

FIG. 1.

Workflow for the estimation process: (a) complex (bio-inspired) dynamical system with sparse sensors; (b) training datased made up of a collection of time series for different configurations, αi [each collection comprises N time series, ηn(t), from sparse sensors and M time series, ξm, corresponding to the variables to estimate]; (c) phase-space definition of dimension equal to N + M, and construction of the transition networks; (d) collection of a testing dataset of input time series, η̃n; and (e) estimation of the M unknown signals, ξm, using the real-time input in panel (d) and a pre-trained transition networks of panel (c).

Close modal

Data-driven approaches overcome the need to define explicit aerodynamic models, which are usually problem dependent—hence, less generalizable—and heavily rely on the knowledge of some system's parameters (e.g., the instantaneous angle of attack).10–12 In fact, to accurately perform real-time signal reconstruction and control during unsteady conditions, these parameters are not typically available. Nevertheless, data-driven methods and explicit (mathematical) models could be combined to advance the methodological capabilities through an integrated framework.13 

With the aim to perform (aerodynamic) load estimation, data-driven methods have been shown to be a valid option, although they usually tend to perform well only within a limited range of unsteady boundary conditions.6–9 The need for large training datasets, as well as the availability of very sparse pressure distributions, represent a challenge for data-driven methods attempting to characterize and predict aerodynamic loads.7 Moreover, strong nonlinear effects emerge under realistic conditions (i.e., for highly separated, unsteady flows at high Reynolds numbers), contributing to the challenge in the load-estimation process.

In order to fully account for the effects deriving from such realistic conditions, and to exploit only a sparse set of sensors—hence, without the knowledge of additional parameters such as the instantaneous angle of attack—existing data-driven methods are continuously improved, as well as novel approaches are proposed. Among other techniques, complex networks represent a powerful and versatile tool for time-series analysis14 that have been recently employed to study fluid flows,15 including vortical flows,16,17 turbulent-combustor dynamics,18–20 as well as mixing in wall-bounded turbulence.21,22 In this context, transition networks—thanks to their connection with Markov models14—have been successfully employed for time series reconstruction,23–26 as well as for reduced-order modeling27–29 and control.30 Specifically, Fernex et al.26 have recently shown that cluster-based transition networks can be used as an effective data-driven tool to model complex nonlinear dynamical systems (including turbulence) without any prior knowledge.

Despite the recent progress, cluster-based transition networks are still predominantly used to reconstruct data resulting from accurate, numerically obtained data. In the reconstruction process, the newly generated time series are only expected to be globally similar (i.e., sharing similar statistical features) to the reference time series, which is explicitly included in the training dataset.26,28 In the present study, instead, we apply transition networks to perform signal estimation in real-time and with experimentally obtained data. Here, we generate new signals that, based on sparse (sensor) input, estimate the instantaneous aerodynamic state of an aerodynamic body with good local accuracy.

Hereby, three main challenges arise from real-world (experimental) implementations: (i) limited amount of training data (e.g., range of boundary conditions); (ii) sparse data (limited amount of sensors); and (iii) realistic (noisy) data. To tackle these issues, we present an algorithm that—relying on cluster-based transition networks—is able to perform signal estimation in highly separated experimental flows where sparse sensors are available (Fig. 1). In particular, the network-based strategy can exploit sparse datasets as well as the system's dynamics in the recent past for signal estimation, thus mitigating the need to collect large datasets typical of data-driven approaches.

The methodological steps required to build transition networks are reported in Sec. II, where the new network-based strategy allowing for signal estimation is described (Sec. II C). Our strategy is tested on a simple yet challenging, experimental test case of an accelerated elliptical plate (Sec. III), captured by only two differential pressure sensors. The results of the load estimates are presented and discussed in Sec. IV. Conclusions and future outlook are eventually drawn in Sec. V.

This section describes the signal estimation process, which exploits the features of a transition network built on an experimental training dataset, and a testing dataset of (input) sparse pressure measurements enforcing a constraint to the estimation process of unknown load signals. As shown in Fig. 1, the overall method is characterized by four main steps: (i) the collection of a training dataset [Fig. 1(b)]; (ii) the definition of a phase space from training data and the construction of transition networks [Fig. 1(c)]; (iii) the measurement of (real-time) input data [Fig. 1(d)]; and (iv) the estimation of the load signal [Fig. 1(e)]. We note that, while the procedural steps leading to the construction of transition networks (see Secs. II A and II B) are mainly based on standard practices in the literature,26–28 the novel strategy adopted here to estimate the load signal is provided in Sec. II C.

A training dataset consists of N synchronized time series from sparse sensors (here pressure probes) and M signals corresponding to the variables that have to be estimated (here aerodynamic load time series). In general, the training dataset can comprise multiple collections of synchronized time series [see Fig. 1(b)], where each collection belongs to a different configuration parameter value, αi. The configuration parameter can be, e.g., the Reynolds number, a boundary condition, or a geometrical configuration. In this work, α represents different angles of attack (see Sec. III).

In this study, we consider N =2 pressure signals and M =1 loads (as described in Sec. III), so that a 3D phase space can be obtained. In a phase space, each variable of the training dataset corresponds to a direction, ηn or ξm with n=1,,N and m=1,,M. Figure 2(a) shows a 3D phase space with directions η1, η2, and ξ, and an exemplifying trajectory depicted as a blue dotted arrow. The rationale behind the phase-space construction is to provide a geometrical representation of a multivariate time series, where each set of values {η1(t),η2(t),ξ(t)} at a given time, t, indicates a unique dynamical state of the (flow) system through a unique point in the phase space. By mapping signal data at different times into the phase space, an oriented trajectory can then be formed whose direction is in increasing time.

FIG. 2.

(a) An example of 3D phase space with a trajectory (blue dotted line) and clusters centroids (red dots). (b) A zoomed-in inset in the trajectory, displaying three different Voronoi cells associated with three centroids. The transition network diagram refers to the inset, where green arrows indicate the direction of transition while numbers refer to transition probability values. (c) Schematic for the computation of transition times. The Voronoi cell of S2 is highlighted through gray faces, while the remaining two cells for nodes S1 and S3 are highlight via gray edges.

FIG. 2.

(a) An example of 3D phase space with a trajectory (blue dotted line) and clusters centroids (red dots). (b) A zoomed-in inset in the trajectory, displaying three different Voronoi cells associated with three centroids. The transition network diagram refers to the inset, where green arrows indicate the direction of transition while numbers refer to transition probability values. (c) Schematic for the computation of transition times. The Voronoi cell of S2 is highlighted through gray faces, while the remaining two cells for nodes S1 and S3 are highlight via gray edges.

Close modal

Trajectory points in the phase space are grouped by means of the k-means algorithm.26,28,31,32 Phase-space clustering is usually performed to gain a simplification of the trajectories in the phase space, thus providing a reduced-order representation of the system.26–28 Specifically, the k-means algorithm partitions the phase space into Voronoi cells represented by their cell centroids, Sα, whose entries are the centroid's coordinates in the phase space. The superscript α here indicates that the clustering is applied to the trajectory corresponding to the configuration α of the training dataset. As a result, a cluster-based representation of the training dataset is obtained, where each centroid in the phase space represents a coarse-grained aerodynamic state. Following this idea, in this work each centroid represents a specific set of pressure and load values resulting from a specific (unsteady) condition, which in turn is dependent on the angle of attack and its time history).

In the example of Fig. 2(a), cluster centroids are depicted as red dots, capturing the essential features of the blue-dotted trajectory. The inset in Fig. 2(b) presents three Voronoi cells associated with three centroids, where red straight lines indicate the cell edges. Note that, although the k-means performs an unsupervised clustering, it requires an a priori definition of the number of clusters (i.e., centroids), Ncl, which should be large enough to capture the essential geometrical features of the phase-space trajectories. After Ncl is fixed, the k-means algorithm is applied to each trajectory corresponding to each configuration α.

Transition networks are constructed from clustered trajectories, where cluster centroids are assigned to network nodes. Accordingly, a univocal correspondence exists between Voronoi cells, their centroids, and network nodes, all indicated through the symbol Sα. Network links are weighed on the probability of (temporal) transition between two nodes, thus capturing the temporal dynamics of the complex system. In particular, the probability of transition from node Siα to node Sjα is given by

Pi,jα=N(Siα,Sjα)Nall(Siα),i,j=1,,Ncl,
(1)

where N(Siα,Sjα) is the number of times that trajectories directly transit from node Siα to node Sjα, while Nall(Siα) is the total number of times that trajectories exit from node Siα. In general, Pi,j is not symmetric, i.e., Pi,jPj,iij, so that network links can be illustrated by means of arrows indicating the direction of transition.33 For example, the blue trajectory in the inset of Fig. 2(b) uniquely transits from S1 to S3, so that P1,3=1, as reported in the transition network diagram where links are depicted as green arrows. Moreover, j=1NclPi,j=1 (by definition of probability), and Pi,i=0 (by definition of direct transition between nodes28) for any i=1,,Ncl.

To fully characterize the transition properties of the network, transition times, Ti,jα, are also defined as the average amount of time needed for the transition from a node Siα to a node Sjα.28Figure 2(c) shows a 3D sketch to illustrate the computation of transition times for a given reference cell identified by node S2, where three trajectories (or three intervals of the same trajectory) are illustrated as colored dotted arrows. For node S2, the transition times are T2,3=1.5 and T2,1=3 since the trajectories take, on average, 1.5 and 3 time steps to transit from node S2 to S3 and S1, respectively. As per matrix P, the transition times matrix, T, is also generally asymmetric (Ti,jαTj,iα).

The features of the transition-probability matrix, P, and the transition-time matrix, T, can then be used to generate a new set of signals. Owing to the versatility in the transition network construction, in this work we present a strategy to perform signal estimation referred to as the weighted-average-based (WAB) transition network. The methodology proposed here performs a weighted average among different states in the phase space to create a trajectory of newly generated nodes, as well as an optimization procedure that minimizes the difference between the estimated and input (measured) pressure. In particular, here we assume that input values from N time series are known during the signal generation. For example, input values can originate from a set of N sparse pressure sensors collecting data in real-time, thereby supporting the time-series estimation [Fig. 1(d)].

The time series from the testing dataset are hereafter indicated via ̃ notation (i.e., {η1̃,η2̃,ξ̃}), while the newly generated signals are hereafter indicated via ̂ notation (i.e., {η̂1,η̂2,ξ̂}). We also note that an estimated time vector, t̂, can also be defined since, in general, T entries do not exactly correspond to the time step Δt (a consequence of the clustering operation).

The WAB methodology is described here by highlighting its key features and then sketched in Fig. 3, while procedural details (due to their elaborated nature) are extensively reported in  Appendix A:

  1. The first step of WAB is its initialization [Fig. 3(a)]. A load estimation is obtained at the first time, t0, by employing a nearest-neighbor approach, because a transition approach requires at least two times. Nnn nearest nodes [see cyan filled circles in Fig. 3(a)] are identified for each trajectory (corresponding to each α) with respect to the measured input pressure at t0. Hence, the load value at t0 is evaluated as the distance-based weighted average of the load values coming from each of the Nnn nodes in each trajectory;

  2. The transition probabilities of the networks are then used to continue estimating the load signal at a generic th>t0 [Fig. 3(b)]. In particular, Nnn nearest nodes, with respect to the pressure input and the previously estimated load, {η̃1,η̃2,ξ̂}(th1), are first identified for each trajectory. By so doing, a set of weights, wi,distα, can be defined that are inversely proportional to the distance between the Nnn closest nodes and {η̃1,η̃2,ξ̂}(th1) [see cyan dashed lines in Fig. 3(b)]. The transition matrix is then exploited to identify the transition target nodes of each of the Nnn (source) nodes, following the criterion of maximum transition probability [see magenta circles and arrows in Fig. 3(b)];

  3. To not only account for the present system state, but also for the recent past of the pressure input, we define a second set of α-specific weights, woptα [see Fig. 3(c)]. The weights, woptα, are generated through an optimization process that minimizes the error between the estimated pressure η̂ and the input (measured) pressure η̃ over the recent time period Δt. As such, the configurations, α, with a similar pressure history as the input data, will have a higher impact on the load estimate; and

  4. The load ξ̂(th) is eventually estimated as a weighted average of the loads from the identified target nodes [see red-filled dot in Fig. 3(d)],
    ξ̂(th)=αwoptαiwi,distαξ(Siα)αwoptαiwi,distα.
    (2)
FIG. 3.

Schematic of the WAB method for a 3D phase space, including three trained trajectories depicted as orange, green, and blue lines. (a) Initial estimate of the load by weighted averaging using the 2D distance dη. (b) Identification of transition targets. (c) Weight optimization (optional) based on the recent pressure history. (d) New estimation using transition targets (highlighted in magenta) for the weighted average. Here, wiα is calculated via Eq. (A3) using the weights wdist,iα and woptα, obtained in panels (b) and (c), respectively. As in this example the history of α3 is significantly different from the measured data η̂, the α-dependent weight is small (woptα31). Thus, the impact of the corresponding transition network [α3, transparent in panel (d)] on the estimated state is negligible.

FIG. 3.

Schematic of the WAB method for a 3D phase space, including three trained trajectories depicted as orange, green, and blue lines. (a) Initial estimate of the load by weighted averaging using the 2D distance dη. (b) Identification of transition targets. (c) Weight optimization (optional) based on the recent pressure history. (d) New estimation using transition targets (highlighted in magenta) for the weighted average. Here, wiα is calculated via Eq. (A3) using the weights wdist,iα and woptα, obtained in panels (b) and (c), respectively. As in this example the history of α3 is significantly different from the measured data η̂, the α-dependent weight is small (woptα31). Thus, the impact of the corresponding transition network [α3, transparent in panel (d)] on the estimated state is negligible.

Close modal

As such, the weighted average accounts for both the phase-space distances (via wi,distα) and the temporal dynamics of each trajectory (via woptα).

To conclude, an estimated time, t̂h, is also computed by using transition times T instead of load values in the weighted average.

As a test case for the transition-network frameworks presented in Sec. II, realistic experimental flow data were captured in a highly separated and unsteady flow at a high Reynolds number. In particular, the flow around an accelerating elliptical plate was characterized via pressure and load measurements, and the same experimental setup was used to obtain both training and testing datasets, as described in Sec. II.

The experiments were performed in a fully enclosed, water-filled (viscosity ν), 15-m-long towing-tank facility with 1m×1m cross section [Fig. 4(a)]. The model consisted of an elliptical plate [Fig. 4(a)], with principal axes b=0.3 m and w=0.15 m and a cross-sectional area A=πbw. The model was connected to the traverse above the towing tank via a horizontal sting with diameter 0.08b and length 2b, and a vertical symmetric profile of thickness 0.08b. The plate was towed from rest with the plate velocity Up accelerating at a rate of 0.4 m/s2 until hitting its final velocity U=1 m/s, resulting in a terminal Reynolds number of Re=Ub/ν=194000. The plate velocity (Up=U) was then kept constant over a distance of 40Dh before it was decelerated to rest, where Dh is the hydraulic diameter. The same kinematics were tested for the plate mounted at various angles of attack α [as defined in Fig. 4(b)], in the range of 45°α130°.

FIG. 4.

(a) Towing-tank facility; (b) front and side views of the elliptical plate model with force-torque sensor, an angle adapter to determine α and two differential pressure sensors. Normalized plots of (c) the differential pressures, and (d) the plate-normal force as a function of the normalized time t*=Ut/Dh and for various angles of attack. The acceleration stage (red), steady stage (blue), and the deceleration stage (green) are highlighted. (e) 3D phase space built on CN, Cp1, and Cp2. Phase-averaged data (colored), and single-run data (gray) are both shown. For α=45° (red trajectory), Ncl = 60 nodes, Sα (corresponding to cluster centroids), are also displayed.

FIG. 4.

(a) Towing-tank facility; (b) front and side views of the elliptical plate model with force-torque sensor, an angle adapter to determine α and two differential pressure sensors. Normalized plots of (c) the differential pressures, and (d) the plate-normal force as a function of the normalized time t*=Ut/Dh and for various angles of attack. The acceleration stage (red), steady stage (blue), and the deceleration stage (green) are highlighted. (e) 3D phase space built on CN, Cp1, and Cp2. Phase-averaged data (colored), and single-run data (gray) are both shown. For α=45° (red trajectory), Ncl = 60 nodes, Sα (corresponding to cluster centroids), are also displayed.

Close modal

Two Omega differential pressure transducers captured the instantaneous differential pressure Δp between the two sides of the plate. Figure 4(b) shows the positions of the two pressure transducers at y/b=0.5 (Δp1) and y/b=0.75 (Δp2), respectively. The pressure sensors measure a range of ±6895 Pa and have a response time of 103 s, and an accuracy of ±0.25% the full-scale best fit straight line (FS BFSL) with hysteresis and repeatability of 0.2% FS. In order to measure forces and moments on the plate, an ATI Nano 25 six-axis force-torque sensor was mounted between the plate and the horizontal sting. The transducer has a resolution of 0.125 N. Pressures and forces were recorded at a sampling frequency of 1000 Hz.

Figures 4(c) and 4(d) present the temporal evolution of the normalized pressures Cpi=2Δpi/(ρU2), and the plate-normal load CN=2FN/(ρAU2), for the six angles of attack α={45°,60°,75°,90°,110°,130°}. Here, i =1, 2 refers to sensor position at y/b=0.5 and y/b=0.75, respectively. Inertial effects due to the plate's mass are within 1% of the average loads and, therefore, they are not accounted for. The pressures and loads were phase-averaged over 10 runs and temporally filtered with a Savitzky–Golay filter.34Figure 4(e) visualizes the same data of Figs. 4(c) and 4(d) in 3D phase space, whose directions are η1=Cp1,η2=Cp2, and ξ=CN. Single-run data are also illustrated in Fig. 4(e) as gray trajectories, in addition to the phase-averaged data (colored).

For smaller α values, the spacing between different trajectories is notably visible [e.g., the red and yellow trajectories in Fig. 4(e)], while similar pressures are observed for α>75°, making it difficult in this phase space to distinguish between the trajectories of different α. As such, the present dataset is particularly challenging with regard to accurate load estimates using transition networks. Specifically, ambiguous states are likely to appear, namely, points in the phase space with similar pressure values but different loads.

The transition networks were built using a training dataset made up of the pressure data (η1=Cp1,η2=Cp2) and the plate-normal load (ξ=CN), as well as following the description provided in Sec. II. The order of the experimental data was then reduced by clustering the phase-space trajectories for each α into Ncl = 300 centroids, Sα. As a representative example, centroids are shown in Fig. 4(e) for the (phase-averaged) trajectory corresponding to the configuration α=45°. In general, small values of Ncl serve to reduce the computational effort of the method. However, if Ncl is too small, the dynamics of a trajectory in the phase space cannot be properly resolved, thus leading to higher estimation errors. In the present study, Ncl = 300 (with 12 300 time-series instants, i.e., trajectory points) provides a good balance between estimation accuracy of the load ĈN and computational effort. A parametric analysis on the effects of Ncl on the results is provided in  Appendix B.

Once the transition networks are established, a real-time estimate of the plate-normal load ξ̂=ĈN can be obtained by utilizing the pressure sensors' (real-time) input (η̃1=C̃p1 and η̃2=C̃p2) in combination with the WAB procedure (Sec. II C). Specifically, the number of nearest-neighbors Nnn [cyan nodes in Figs. 3(a) and 3(b)] was set equal to 10. Although Nnn is usually set to 3 or 4,26Nnn = 10 leads to smoother load estimates without substantial changes in the overall results.

Furthermore, pressure-error minimization was performed over a temporal window Δt*6 (i.e., a traveled distance at most equal to six times the plate hydraulic diameter). The value of Δt* corresponds to half of the acceleration time t*12 [see Fig. 4(d)] and allows us to sufficiently capture the temporal features of the acceleration and deceleration phases. Larger Δt* values do not lead to substantial changes in the results.

In this section, we present and discuss the results stemming from the application of the WAB approach when the transition networks are used to estimate loads for omitted flow configurations, i.e., time series corresponding to α values that were not available in the training data. To mimic realistic estimation conditions, a randomly selected single run [gray lines in Fig. 4(e)] is used as a testing time series, rather than the phase-averaged signals [colored trajectories in Fig. 4(e)]. In this way, we account for single-run noise, as phase-averaged signals are less noisy.

Figure 5 shows the resulting load estimates ĈN compared to the measured loads C̃N, for different omitted α values (reported in each panel's title). To highlight the impact of including the system dynamics in the estimation process, we compare the results with and without utilizing the α-specific weights woptα obtained via the pressure optimization strategy presented in Fig. 3(c) (dark blue and black lines, respectively). In general, the WAB approach is able to reproduce well the shape of the measured loads C̃N.

FIG. 5.

Load estimates for various omitted α values reported in each panel's title: (a) α=45°; (b) α=90°; (c) α=110°; and (d) α=130°. The input pressure data were randomly selected among all single runs of the respective α. Estimated loads, ĈN, with or without pressure-error minimization, are depicted as dark blue and black lines, respectively. As a reference, the measured loads C̃N are shown using the same colors as in Fig. 4(e).

FIG. 5.

Load estimates for various omitted α values reported in each panel's title: (a) α=45°; (b) α=90°; (c) α=110°; and (d) α=130°. The input pressure data were randomly selected among all single runs of the respective α. Estimated loads, ĈN, with or without pressure-error minimization, are depicted as dark blue and black lines, respectively. As a reference, the measured loads C̃N are shown using the same colors as in Fig. 4(e).

Close modal

To quantify the estimation performance in more detail, Fig. 6 presents the normalized error, E=Eabs/C̃N, of the estimates ĈN with respect to the test data C̃N, where Eabs=|ĈNC̃N| is the absolute error while indicates the time average. The accuracy of the estimate varies depending on the flow stage [i.e., acceleration, steady-state, deceleration; see Figs. 4(c) and 4(d)], the estimation strategy (with or without pressure optimization), as well as on the omitted α. The magnitude of E remains relatively small (less than 20%) throughout the whole estimation process, even if optimization is not performed (see also Fig. 7 for a comprehensive assessment).

FIG. 6.

Normalized absolute errors, E, for the load estimates shown in Fig. 5: (a) α=45°; (b) α=90°; (c) α=110°; and (d) α=130°.

FIG. 6.

Normalized absolute errors, E, for the load estimates shown in Fig. 5: (a) α=45°; (b) α=90°; (c) α=110°; and (d) α=130°.

Close modal
FIG. 7.

Effect on the number nodes (i.e., cluster centroids) on the estimation performance. The legend indicates the estimated case: (a) α=45°,60°,75°; and (b) α=90°,110°,130°. The black vertical line corresponds to Ncl = 300 as in Sec. IV.

FIG. 7.

Effect on the number nodes (i.e., cluster centroids) on the estimation performance. The legend indicates the estimated case: (a) α=45°,60°,75°; and (b) α=90°,110°,130°. The black vertical line corresponds to Ncl = 300 as in Sec. IV.

Close modal

Note that, in contrast to previous studies on signal reconstruction,26,28 we use the (normalized) absolute error E instead of statistical quantities, such as autocorrelation to assess the estimation quality. A discussion on the relevance of statistical similarity in the context of signal estimation is provided in  Appendix C.

Estimating the aerodynamic state in real-world applications, such as our accelerated flat plate, imposes several challenges as outlined in Sec. I: (i) a limited amount of training data; (ii) a limited amount of sensors; and (iii) realistic (noisy) data. In the following, we use the present results to discuss how the network-based estimation algorithm presented in Sec. II C tackles those challenges.

Training a data-driven algorithm in an experimental setting comes with significant effort. In fact, experimental campaigns are often time intensive and involve costly facilities. To potentially reduce the required amount of training data needed for accurate load estimates, the WAB approach estimates unknown dynamics by combining information from different configurations, α, via a weighted average. A similar approach was proposed by Fernex et al.,26 who evaluated the weighted average between two configurations that were identified manually and a priori. In contrast, the WAB approach used here takes all available configurations, α, of the training dataset into account, since one does not know a priori which α value(s) are suitable for the present estimation.

By omitting a configuration α from the training data, and then estimating the aerodynamic loads of the omitted α, we can assess the capabilities of the WAB approach to estimate an unknown signal. As shown in Fig. 5, the WAB approach successfully estimates loads from α that were omitted in the training data. The estimates are particularly accurate for α=90° and α=110° [see Figs. 5(b) and 5(c) and Figs. 6(b) and 6(c)]. However, the performance of the WAB approach deteriorates if α=130° or α=45° are omitted in the training dataset and then estimated [see panels (a) and (d), respectively, in Figs. 5 and 6].

Taking α=130° as a representative case, this behavior can be explained by the fact that the trajectory for α=130° is not fully surrounded by other trajectories in the phase space [see Fig. 4(e)], but is only close to the trajectory for α=110° [orange line in Fig. 4(e)]. From the point of view of WAB, the weighted average to estimate the load for α=130° is performed using load data that are always lower than the expected C̃N, so that the average is unavoidably driven by lower CN values, thus leading to higher estimation errors. In contrast, the trajectories of α=90° and α=110° are surrounded by states of various α [as shown in the phase-space representation of Fig. 4(e)], so that the WAB strategy can better interpolate to the estimated force.

This limitation of the WAB approach is a common feature of interpolation-based techniques, and can be overcome by properly collecting training data [Fig. 1(b)] so that all expected peripheral boundary conditions are accounted for. For example, in the present experimental dataset, the training should contain the maximum and minimum α that can be experienced by the system.

Physically implementing a dense network of pressure sensors on an aerodynamic body requires high design and production costs. However, when the amount of sensors is significantly reduced, ambiguous states will likely occur (i.e., same pressure input values, but different loads; see Sec. III). To accurately estimate the aerodynamic loads with sparse data, the present WAB approach concurrently relies on the instantaneous sensor data and the recent history of the system. Namely, the information from the previous state estimate at th1 and the instantaneous sensor input are combined to determine the node-specific weights wi,distα [Fig. 3(b)]. In addition, the recent pressure history is used to find a set of α-specific weights woptα that minimize the error of the pressure estimate within the recent past [window Δt*, Fig. 3(c)]. As such, both weights, wi,distα and woptα, contribute to the WAB's performance in systems with sparse sensors.

The positive effects of using woptα, obtained by the pressure-error minimization, is apparent when comparing the estimated loads with and without optimization in Figs. 5 and 6. The estimation performance is consistently better when woptα is used. This is particularly evident for α=90° and α=110°, in which optimized WAB [black lines in Figs. 5(b) and 5(c)] is able to capture the initial load bump occurring at the end of the acceleration phase (t*12). This local increase is particularly challenging to be estimated as a result of the strong unsteadiness and nonlinearity in the system, thus highlighting the potential of optimized WAB in performing load estimation effectively.

As shown in Sec. IV B, including the system dynamics in the estimation process can lead to better cluster-based modeling and estimation performance. While our WAB approach relies on a pressure-error optimization, alternative approaches were suggested to account for system dynamics during the state estimation. For instance, Nair et al.30 accounted for the system dynamics by adding the temporal derivative of their (numerical) input data as an additional axis of the phase space, thus obtaining a better aerodynamic state identification. However, under realistic (experimental) conditions, training and input data display noise levels as a result of several (systematic or randomly appearing) factors affecting the measurements. In particular, the noise level of experimental pressure data obtained in separated, high-Reynolds number flows is typically very high, leading to inaccurate evaluations of temporal gradients. Therefore, in such cases, temporal gradients do not generally represent a suitable choice to be included in the phase space. Accordingly, our WAB approach has been conceived to rely only on absolute values of the sensor input and not on their temporal gradients.

Furthermore, to mitigate the impact of experimental noise on estimation, we used phase-averaged data as a training dataset [see Figs. 4(c)–4(e)]. Nevertheless, the estimation capabilities were tested on (randomly selected) single-run time series, which differ from the phase-averaged signals, and provide additional challenges to the estimation accuracy. In spite of these challenges, the WAB approach still proved to be accurate even in the presence of noisy input data.

In this study, we extend the application of transition networks from signal reconstruction to load estimation with real-time input. In particular, we generate new signals that were not included in the training dataset, utilizing the input of N =2 sparse sensors. A weighted average-based (WAB) network strategy is proposed and tested under realistic conditions of unsteady flows. In particular, the network-based approach is tested on an experimental dataset with pressure and load data from an accelerating elliptical plate at various angles of attack. The WAB strategy exploits the features of transition networks (which comprise the definition of a phase space and a clustering algorithm) and a real-time input of sparse pressure signals. Furthermore, an optimization process that minimizes the difference between estimated and measured (input) pressure is implemented.

The potential and limitations of the WAB approach are discussed for estimates corresponding to different (omitted) angles of attack. The results indicate that transition networks can estimate configurations that were unknown during the training stage, with global estimation errors below 20%, and for some cases even below 10%. Moreover, the pressure optimization approach is able to further refine the estimation outcomes by also capturing characteristic local behaviors of the unsteady load signals, e.g., after the acceleration phases. Therefore, our WAB approach proves to be a robust and accurate tool for signal estimation, even in the presence of sparse input data and with limited training data.

While the current approach represents a first effort to employ transition networks to estimate aerodynamic loads in unsteady and high-Reynolds number flows, several methodological advancements can be implemented to potentially enhance the capabilities of the WAB approach. In fact, transition networks do not represent a black box tool between input and output variables but provide a versatile framework that can be easily modified to account for advanced information on the flow system.

In this regard, the implementation of the system dynamics through pressure-error minimization represents a paradigmatic example. Additional physical insights can be included in the network model by expanding the size of the phase space, where additional axes could represent other measured variables of the system. In this vein, future efforts aim to incorporate external flow measurements (e.g., velocity or vorticity fields) in the transition network model. On this note, the outcomes from simple models could also be incorporated in the algorithms, either as additional axes of the phase space or as additional input data to constrain the estimation process.

With the aim to account for noise in experimental data, different clustering strategies could also be applied (such as fuzzy algorithms), and the probabilistic nature of transition networks could be further exploited, e.g., implementing Bayesian statistics.35 Furthermore, different interpolation schemes can be used, where the weighted average can be performed on a limited set of configurations chosen via an optimization routine.

In conclusion, transition networks show great potential for real-time estimation of unknown variables in fluid dynamics problems, even under challenging flow conditions and sparse training datasets. Therefore, we believe the proposed network-based methodology, owing to its versatility, can be a promising tool for the real-time estimation of realistic flows even when limited by sparse data.

D.E.R. acknowledges support from the Air Force Office of Scientific Research (AFOSR) under Grant No. FA9550-20-1-0086, monitored by Dr. Gregg Abate.

The authors have no conflicts to disclose.

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

The details of the WAB methodology are here reported. We recall that the time series from the testing dataset are indicated via ̃ notation (i.e., {η1̃,η2̃,ξ̃}), while the newly generated signals are hereafter indicated via ̂ notation (i.e., {η̂1,η̂2,ξ̂}). The WAB approach comprises the following four main steps (see Fig. 3):

  1. At the first time, t0, a nearest-neighbor approach is used to estimate ξ̂(t0) because a transition approach requires at least two times. First, a set, Snnα, of Nnn nodes is selected for each configuration α. Specifically, each set of nodes comprises the closest Nnn nodes to the measured pressures of the testing dataset, namely, η̃(t0)={η̃1(t0),η̃2(t0)}. For example, in Fig. 3(a), Nnn = 2 and the closest nodes to the two configurations α1 and α2 are identified by dashed cyan lines, which highlight the 2D distances dη=||η(Snnα)η̃(t0)||2. The load value, ξ̂(t0), is then evaluated as the weighted average of the ξ values of each node in Snnα, namely,
    ξ̂(t0)=αiwiαξ(Siα)αiwiα,
    (A1)

    where SiαSnnα, while wiα=1/dη is a set of distance-dependent weights. A newly estimated node Ŝ(t0)={η1̃(t0),η2̃(t0),ξ̂(t0)} is then obtained, which is illustrated as a filled red dot in Fig. 3(a);

  2. The transition probabilities of the networks are then used to continue estimating the load signal. At a generic time th, the Nnn closest nodes, Snnα, are first identified. In particular, the nodes in Snnα are selected to minimize the Euclidean distance, dη,ξ, between the point {η̃1,η̃2,ξ̂}(th1) [see filled red circle in Fig. 3(b)] and all the nodes belonging to each trajectory. For any node in Snnα, the transition matrix is exploited to identify the transition target nodes following the criterion of maximum transition probability. Target nodes are highlighted by magenta circles in each trajectory of Fig. 3(b), while transitions are shown via magenta arrows.

    Similar to the initialization at t0 [Fig. 3(a)], a set of weights wdist,iα=1/dη,ξ can be defined, implying that closer Snnα nodes will have a higher impact (i.e., a higher weight) on the estimation of the next state Ŝ(th). We note that during initialization, a 2D distance (dη) was used. For th>t0, a load estimate ξ̂(th) exists that can be exploited to calculate a 3D distance (dη,ξ), as shown by cyan dashed lines in Fig. 3(b). In contrast to a 2D distance, the 3D distance provides more robustness against ambiguous estimations. These arise when the η1 and η2 values of the nodes Snnα are similar, but their ξ values are significantly different;

  3. Although the identified nodes Snnα are close to Ŝ(th1) in the phase space, this does not necessarily imply that Snnα nodes belong to a trajectory (i.e., an α case) displaying a similar temporal dynamics as the input data η̃. Therefore, to identify the α cases that best capture the dynamics of the input data, we rely on the recent past of our input data η̃. In particular, we define a new set of α-specific coefficients, woptα, which are generated through a pressure-error minimization strategy [Fig. 3(c)]. The input data, η̃, in a temporal window [thδh,th] are used as reference values and compared with the estimated pressure values η̂, computed as
    η̂=αiwiαη(Siα)αiwiα,
    (A2)
    where
    wiα=wi,distα·woptα.
    (A3)
    A minimization problem is then solved which aims to find the set of weights woptα that minimizes the maximum absolute difference between input and estimated pressure [dashed and solid red lines in Fig. 3(c)] over the chosen temporal window, [thδh,th], namely,
    argminwoptα[maxt[thδh,th][|η̂1,2(t;woptα)η̃1,2(t)|]].
    (A4)

    We note that woptα only depends on α, thus providing a measure of the reliability of each trajectory (corresponding to configurations α) to fulfill the constrain on estimated pressure coming from input (testing) data. For example, the blue trajectory corresponding to α3 in Fig. 3(c) is much less reliable than the remaining two trajectories (for α1 and α2) in providing a good estimation for pressure, thus leading to woptα31. If optimization is not performed, woptα=1 for any α configuration and wiαwi,distα; and

  4. Finally, a weighted average of the ξ values of the target nodes is computed to estimate the load value ξ̂ at time th, thus obtaining the point Ŝ(th) [Fig. 3(d)]. Equation (A1) is still exploited to get Ŝ(th), but the weights wiα=wi,distα·woptα [Eq. (A3)] are defined as the product of the distance-related weight wi,distα=1/dη,ξ [cyan dashed lines in Fig. 3(b)] and the α-specific novel coefficient woptα [Fig. 3(c)] which accounts for the system dynamics.

To conclude the procedure, the estimated time t̂h is also computed as the sum of the previous estimated time, t̂h1, and a weighted-average transition time from Eq. (A1) where transition times T are used instead of ξ.

This Appendix describes the effects the Ncl parameter on the estimation performances of the WAB transition network strategy. We recall that Ncl indicates the number of nodes in the network, which correspond to the centroids of the Voronoi cells obtained from the k-means clustering.

Figure 7 shows the average performance of WAB when Ncl is varied, either with or without pressure-error minimization. Here, a global error is computed as E=Eabs/C̃N. In general, E increases toward small Ncl values because Ncl becomes comparable with the number of nodes used to perform the weighted average, Nnn = 10. As discussed in Sec. IV, the WAB method performs well when intermediate configurations have to be estimated, which is confirmed in Fig. 7 for α={75°,90°,110°} (see blue, green, and orange lines, respectively). In particular, it is evident as the pressure-error minimization (filled-dotted lines) always reduces the overall estimation error for any omitted α. Finally, we highlight that E remains quite constant for Ncl > 300, with values below 10% for the intermediate cases and below 20% for external cases (i.e., α=45°,130°), thus justifying the choice of Ncl = 300 in Sec. IV.

Data-driven approaches have often been used to reconstruct data, so that the newly generated reconstructed time series are expected to share very similar statistical features with respect to the corresponding signal included in the training dataset.26,28 In contrast, in the present study, we aimed to estimate unknown signals, i.e., not included in the training dataset. In general, the newly estimated signals could display very similar statistical features with respect to the expected time-series, while still containing considerable local errors. In other words, although an estimation process could produce a globally (statistically) similar time-series with respect to the expected signal, local errors can be non-negligible.

In our work, this could be a consequence of the fact that each single run will always differ locally from other runs or from phase-average signals that are included in the training dataset. A representative example is provided in Fig. 8: we show as a black line the autocorrelogram of the estimated load ĈN (without performing optimization) for the omitted case α=130°, as a function of the temporal lag Δτ*. For comparison, the autocorrelogram of the measured (reference) load C̃N is also reported as a cyan line. Autocorrelation is chosen here in analogy with previous studies assessing the capabilities of transition networks.26,28 While estimation errors can be locally non-negligible [as illustrated in panel (d) of Figs. 5 and 6], the difference between the autocorrelogram for the estimated and measured loads is instead very small.

FIG. 8.

Autocorrelogram for the measured (cyan line) and estimated (black line) loads for omitted case α=130°. Estimated load is obtained without optimization. The root mean square error between autocorrelations is 0.022.

FIG. 8.

Autocorrelogram for the measured (cyan line) and estimated (black line) loads for omitted case α=130°. Estimated load is obtained without optimization. The root mean square error between autocorrelations is 0.022.

Close modal

Therefore, while statistical tools like the autocorrelation could be effectively used to assess reconstruction performances, they might not always provide a reliable measure of the estimation performances, especially in the context of unsteady load estimation.

1.
J.
Liao
,
D.
Beal
,
G.
Lauder
, and
M.
Triantafyllou
, “
Fish exploiting vortices decrease muscle activity
,”
Science
302
,
1566
1569
(
2003
).
2.
S.
Portugal
,
T.
Hubel
,
J.
Fritz
,
S.
Heese
,
D.
Trobe
,
B.
Voelkl
,
S.
Hailes
,
A. M.
Wilson
, and
J. R.
Usherwood
, “
Upwash exploitation and downwash avoidance by flap phasing in ibis formation flight
,”
Nature
505
,
399
402
(
2014
).
3.
R.
Zbikowski
, “
Sensor-rich feedback control: A new paradigm for flight control inspired by insect agility
,”
IEEE Instrum. Meas. Mag.
7
,
19
26
(
2004
).
4.
S.
Sterbing-D'Angelo
,
M.
Chadha
,
C.
Chiu
,
B.
Falk
,
W.
Xian
,
J.
Barcelo
,
J. M.
Zook
, and
C. F.
Moss
, “
Bat wing sensors support flight control
,”
Proc. Natl. Acad. Sci. U. S. A.
108
,
11291
11296
(
2011
).
5.
A.
Saini
,
S.
Narsipur
, and
A.
Gopalarathnam
, “
Leading-edge flow sensing for detection of vortex shedding from airfoils in unsteady flows
,”
Phys. Fluids
33
,
087105
(
2021
).
6.
M. L.
Provost
,
X.
He
, and
D. R.
Williams
, “
Real-time roll and pitching moment identification with distributed surface pressure sensors on a UCAS wing
,” in
AIAA Aerospace Sciences Meeting
(AIAA,
2018
).
7.
L.
Burelle
,
W.
Yang
,
F.
Kaiser
, and
D. E.
Rival
, “
Exploring the signature of distributed pressure measurements on non-slender delta wings during axial and vertical gusts
,”
Phys. Fluids
32
,
115110
(
2020
).
8.
K. T.
Wood
,
S.
Araujo-Estrada
,
T.
Richardson
, and
S.
Windsor
, “
Distributed pressure sensing-based flight control for small fixed-wing unmanned aerial systems
,”
J. Aircr.
56
,
1951
1960
(
2019
).
9.
W.
Hou
,
D.
Darakananda
, and
J.
Eldredge
, “
Machine-learning-based detection of aerodynamic disturbances using surface pressure measurements
,”
AIAA J.
57
,
5079
5093
(
2019
).
10.
B.
Pamadi
,
J.
Brandon
,
V.
Klein
, and
P.
Murphy
, “
Prediction of unsteady aerodynamic coefficients at high angles of attack
,” in
AIAA Atmospheric Flight Mechanics Conference and Exhibit
(
AIAA
,
2001
), p.
4077
.
11.
P. C.
Murphy
,
V.
Klein
, and
N. T.
Frink
, “
Nonlinear unsteady aerodynamic modeling using wind-tunnel and computational data
,”
J. Aircr.
54
,
659
683
(
2017
).
12.
X.
An
,
D. R.
Williams
,
J. D.
Eldredge
, and
T.
Colonius
, “
Lift coefficient estimation for a rapidly pitching airfoil
,”
Exp. Fluids
62
,
11
(
2021
).
13.
J. L.
Callaham
,
J. V.
Koch
,
B. W.
Brunton
,
J. N.
Kutz
, and
S. L.
Brunton
, “
Learning dominant physical processes with data-driven balance models
,”
Nat. Commun.
12
,
1016
(
2021
).
14.
Y.
Zou
,
R.
Donner
,
N.
Marwan
,
J.
Donges
, and
J.
Kurths
, “
Complex network approaches to nonlinear time series analysis
,”
Phys. Rep.
787
,
1
97
(
2019
).
15.
G.
Iacobello
,
L.
Ridolfi
, and
S.
Scarsoglio
, “
A review on turbulent and vortical flow analyses via complex networks
,”
Physica A
563
,
125476
(
2021
).
16.
K.
Taira
,
A.
Nair
, and
S.
Brunton
, “
Network structure of two-dimensional decaying isotropic turbulence
,”
J. Fluid Mech.
795
,
R2
(
2016
).
17.
M.
Gopalakrishnan Meena
and
K.
Taira
, “
Identifying vortical network connectors for turbulent flow modification
,”
J. Fluid Mech.
915
,
A10
(
2021
).
18.
A.
Krishnan
,
R.
Sujith
,
N.
Marwan
, and
J.
Kurths
, “
Suppression of thermoacoustic instability by targeting the hubs of the turbulent networks in a bluff body stabilized combustor
,”
J. Fluid Mech.
916
,
A20
(
2021
).
19.
T.
Kobayashi
,
S.
Murayama
,
T.
Hachijo
, and
H.
Gotoda
, “
Early detection of thermoacoustic combustion instability using a methodology combining complex networks and machine learning
,”
Phys. Rev. Appl.
11
,
064034
(
2019
).
20.
M.
Murugesan
and
R.
Sujith
, “
Combustion noise is scale-free: Transition from scale-free to order at the onset of thermoacoustic instability
,”
J. Fluid Mech.
772
,
225
245
(
2015
).
21.
G.
Iacobello
,
S.
Scarsoglio
,
J.
Kuerten
, and
L.
Ridolfi
, “
Lagrangian network analysis of turbulent mixing
,”
J. Fluid Mech.
865
,
546
562
(
2019
).
22.
D.
Perrone
,
J.
Kuerten
,
L.
Ridolfi
, and
S.
Scarsoglio
, “
Wall-induced anisotropy effects on turbulent mixing in channel flow: A network-based analysis
,”
Phys. Rev. E
102
,
043109
(
2020
).
23.
A.
Shirazi
,
G.
Jafari
,
J.
Davoudi
,
J.
Peinke
,
M.
Tabar
, and
M.
Sahimi
, “
Mapping stochastic processes onto complex networks
,”
J. Stat. Mech.
2009
,
P07046
.
24.
A.
Campanharo
,
M. I.
Sirer
,
R. D.
Malmgren
,
F. M.
Ramos
, and
L.
Amaral
, “
Duality between time series and networks
,”
PLoS One
6
,
e23378
(
2011
).
25.
M.
McCullough
,
K.
Sakellariou
,
T.
Stemler
, and
M.
Small
, “
Regenerating time series from ordinal networks
,”
Chaos
27
,
035814
(
2017
).
26.
D.
Fernex
,
B. R.
Noack
, and
R.
Semaan
, “
Cluster-based network modeling—From snapshots to complex dynamical systems
,”
Sci. Adv.
7
,
eabf5006
(
2021
).
27.
E.
Kaiser
,
B.
Noack
,
L.
Cordier
,
A.
Spohn
,
M.
Segond
,
M.
Abel
,
G.
Daviller
,
J.
Östh
,
S.
Krajnović
, and
R.
Niven
, “
Cluster-based reduced-order modelling of a mixing layer
,”
J. Fluid Mech.
754
,
365
414
(
2014
).
28.
H.
Li
,
D.
Fernex
,
R.
Semaan
,
J.
Tan
,
M.
Morzyński
, and
B. R.
Noack
, “
Cluster-based network model
,”
J. Fluid Mech.
906
,
A21
(
2021
).
29.
F.
Foroozan
,
V.
Guerrero
,
A.
Ianiro
, and
S.
Discetti
, “
Unsupervised modelling of a transitional boundary layer
,”
J. Fluid Mech.
929
,
A3
(
2021
).
30.
A.
Nair
,
C.
Yeh
,
E.
Kaiser
,
B.
Noack
,
S.
Brunton
, and
K.
Taira
, “
Cluster-based feedback control of turbulent post-stall separated flows
,”
J. Fluid Mech.
875
,
345
375
(
2019
).
31.
S.
Lloyd
, “
Least squares quantization in PCM
,”
IEEE Trans. Inf. Theory
28
,
129
137
(
1982
).
32.
D.
Arthur
and
S.
Vassilvitskii
, “
k-means++: The advantages of careful seeding
,” Report No. 2006-13 (
Stanford
,
2006
).
33.
M.
Newman
,
Networks
(
Oxford University Press
,
2018
).
34.
A.
Savitzky
and
M. J. E.
Golay
, “
Smoothing and differentiation of data by simplified least squares procedures
,”
Anal. Chem.
36
,
1627
1639
(
1964
).
35.
F.
Kaiser
,
G.
Iacobello
, and
D. E.
Rival
, “
Aerodynamic state estimation from sparse sensor data by pairing Bayesian statistics with transition networks
,” in
AIAA Aerospace Sciences Meeting
(AIAA2022).