University of Southern California Dissertations and Theses

Reconstruction and estimation of shear flows using physics-based and data-driven models

(USC Thesis Other)

Reconstruction and estimation of shear flows using physics-based and data-driven models

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content RECONSTRUCTION AND ESTIMATION OF SHEAR FLOWS USING PHYSICS-BASED AND DATA-DRIVEN MODELS by Vamsikrishna Chinta A Dissertation Presented to the FACULTY OF THE USC GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirement for the Degree DOCTOR OF PHILOSOPHY (MECHANICAL ENGINEERING) December 2022 Copyright 2022 Vamsikrishna Chinta Dedicated to my parents for their unconditional support and encouragement that enabled me to complete the journey to the highest rung of education. ii Acknowledgments First and foremost, I would like to thank my adviser, Dr. Mitul Luhar for his constant encouragement, patience, and support. I really enjoyed all the technical discussions with him. I learnt a lot from him especially his advising style: how encouragement can be a far more eective motivator than fear and intimidation. His advising style made me curious to explore other elds and broaden my research interests. My PhD thesis involved working with various collaborators. I am grateful to have made these professional connections, which I believe is also an important part of PhD journey. His promptness in replying to emails and providing feedback on manuscripts was a strong encouragement for me. Lastly, his ever cheerful attitude made working with him a pleasant experience. I would also like to take this opportunity to thank my adviser at Indian Institute of Science, Dr. Jaywant Arakeri. He is denitely one of the people with out whose help I would not have been able to reach the PhD program. Working with him kindled my interest in uid mechanics, boosted my self condence to pursue a PhD, and helped me overcome many diculties along the way. My long endeavor for research and graduate with a PhD materialized because of these two amazing advisers. My parents have always been very supportive of my graduate studies. Their unwavering support during the pandemic and during my master's helped me endure during some very tough times. A special thanks to my friend and colleague, Chan-Ye (Chris) Ohh. I had the oppor- tunity of working with her on the collaborative strat wakes project with Dr. Spedding. I thoroughly enjoyed working with them. Chris patiently answered several questions I had (O(10 2 )) via emails during the pandemic lockdown and helped me out with a video call to clear my doubts before the conference presentation of SoCal uids 2021. Her cheerful atti- tude also made her a very pleasant oce roommate. During my time as a PhD student, I had iii the privilege of working with two collaborators of Dr. Luhar: Dr. Georey Spedding, and Dr. Maziar Hemati. Both of them have been very kind, supportive, and encouraging. They also provided valuable feedback regarding the project, conference presentations, as well as manuscripts. I would also like to thank Dr. Ketan Savla for many fruitful discussions as well as serving on my thesis committee. Andrew deserves a special mention here. He sel essly oered rides to the conference venue or airport. He was a go to computer expert, a very decent person, and a passionate researcher. Another person who was equally passionate and nice discuss about career or uid mechanics was Anup Kanale. I thoroughly enjoyed learning about low-Reynold's number uid mechanics. I still remember the time when he presented the paper \Life at low reynolds number" by EM Purcell. That discussion made me feel that uid mechanics or physics in general is so vast and so interesting that it is unfortunate that we get to do research in a very narrow topic of this vast ocean. A special thanks to Michael Kruger for sharing the template of this thesis that is consistent with the formatting style re- quired by the USC graduate school. For reference to other graduating students, the template is available at https://github.com/michK. I would also like to thank all the other uids supergroup members: Shilpa, Christoph, Emma, Morgan, Madeleine, JP, Vanessa, Chase, Mark, Idan, and Sultan. My stay here at USC was made enjoyable by my friends and roommates: Melih, Han ching, Aneesh, Rodolfo, Gaurav, Mukul, Ankit, and Omkar. I dearly miss Aneesh whose journey at USC was cut short due to unfortunate circumstances. Although this is not an exhaustive list people who helped me during my PhD, I would like to thank all the higher order friendships. The support and resources provided by USC are gratefully acknowledged. Lastly, I am grateful for the nancial support provided by the graduate school through the Provost fellowship. This thesis contains work supported by AFOSR grant FA9550-19-1- 7027 and ONR grant N00014-22-1-2029 and is gratefully acknowledged. The data generated by Chan-Ye Ohh and Dr. Spedding from ONR Grant No. N00014-20-1-2584 under Dr. Peter Chang is utilized in this thesis and is gratefully acknowledged. iv Contents Dedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Acknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xviii 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Data Assimilation for Fluid Flows . . . . . . . . . . . . . . . . . . . . . . . . 2 1.3 Thesis Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 2 Temporal reconstruction of turbulent channel ows from non-time-resolved PIV measurements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1.1 Problem Description . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1.2 Previous Reconstruction Eorts . . . . . . . . . . . . . . . . . . . . . 7 2.1.3 Rapid Distortion Theory and Taylor's Hypothesis . . . . . . . . . . . 9 2.1.4 Contribution and Outline . . . . . . . . . . . . . . . . . . . . . . . . 10 2.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.2.1 Rapid Distortion Theory and Taylor's Hypothesis . . . . . . . . . . . 12 2.2.2 Fusion of Forward and Backward Estimates . . . . . . . . . . . . . . 14 2.2.3 Numerical Evaluation and Error Quantication . . . . . . . . . . . . 17 2.3 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 2.3.1 Reconstruction Accuracy for Benchmark Test Case . . . . . . . . . . 22 2.3.2 Eect of Measurement Spatiotemporal Resolution . . . . . . . . . . 29 2.3.3 Reconstructed Statistics and Spectra . . . . . . . . . . . . . . . . . . 34 2.3.4 Eect of Noise on Reconstruction Accuracy . . . . . . . . . . . . . . 39 2.4 Chapter Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3 Spatio-temporal resolvent-based estimation of turbulent channel ows . . . . . . . 44 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.2.1 Resolvent Mode Computation and Database . . . . . . . . . . . . . . 47 v 3.2.2 Problem Formulation . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.2.3 Identifying Dominant Resolvent Modes Using FROLS . . . . . . . . . 55 3.2.4 Enforcing Statistical Constraints . . . . . . . . . . . . . . . . . . . . 59 3.2.5 Numerical Evaluation . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 3.2.6 Cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 3.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 3.3.1 Temporal Reconstructions . . . . . . . . . . . . . . . . . . . . . . . . 63 3.3.2 Estimation of Unobserved Variables . . . . . . . . . . . . . . . . . . . 66 3.4 Chapter conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 4 Reconstruction and estimation for stratied wakes . . . . . . . . . . . . . . . . . . 74 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 4.1.1 Motivation and Problem Statement . . . . . . . . . . . . . . . . . . . 74 4.1.2 Previous Data-Driven Classication and Reconstruction Eorts . . . 79 4.1.3 Contribution and Outline . . . . . . . . . . . . . . . . . . . . . . . . 81 4.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 4.2.1 Numerical Simulations and DMD Database . . . . . . . . . . . . . . . 83 4.2.2 Laboratory Experiments . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.2.3 Problem Formulation . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 4.2.4 Identifying Dominant DMD Modes Using FROLS . . . . . . . . . . . 93 4.2.5 Estimating Regime and Quantifying Accuracy and Condence . . . . 95 4.3 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 4.3.1 2D2C PIV Snapshots from DNS . . . . . . . . . . . . . . . . . . . . . 99 4.3.2 2D2C PIV Snapshots from Experiments . . . . . . . . . . . . . . . . 107 4.3.3 Point Measurements . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 4.4 Chapter Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 5.1 Implications and Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . 120 Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 vi List of Tables 2.1 DNS datasets acquired from the JHTDB [39]. The time interval between individual snapshots in the dataset is t + and the the total number of snap- shots is N. For datasets 1 and 3, multiple snapshot sequences are obtained for ensemble averaging purposes. Dataset 1 includes 48 dierent realizations obtained at 6 dierent spatial locations and over 8 dierent time windows. Dataset 3 includes 12 realizations obtained at dierent spatial locations. The last column lists the section(s) of the chapter in which results corresponding to each dataset appear. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.2 Description of the dierent reconstruction techniques used in this study. In general, superscripts (+,, or) represent time evolution while subscripts (t or tx) represent the fusion scheme. The () notation denotes characteristics- based evolution of the ow eld with TH. . . . . . . . . . . . . . . . . . . . . 21 3.1 Cases discussed in this study. In the second column, represents resolvent analysis with the molecular viscosity while T is for an eddy viscosity model. N is the number of snapshots used for training. The fourth column indicates the planes and velocity data available for training. is the hyperparameter used in constrained optimization, see Eq. 3.27 . . . . . . . . . . . . . . . . . 62 vii 3.2 Reconstruction errors in the velocity eld and statistical quantities for the trained and watched variables. For each case, the top row shows reconstruc- tion error after the FROLS algorithm is used to identify and calibrate 200 resolvent mode amplitudes based on the training snapshots. The bottom row (in parentheses) shows reconstruction error after the optimization algorithm modies the amplitudes to enforce statistical constraints. . . . . . . . . . . . 68 4.1 The accuracy (, Eq. 4.17) and condence (, Eq. 4.18) metrics, averaged over all 24 (Re;Fr) cases, using the full ow eld and the mean-subtracted ow eld. Results are shown for dierent input data and database sizes. The numbers in parentheses (10 or 30) indicate the number of 2D snapshots used in the input dataset. For the cases in which 10 snapshots are used, the accuracy and condence metrics are averaged over 3 dierent realizations. . . . . . . . 99 viii List of Figures 2.1 Left: Non-time resolved PIV snapshots stacked along the time axis. Right: Time-resolved reconstruction of the velocity eld from the PIV snapshots. . . 5 2.2 Schematic showing the combined Region of In uence (ROI) and Domain of Dependence (DOD) of two consecutive snapshots in thext plane at a given wall-normal location. Characteristics with slope determined by the local mean velocity (dt=dx = 1=U(y)) are shown by the solid black lines at the yellow- green and green-blue interfaces. The green region shows where the ROI of Snapshot 1 (S1) and the DOD of Snapshot 2 (S2) coincide. In this region both the forward and backward estimates are used for fusion. . . . . . . . . 16 2.3 Time variation of integrated reconstruction error for the benchmark case with T + 16 and x + = y + 16. The error is computed using (2.12) and averaged over 48 dierent realizations. The solid, dashed, and dotted lines represent the mean and the color bands are spaced at 1 standard deviation from the mean line. Table 2.2 provides a description of the dierent recon- struction techniques used in this gure. . . . . . . . . . . . . . . . . . . . . 23 2.4 Wall-normal variation in reconstruction error as a function of time for the (a) RDT tx (b) TH tx and (c) TH tx models. The error is calculated using (2.13) and averaged over 48 dierent realizations. . . . . . . . . . . . . . . . . . . 26 ix 2.5 Snapshots of the velocity eld from DNS (a,b), the RDT tx reconstruction (c,d), and the TH tx reconstruction (e,f) in the middle of the time horizon, when error is maximum. Proles of streamwise and wall-normal velocity at y + 200 are plotted in panels (g) and (h). The bold solid lines show the reconstructed velocity eld for RDT tx , dotted lines are for TH tx , and the ne solid lines show the DNS velocity eld. The locations of the shocks in the weighting functions for the reconstruction are highlighted in red. . . . . . . . 27 2.6 Maximum reconstruction error ( max ) as a function of the grid resolution (x + ) and time horizon T + for (a) RDT tx , (b) TH tx , and (c) TH tx . The contour lines are shown at intervals of 0.05, and the benchmark case is shown as the black cross. Figure 2.7 shows the variation in integrated error as a function of time for the case marked with a black square (x + = 4;T + = 32). 30 2.7 Variation in integrated error (2.12) as a function of time for the RDT, TH, and TH reconstructions at the grid resolution and time horizon marked with a black square (x + = 4;T + = 32) in gure 2.6. . . . . . . . . . . . . . . . 31 2.8 Comparison of wall-normal proles of (a) u 2 , (b) v 2 , and (c) uv from DNS with theTH tx and TH tx reconstructions. The solid black lines show statistics computed directly from DNS data. The dashed lines show statistics for TH tx and the dotted lines show statistics for TH tx . Reconstructions were carried out with the DNS data sub-sampled at intervals of T + 16. . . . . . . . . . 34 2.9 Premultiplied power spectral density for streamwise velocity fE uu =u 2 , (a) - (c), and wall-normal velocities fE vv =u 2 , (d) - (f). Panels (a) and (d) are obtained from DNS data. Panels (b) and (e) show the power spectra from the PIV-like data, i.e., DNS data sub-sampled at a frame rate correspond- ing to T + 16. Panels (c) and (f) show the power spectra obtained from the characteristics-based TH tx reconstruction. Contour lines are shown at intervals of 0.2 for panels (a)-(c), and 0.04 for panels (d)-(f). . . . . . . . . 37 x 2.10 Error in reconstructed power spectral densities for (a) streamwise velocity, and (b) wall-normal velocity. The error is computed by taking the dierence between power spectra obtained from DNS and the TH tx reconstruction. The horizontal white line corresponds to the Nyquist limit for the PIV-like data sampled atT + 16. Contour lines are shown at intervals of 0.05 for (a), and 0.01 for (b). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 2.11 Variation of integrated error as a function of time with 4 dierent initial noise levels for (a) TH tx and (b) TH tx . The solid lines show results averaged over 48 dierent realizations. A shaded gray band representing one standard deviation above and below the average value from the 48 realizations is shown for the noisiest case, i.e., SNR = 5. . . . . . . . . . . . . . . . . . . . . . . 40 3.1 Schematic of the resolvent framework. The resolvent operator can be assem- bled by dierent k i.e., dierent wavenumber and frequency combinations. Resolvent modes then can be obtained via a gain-based decomposition of the resolvent operator. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 3.2 Subplots (a)-(c) show the dierent types of measurement planes we use for this study. (a) shows a single measurement plane in xy. (b) shows the measurement planes in two orthogonal xy and yz planes. (c) shows the measurement planes in two orthogonal xy and xz planes at y + = 200. Subplot (d) shows the target 3D velocity eld we aim to reconstruct. x, y, and z are the streamwise, wall-normal, and spanwise directions respectively. 52 3.3 Schematic showing successive steps in the FROLS algorithm. Subplot (a) shows the database of basis functions in blue. The training vector (input data) is shown in red and the closest match to the training vector is shown in green. Subplot (b) shows the remaining basis functions after the projection along the green vector in (a) is subtracted. . . . . . . . . . . . . . . . . . . . 57 xi 3.4 (a) Overview of the sparse regression technique. The FROLS algorithm iden- ties and calibrates the dominant resolvent modes from a large database that best represent the measurements. The weight vector is sparse. (b) Reconstruc- tions with 200 resolvent modes in comparison with the DNS ground `truth'. . 59 3.5 Flow elds of stream wise and wall-normal uctuations of velocity. subplots (a) and (b) are for DNS and subplots (c)-(h) are for ow elds reconstructed by increasing the number of modes. . . . . . . . . . . . . . . . . . . . . . . . 65 3.6 Error computed using equation 3.26 is plotted with increasing number of resolvent modes. The error decreases steeply in the beginning but later it asymptotes. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 3.7 Reconstructed statistics from FROLS (blue line), after constrianed optimiza- tion (red dots), and the statistics from DNS (red dashed line) are shown for three statistical quantities for Case 1. The statistics from FROLS are several orders of magnitude larger and are shown on the left axis (blue). The statistics after constrained optimization match the DNS statistics closely. . . . . . . . 66 3.8 The error of reconstruction from FROLS, after constrained optimization, and linear interpolation between the snapshots is shown in gure on the left. The gure on the right shown the reconstructed ow eld foru andv components along with the DNS truth. The error after optimization is slightly higher than the FROLS reconstructions and the ow elds match qualitativley and quantitatively with the the DNS. . . . . . . . . . . . . . . . . . . . . . . . . 67 3.9 Error evolution for trained variables in the xy plane (a) and yz plane (b). Solid lines represent reconstructions from FROLS and dashed lines from constrained optimization. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 xii 3.10 Reconstructed statistics from FROLS for Case 2 and constrained optimization for Cases 2 and 3. The blue scale on the left corresponds to the FROLS recon- structions with the orange scale on the right corresponds to the optimization results and DNS proles. Note the substantial dierence in magnitude from the FROLS proles to the optimization and DNS proles. The statistics from eddy viscosity model are smoother and more closely approximate the DNS proles. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 3.11 Reconstructed velocity and pressure elds after optimization for Case 2 are compared with DNS data. Subplots (a)-(b) show reconstructions for the ob- served (trained) components of velocity in the xy plane. Subplots (c)-(d) show reconstruction for the unobserved (watched) velocity and pressure elds. Subplots (e)-(f) show reconstructions at a dierentxy plane (z=h = 0:125) for the trained components (u;v) . . . . . . . . . . . . . . . . . . . . . . . . 71 4.1 Regime space in terms of Reynolds number (Re) and Froude number (Fr). Filled circles show the cases used in this study. Physical regimes identied in previous studies are denoted by the color of the circles [57, 73] and the back- ground shading [55, 24]. Regimes are color-coded as: vortex street (purple), symmetric non-oscillation (red), asymmetric non-oscillation (yellow), planar oscillation (green), and spiral mode (blue). Overlapping shaded regions re ect dierences in regime boundaries identied in previous studies. . . . . . . . . 76 xiii 4.2 Snapshots of normalized streamwise velocity (u) from DNS (left column) and cross-stream velocity components (v, w) for the dominant DMD mode (right column) in dierent regimes: (a) vortex street at (Re;Fr) = (1000; 0:5) , (b) symmetric non-oscillation at (Re;Fr) = (500; 1), (c) asymmetric non- oscillation at (Re;Fr) = (300; 4), (d) planar oscillation at (Re;Fr) = (500; 4), and (e) spiral mode at (Re;Fr) = (1000; 8). Shading for the DMD modes shows positive (red) and negative (blue) values for v and w, normalized by 0:75max(jvj;jwj). The Strouhal number, St = fD=U o , for the oscillation frequency f of the DMD modes is provided. . . . . . . . . . . . . . . . . . . 77 4.3 Meshes used for numerical simulations (a) at Re =f200; 300; 500g and (b) at Re = 1000. The mesh used for the lower Re cases has four renement levels (2D [25] D). The mesh used for the higher Re case has six renement levels (2D [27] ). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 4.4 (a) Schematic showing experimental setup. The annotated scale is based on a typical sphere diameter ofD = 4 cm. (b) Schematic showing temporal data reconguration where the colored box A is moving in reference to the sphere. 86 4.5 (a) Streamwise velocity contours normalized by freestream velocityU o . Three point sensors are shown in blue, green, and red. These sensors are displaced from a common central location, (x;y;z) = (24R;0:25R;0:25R), by a dis- tance of 1:25R in the streamwise, lateral, and vertical directions, respectively (n.b., the green sensor is not in the same xz plane as the blue and red sensors). (b) Streamwise, vertical, and spanwise velocity signals are plotted against time for all the three sensors. The signals are color coded by sensor. 88 4.6 (a) Overview of the sparse regression technique. The FROLS algorithm iden- ties and calibrates the dominant DMD modes from a large database that best represent the input measurements. The weight vector is sparse. (b) Reconstructed snapshots in comparison with the the input ow elds. . . . . 94 xiv 4.7 Component-wise alignment between mean ow elds is evaluated between each of the (Re;Fr) cases in gure 4.1. There are a total of 24 cases, yielding a 24 24 grid. (a) Alignment between streamwise mean velocity elds,ju i u j j=(ku i kku j k), and (b) Alignment between vertical mean velocity elds,jw i w j j=(kw i kkw j k). Here, the () notation implies an element-wise product and sum. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 4.8 Comparisons between streamwise (a, b) and vertical velocities (c, d) over the testing data for (Re;Fr) = (500; 4). (e) Evolution of error with number of iterations for training and validation data sets. The minimum point for the validation error (i.e., the optimal iteration) is highlighted with a red circle. (f) The evolution of the predictions Re w and Fr w is shown until the optimal iteration. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 4.9 Plots showing prediction accuracy and condence averaged over three dierent sets of input data for all (Re;Fr) cases in the database. Regime predictions made using full ow elds as input data are shown in (a,c) and predictions generated using mean-subtracted ow elds are shown in (b,d). Large boxes represent prediction accuracy (Eq. 4.17) while smaller inset boxes represent condence (Eq. 4.18). Panels (a,b) show the accuracy computed using both Re and Fr, while panels (c,d) show the accuracy computed using Fr alone. 102 4.10 Comparison of reconstructions between streamwise (a, b) and vertical veloc- ities (c, d) over the testing data for (Re;Fr) = (500; 4) using a sub-sampled DMD database. (e) Evolution of error with number of iterations for train- ing and validation data sets. The minimum point for the validation error (i.e., the optimal iteration) is highlighted with a red circle. (f) The evolu- tion in Re w and Fr w until the optimum iteration. The nal prediction is (Re w ;Fr w ) = (664; 3:9) with a condence of 94%. . . . . . . . . . . . . . . . 105 xv 4.11 Plots showing prediction accuracy and condence averaged over three dierent sets of input data. These predictions are made using a database only consisting of DMD modes from the cases marked with a cross (). Regime predictions made using full ow elds as input data are shown in (a) and predictions generated using mean-subtracted ow elds are shown in (b). Large boxes represent prediction accuracy (Eq. 4.17) while smaller inset boxes represent condence (Eq. 4.18). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 4.12 Comparisons between measurements and reconstructions for the streamwise (a, b) and vertical velocity (c, d) over the testing data for the experimental dataset at (Re;Fr) = (533; 3:7). Note the diering spatial extent of the wakes for the experimental and DNS data. (e) Evolution in error for the training and validation data sets. The minimum point for the validation error (i.e., the optimal iteration) is highlighted with a red circle. (f) The evolution of Re w and Fr w until the optimal iteration. The nal prediction is (Re w ;Fr w ) = (349; 7:5) with a condence of 83% per Eq. 4.18. . . . . . . . . 108 4.13 Prediction accuracy (large boxes, Eq. 4.17) and condence (smaller inset boxes, Eq. 4.18) averaged over three snapshot sequences from experiments, using the full DMD mode database. As before, panel (a) shows results for the full ow eld while panel (b) shows results obtained using the mean-subtracted ow eld. Experimental data for theFr = 16 cases are unavailable. However, the corresponding DMD modes from DNS have been included in the database. 109 4.14 Prediction accuracy (large boxes, Eq. 4.17) and condence(smaller inset boxes, Eq. 4.18) using experimental measurements of only the vertical componentw, averaged over 3 dierent sets of input data. As before, (a) shows results obtained using the full ow elds while (b) shows results obtained using the mean-subtracted ow elds. . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 xvi 4.15 Comparisons between the predictor (^ z) and the DNS (z) over the training (a) and testing (b) datasets. These reconstructions correspond to sensor location 1. The evolution of reconstruction error with the number of iterations for the training and validation sets is shown in (c). The corresponding variation of Re w and Fr w is shown in (d). The input data corresponds to (Re;Fr) = (500; 4) and the regime prediction corresponds to (Re w ;Fr w ) = (500; 4) with a condence of 97%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 4.16 Prediction accuracy (large boxes, Eq. 4.17) and condence (smaller inset boxes, Eq. 4.18) with point measurements as input data using the full DMD mode database. Subplots (a), (b), and (c) correspond to sensor locations 1, 2, and 3 respectively. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113 4.17 Error in predicted sensor location normalized byR. Subplots (a), (b), and (c) correspond to sensor locations 1, 2, and 3 respectively. Note that the error is capped at 1R and red shading indicates lower error. . . . . . . . . . . . . . . 113 xvii Abstract The ability to generate spatially and temporally resolved velocity measurements in turbulent ows is essential for many applications, ranging from numerical weather predictions to ac- tive ow control. However, even state-of-the-art measurement techniques are often limited in their ability to fully resolve the spatio-temporal uctuations inherent in turbulent ows due to hardware constraints. Therefore, the ability to infer dynamic parameters and generate higher-resolution reconstructions of turbulent ows from limited measurements can benet many applications. This thesis develops dynamical models that can make accurate predic- tions of wall-bounded turbulent ows over short time horizons from limited measurements. These models are derived from the governing equations that ensures physical consistency of the reconstructed ow eld. The advantages and disadvantages of each model are discussed, and can potentially be used for data assimilation. The dynamical models developed here can be used to estimate ow eld in time, space, and possibly unmeasured variables. We also introduce a novel method to reconstruct the ow eld for a stratied wake to estimate its dynamic regime. Stratied wakes of blu bodies are known to exhibit a rich variety of three dimensional structures that enables us to classify these ow structures into dierent dynamical regimes. Previous studies had access to 3D velocity elds or 2D velocity elds in two orthogonal planes to classify these structures in to ve dierent classes. Here, we develop a purely data-driven method to reconstruct the input measurements, which in turn enables direct estimation of the dynamic parameters, Reynolds number (Re) and Froude number (Fr). In the rst part of this thesis, we develop models grounded in Rapid Distortion The- ory (RDT) and Taylor's Hypothesis (TH) to reconstruct the time evolution of a turbulent ow eld in the intermediate period between consecutive Particle Image Velocimetry (PIV) snapshots obtained using a non-time resolved system. The linear governing equations are xviii evolved forwards and backwards in time using the PIV snapshots as initial conditions. The ow eld in the intervening period is then reconstructed by taking a weighted sum of the forward and backward estimates. This spatiotemporal weighting function is designed to ac- count for the advective nature of the RDT and TH equations. Reconstruction accuracy is evaluated as a function of spatial resolution and reconstruction time horizon using Direct Numerical Simulation (DNS) data for turbulent channel ow from the Johns Hopkins Tur- bulence Database (JHTDB). This method reconstructs single-point turbulence statistics well and resolves velocity spectra at frequencies higher than the temporal Nyquist limit of the acquisition system. Reconstructions obtained using a characteristics-based evolution of the ow eld under TH prove to be more accurate compared to reconstructions obtained from numerical integration of the discretized forms of RDT and TH. The eect of measurement noise on reconstruction error is also evaluated. The second part of this thesis develops a resolvent-based framework that can be used to reconstruct the velocity eld for turbulent channel ow from limited measurements. For re- construction, resolvent modes computed at dierent frequency and wavenumbers are pooled into a large database. The Forward Regression with Orthogonal Least Squares (FROLS) algorithm is used to sequentially identify modes that best represent the input data. A com- posite loss function is dened with the error for the second-order statistics imposed as a soft constraint to the reconstruction error for velocity eld. The amplitude and phase of the identied modes are then optimized to minimize this composite loss function. After cali- bration, the velocity eld can be reconstructed at arbitrary spatiotemporal resolution using the weighted resolvent modes. The resulting model enables estimation of out-of-plane com- ponents of velocity and pressure as well as reconstruction in other spanwise locations. For proof-of-concept tests of this method, DNS data for turbulent channel ow from JHTDB are used. Results indicate that using measurements from multiple planes results in reduced error for the untrained variables (i.e., non-measured velocity components, pressure, out-of-plane data). Using an eddy viscosity model to compute resolvent modes results in smoother sta- xix tistical proles, but the velocity eld reconstructions do not reproduce small-scale features. In the last part of this thesis, we attempt to identify the dynamic regime from lim- ited measurement data in a stratied wake with (nominally) unknown Re and Fr. A large database of candidate basis functions is compiled by pooling the Dynamic Mode Decom- position (DMD) modes obtained in prior DNS. A sparse model is built using the FROLS algorithm, which sequentially identies DMD modes that best represent the data and cali- brates their amplitude and phase. After calibration, the velocity eld can be reconstructed using a weighted combination of the dominant DMD modes. The dynamic regime for the measurements is estimated via a projection-weighted average ofRe andFr corresponding to the identied modes. Regime identication is carried out from a limited number of 2D ve- locity snapshots from numerical and experimental datasets, as well as 3 point measurements in the wake of the body. A metric to assess condence is introduced based on the observed predictive capability. This approach holds promise for the implementation of data-driven uid pattern classiers. xx Chapter 1 Introduction 1.1 Motivation The ability to reconstruct uid ows and estimate key dynamic parameters is important in many applications of scientic and engineering interest. Such eorts are especially chal- lenging for turbulent ows, which are characterized by features with a wide range of spatial and temporal scales. Despite signicant advances in measurement hardware as well as post- processing and analysis algorithms, current measurement techniques are limited in their ability to fully resolve turbulent ows. This motivates the need to develop methods capable of reconstructing turbulent ows from limited and possibly noisy measurements. In many real-world situations, acquiring the full state of the system is prohibitively expensive and impractical. Often it is only possible to sample a subset of the state variables at low resolu- tion and rely on other techniques to reconstruct these ows. In these situations, estimates from one or more dynamical models can be fused with the measurements (e.g., using the well known Kalman Filter) to reconstruct uid ows at higher resolution. This data assimilation process is used in many applications and has its roots in meteorology. Successful data assim- ilation and ow reconstruction for turbulent ows has the potential to benet applications ranging from weather prediction and climate modeling to active ow control 1 1.2 Data Assimilation for Fluid Flows Data assimilation is the science of combining dierent sources of information (e.g., measure- ments and model predictions) to estimate the state of a system as it evolves in time. In the context of turbulent ows, data assimilation techniques can be used to synergistically combine limited measurements with model predictions to generate ow eld estimates with higher delity. The reliance on models is a central feature of turbulent ow estimation due to impracticality of Direct Numerical Simulation (DNS) in all conditions[65]. Models can be physics based [46, 47] or data driven [12, 45, 87]; a subset of models involve Galerkin projec- tion onto basis functions obtained from prior data or via analysis of the governing equations [7, 102, 83]. On the other hand, experimental data are essential for model development and testing. However, generating spatio-temporally resolved measurements for turbulent ows is challenging despite many recent advances in measurement techniques [2, 106]. For this reason, previous studies have attempted turbulent ow reconstruction using a variety of tech- niques such as projection onto basis functions. Tu et al. [98] and Discetti, Raiola, and Ianiro [31] used basis functions derived from data-driven techniques, while G omez et al. [37] and Beneddine et al. [11] used resolvent modes (physics-based) as the basis functions. Proper Or- thogonal Decomposition (POD) and Dynamic Mode Decomposition (DMD) are data-driven techniques, where the amplitudes and modal basis functions are computed from the data. In contrast, resolvent modes are computed directly from the governing equations. Thus, when using resolvent modes, the key ow features are identied directly from the Navier-Stokes equations. Their modal amplitudes, however, should be calibrated using the measurements from experiments. For completeness, we note that there are several machine learning tech- niques used for super-resolution reconstruction [35, 33, 34], and ow eld reconstruction from wall-based measurements [40, 41]. 2 1.3 Thesis Outline In this thesis, we focus on reconstruction and estimation of ow variables for wall-bounded turbulent ows in chapters 2 and 3. We use two dierent approaches: in chapter 2, we use a dynamical model | a simplied partial dierential equation derived from the governing Navier-Stokes equations | to reconstruct the ow eld. The use of this dynamical model is to be able to predict ow elds with reasonable accuracy in the time between two consecutive measurements. In contrast, in chapter 3, we project the acquired velocity eld on to resolvent modes which serve as the basis functions. These resolvent modes are three dimensional travelling wave ow structures. This enables reconstruction in space and time and also estimation of unmeasured variables. In both these chapters, we use models that are derived from the governing equations so that the reconstructed ow eld will be physically consistent. This is in contrast to previous studies that use purely data-driven approaches [98, 13, 31], where physical consistency is not guaranteed. In chapter 4, we reconstruct the ow eld in order to estimate the dynamic regime for stratied wakes. This is in contrast to chapters 2 and 3, where we do not try to estimate the ow variables and reconstruct higher spatio- temporal ow elds, but only estimate the regime of the ow. In chapter 2 we reconstruct the time evolution of turbulent channel ow between two consecutive non-time-resolved 2D2C PIV snapshots. The PDE of the dynamical model are evolved forwards and backwards in time using the PIV snapshots as the initial conditions. The forward and backward estimates are fused using weights. The accuracy of this method is evaluated by a direct comparison between the DNS data for turbulent channel ow at friction Reynolds numberRe = 1000 obtained from the Johns Hopkins Turbulence Database (JHTDB) [39]. This chapter has been published as Krishna et al. [51]. In chapter 3 we use physics-based model|projection of the velocity eld onto resolvent modes|to reconstruct the velocity eld from limited non-time-resolved measurements for turbulent channel ow. Resolvent modes are generated via a gain-based decomposition of 3 the governing equations, ensuring physical consistency. For reconstruction, a large database of resolvent modes is created. The Forward Regression with Orthogonal Least Squares (FROLS) algorithm is then used to identify the dominant resolvent modes. A composite loss function is dened with the loss from the second order statistics imposed as a soft constraint to the loss from velocity eld. The amplitude and phase of the identied modes are then optimized to satisfy the input velocity eld and the statistics in a least squares sense. After calibration, the velocity eld can be reconstructed at arbitrary spatiotemporal resolution using the weighted resolvent modes. The weighted resolvent modes also enable estimation of out-of-plane components of velocity and pressure as well as in other spanwise locations. For proof-of-concept tests of this method, we again use DNS data of turbulent channel ow from JHTDB [39]. Part of this chapter has been published in Chinta and Luhar [21] Finally, in chapter 4 we attempt regime identication and ow eld reconstruction for stratied wakes with (nominally) unknown Reynolds number (Re) and Froude number (Fr). We start with the assumption that a large library of wake patterns from known regimes is available. We then use the sparse regression algorithm (FROLS) to identify a limited number of wake patterns from the library that best match the measurements. The patterns identied by this algorithm are then used to infer the regime for the measurements. Regime identication is carried out from a limited number of 2D2C velocity led snapshots from numerical and experimental datasets, as well as 3 point measurements in the wake of the body. This chapter has been published as Chinta et al. [22] 4 Chapter 2 Temporal reconstruction of turbulent channel ows from non-time-resolved PIV measurements 2.1 Introduction 2.1.1 Problem Description The last three decades have seen rapid advances in the development of high-power and high-repetition rate lasers, high-speed digital cameras capable of megapixel resolution, and computing power. These hardware advances, together with improvements in the speed and Streamwise Wall-normal 2h 2h Figure 2.1: Left: Non-time resolved PIV snapshots stacked along the time axis. Right: Time-resolved reconstruction of the velocity eld from the PIV snapshots. 5 accuracy of Particle Image Velocimetry (PIV) analysis algorithms, have led to a step change in our ability to make non-intrusive eld measurements at high spatial and temporal reso- lution [106, 2]. Despite these advances, laboratory PIV systems are often limited in their ability to fully resolve the broadband spatiotemporal uctuations inherent in turbulent ows due to hardware limitations or cost constraints. For example, even state-of-the-art PIV sys- tems with kHz-capable cameras and lasers may not yield complete temporal resolution in certain conditions. These limitations in hardware motivate the need to reconstruct the ow eld from limited and noisy measurements. In this chapter, we develop models that can reconstruct the time evolution of wall- bounded turbulent ows in the intervening period between two PIV snapshots from a non- time resolved system. Figure 2.1 provides a schematic view of the problem being addressed. The left panel in Figure 2.1 shows PIV snapshots stacked along the time axis. Consistent with typical planar PIV systems, we assume that only two-dimensional, two-component (2D-2C) snapshots are available. The sampling time between PIV measurements is T . As shown in the right panel, the goal is to reconstruct the evolution of the ow eld between two consecutive snapshots with high temporal resolution, i.e., to generate predictions for the snapshots shown as translucent planes. We focus on turbulent channel ow and use direct numerical simulation (DNS) data from the Johns Hopkins Turbulence Database [39] to develop and test the models used for ow reconstruction. However, we expect these techniques to be equally applicable to other wall-bounded ows (e.g., pipe and boundary layer ows) and be directly transferable to physical PIV systems. Most prior eorts on turbulent ow reconstruction have relied on data-driven ap- proaches, in which sensor data are used in conjunction with static correlation maps or projected onto basis functions obtained from eld measurements [e.g., 3, 15, 94, 98, 13, 31]. In contrast, the present eort employs simplied models based on Rapid Distortion Theory (RDT) and Taylor's Hypothesis (TH) that are derived from the governing Navier-Stokes equations. A brief review of previous reconstruction eorts and models grounded in 6 RDT/TH is provided below. 2.1.2 Previous Reconstruction Eorts From a signal processing point of view, time resolution issues in measurements can be allevi- ated through the use of compressed sensing, particularly when the signal of interest is sparse in frequency space [e.g., 19, 9, 100]. If such narrow-banded signals are sampled randomly in time, then ` 1 minimization techniques can be used to reconstruct their time evolution even with sub-Nyquist average sampling rates. However, such techniques are less reliable for signals that are not sparse in frequency space, as is the case with broad-banded turbulent ows. Moreover, it is not typically possible to generate PIV measurements that are sampled randomly in time. Previous studies have compensated for the limited time resolution of PIV systems through the use of two complementary instruments. In this multi-sensor fusion approach, the high-spatial and low-temporal resolution velocity elds from PIV are typically fused with time-resolved point measurements from instruments such as hot-wire anemometers (HWA) or pressure sensors [98, 13, 31]. Broadly, such sensor fusion and ow eld reconstruction tech- niques can either involve a static approach, such that the estimator relies purely on statistics compiled from prior data, or a dynamic approach, in which an underlying evolution model is included in the estimation procedure. Stochastic estimation with static maps has been employed extensively in the turbulence community [e.g., 3, 26, 42, 69], particularly in conjunction with proper orthogonal decom- position (POD) [15, 94]. Variants of Linear Stochastic Estimation (LSE) have been used in many dierent contexts. For instance, POD-based LSE has been used to educe coherent structure in turbulent ows [4, 15, 42]; multi-sensor stochastic estimation has been used to estimate velocity elds from wall-based pressure and shear stress measurements [68, 70, 32]; and LSE incorporating time delays has been used to predict the future evolution of ow elds for purposes of control [94, 101]. Along similar lines, Sasaki et al. [85] have recently used 7 Large Eddy Simulation data to identify linear and nonlinear transfer functions that enable the estimation of streamwise velocity uctuations in turbulent boundary layers based on measurements at dierent wall-normal locations. Of relevance to the present study, LSE has also been employed to fuse information from fast-time point measurements and slow-time eld data for turbulent channel ow. Specically, Discetti, Raiola, and Ianiro [31] evaluated the correlation between fast-time point measurements and the time variation in POD modes obtained from the eld measurements, at the instances in which these datasets were syn- chronized. This static correlation map was then used to infer the fast-time evolution of the POD modes from the point measurements. Dynamic estimators seek to incorporate information from underlying dynamic models as well as sensor measurements to estimate the state of high-dimensional systems such as turbulent ows. Standard techniques for linear dynamic estimation include Kalman ltering and Kalman smoothing. Kalman ltering estimates the current state of the system based on previous observations, while Kalman smoothing also allows for renement of previous estimates in light of later observations. Tu et al. [98] employed a Kalman smoother to successfully fuse information from fast-time point measurements and slow-time eld mea- surements to reconstruct the velocity eld in the wake of a thick at plate with an elliptical leading edge at low Reynolds number. The dynamic model employed in this study involved projection onto POD modes computed from the eld measurements. Tu et al. [98] found that a dynamic model comprising the following two elements was sucient for accurate ow eld reconstruction: (i) the rst two POD modes oscillating stably at the shedding frequency measured from the probe signal, with amplitudes estimated from the point measurements, and (ii) the remaining POD mode coecients advanced in time via LSE. For the purposes of more general turbulent ow reconstruction, dynamic models like the one developed by Tu et al. [98] have two important limitations. First, the assumption of stable, low dimensional oscillatory dynamics only holds for a small class of narrow-banded ows (e.g., low Reynolds number blu body wakes). Second, projection onto POD modes limits any model predic- 8 tions to the subspace spanned by the data, and does not guarantee that the resulting ow elds will be physically sound. This limitation is especially problematic in the context of the broadband turbulent ows to be considered in this study. Indeed, the necessary data for extracting the appropriate POD-basis (or a reliable data-driven model by other means) are often unavailable, even when using state-of-the-art instrumentation for ow diagnostics. To circumvent these issues relating to data availability, recent studies have attempted ow reconstruction by projecting the velocity eld onto so-called resolvent modes that are obtained via a gain-based decomposition of the governing equations [61, 60]. More details about this technique, how other studies such as G omez et al. [37] and Beneddine et al. [11] have used this technique to reconstruct ows, their short comings, and propose a new method to attempt reconstructing turbulent ows and estimate unobserved are given in Chapter 3 2.1.3 Rapid Distortion Theory and Taylor's Hypothesis As noted earlier, in this chapter we use models grounded in Rapid Distortion Theory (RDT) and Taylor's hypothesis (TH) to reconstruct the ow eld between two consecutive PIV snapshots. In essence, RDT assumes that if a turbulent ow eld is subjected to substantial distortion by the mean shear ow, the higher-order nonlinear interactions can be neglected when predicting the early response. Scaling arguments show that RDT is formally correct if the time-horizon of prediction is much shorter than a typical eddy turnover time [46]. This makes RDT a natural choice for the present eort, in which the velocity eld must be evolved from a known initial state over a short time horizon until the next PIV snapshot becomes available. RDT has been used extensively in both theoretical and experimental turbulence research [86, 46]. Moreover, there are strong connections between RDT and the resolvent analysis framework mentioned in the previous section [60] since both approaches emphasize linear dynamics. Under Taylor's Hypothesis (TH), the RDT equations are simplied further to retain just the time derivative of the uctuations and the mean ow advection term. In other words, TH 9 assumes that the turbulent ow eld is `frozen' and advects downstream with the mean ow [93]. The mean velocity can therefore be used to convert temporal information to spatial in- formation, and vice versa. TH has been used extensively in the turbulence community [109, 29, 28, 64, 112], primarily to infer spatial structure from time-resolved point measurements (e.g., from hot-wire anemometers). Experimental estimates of streamwise wavenumber (k x ) spectra are often converted from frequency (!) spectra under the assumption that the re- sulting convection velocity (c = !=k x ) is equal to the local mean velocity [88]. This is a good assumption in the outer region of the ow, where the convection velocity of the tur- bulent uctuations is known to be close to the local mean velocity. However, experiments and numerical simulations both show that this assumption leads to an underestimate of the convection velocity below the buer region of the ow. Specically, the convection velocity remains at c + 10 below y + 15, even as the mean velocity goes to zero at the wall [see e.g., 50, 52, 78]. Here,y is the wall-normal coordinate (y = 0 at the wall) and a superscript + denotes normalization with respect to the friction velocity and kinematic viscosity. Thus, the use of the local mean velocity in TH leads to underestimation of streamwise length scales in the near-wall region. Moreover, wall-bounded turbulent ows are typically characterized by a broad spectrum of frequencies for a given wavenumber (and vice versa). Hence, using the same convection velocity for all wavenumber-frequency combinations at a given wall-normal location can also lead to spurious peaks in the power spectral density [28, 64]. Since TH has been used to infer spatial structure from time-resolved measurements with some success in turbulent ows, it should also be possible to invoke this hypothesis to solve the inverse problem: to infer the time evolution of a ow eld from measurements of spatial structure. The reconstruction framework described in this chapter develops this idea further. 2.1.4 Contribution and Outline To reconstruct the time evolution of the ow eld between consecutive 2D-2C PIV snap- shots, we use the simplied equations obtained under RDT and TH to generate predictions 10 forwards in time from the initial snapshot and backwards in time from the subsequent snap- shot. A weighted sum of these forward and backward estimates is used to reconstruct the ow eld in the intervening period. The accuracy of the reconstructed ow elds obtained using these models is assessed using DNS data for turbulent channel ow at friction Reynolds number Re = 1000 obtained from the Johns Hopkins Turbulence Database (JHTDB)[39]. Results show that this reconstruction framework signicantly outperforms direct interpola- tion, i.e., the reconstructed velocity elds deviate much less from the DNS data compared to interpolated velocity elds. Moreover, frequency spectra computed from the reconstructed velocity elds are found to closely resemble spectra obtained from the DNS data, even at frequencies higher than the Nyquist limit of the PIV-like data. Note that the reconstruction framework developed here only makes use of 2D-2C velocity snapshots and simplied models grounded in the Navier-Stokes equations. Unlike previous eorts [98, 31], no additional time- resolved point measurements are used to improve temporal resolution. Instead, the observed improvement in temporal resolution stems from the use of RDT and TH to reconstruct the time evolution of the ow eld from the spatial information present in the snapshots. The remainder of this chapter is structured as follows. In Section 2.2, we provide a brief review of the governing equations obtained under RDT and TH. We also describe the methods used for generating the forward- and backward-time predictions, the weighting functions used to fuse these forward- and backward-time reconstructions of the ow eld, and the metrics used for error quantication. In Section 2.3, we assess reconstruction accuracy for the various models developed in Section 2.2 (e.g., RDT vs. TH, dierent weighting functions). We also evaluate the eects of measurement spatiotemporal resolution and noise on reconstruction accuracy, and compare turbulence statistics and spectra obtained from the reconstructed ow elds against DNS results. Finally, we present concluding remarks in Section 2.4. 11 2.2 Methods 2.2.1 Rapid Distortion Theory and Taylor's Hypothesis Under Rapid Distortion Theory, the Navier-Stokes equations are linearized about the mean prole [10, 86, 46] to yield the following momentum equation and continuity constraint: @u @t + Uru + urU =rp + 1 Re r 2 u + (NL); (2.1) and r u = 0: (2.2) In the expressions above, U = (U(y); 0; 0) represents the mean prole, u = (u;v;w) denotes the turbulent velocity uctuations, p represents pressure uctuations, and (NL) represents the (neglected) nonlinear terms. A standard Cartesian coordinate system is used, in which x is the streamwise direction,y is the wall-normal direction, andz is the spanwise direction; t is time. Scaling arguments show that the nonlinear terms can be neglected in turbulent shear ows for time horizons that are shorter than the typical eddy turnover time [86, 46]. This makes RDT an appropriate choice for the present problem requiring temporal reconstruc- tion between sequential PIV snapshots. However, even with the substantial simplication aorded by linearization, reconstruction based on the full RDT equations is dicult in prac- tice. This is because most common PIV systems are only capable of generating planar 2D-2C eld measurements. Assuming these PIV measurements are carried out in the (x;y) plane to yield velocity components (u;v), additional simplifying assumptions are needed to account for the out-of-plane ow and pressure gradient terms. Here, we simply neglect these terms to yield the following coupled advection-diusion equations for streamwise and wall-normal 12 velocity: @u @t +U @u @x | {z } advection term = 1 Re @ 2 u @x 2 + @ 2 u @y 2 | {z } diusion term v @U @y |{z} coupling term ; (2.3) and @v @t +U @v @x | {z } advection term = 1 Re @ 2 v @x 2 + @ 2 v @y 2 | {z } diusion term : (2.4) In other words, the continuity constraint is not enforced. This ad hoc simplication is not rigorously justied. However, solving two-dimensional versions of (2.1) and (2.2) also introduces additional modeling assumptions that are not rigorously justied and requires solution of the pressure Poisson equation. Moreover, the initial and nal velocity elds do not satisfy the two-dimensional continuity equation. So, imposing this constraint for the intermediate reconstructions could lead to additional numerical errors. Since the goal here is to generate simple models that can be used for ow reconstruction, we proceed with (2.3)- (2.4). These simplied, linearized versions of the Navier-Stokes equations are referred to as RDT for the remainder of this chapter to acknowledge their conceptual origin (though we recognize that this terminology is not entirely accurate). Under Taylor's frozen turbulence hypothesis, the equations above are further simplied by assuming that the advection term is dominant. This yields @u @t +U @u @x = 0: (2.5) As noted earlier, TH essentially assumes that the turbulent velocity eld is `frozen' in form and advects downstream with the mean ow. So, the local mean velocity can be used to convert between temporal variations and spatial variations. Below, we employ the linear models grounded in RDT (2.3)-(2.4) and TH (2.5) to reconstruct the time evolution of a turbulent channel ow in the time interval (T ) between two planar 2D-2C eld measurement snapshots (e.g., from PIV). Note that the only input required for these reconstructions is the mean velocity prole appearing in (2.3)-(2.5), which 13 can be obtained from the snapshots themselves. Keep in mind that reconstruction can proceed both forwards and backwards in time. In other words, the equations above can be evolved forwards in time using the rst snapshot as the initial condition, as well as backwards in time using the second snapshot as the initial condition. These forward- and backward-time predictions can be generated via numerical integration of appropriately discretized versions of (2.3)-(2.5). In addition, the advection equation in (2.5) can also be evolved backwards or forwards in time using the method of characteristics. When using discretized versions of the governing equations, the backward- time integration uses the transformation = Tt, such that the new time variable has value = 0 for the nal snapshot at t = T , and value = T for the initial snapshot at t = 0. Note that this transformation switches the sign of the time derivative term in equations (2.3)-(2.5). This is equivalent to solving the RDT equations with negative convection velocity, mean shear, and viscosity, and solving the TH equations with negative convection velocity. An appropriately weighted combination of these forward- and backward- time estimates has the potential to improve reconstruction accuracy. Below, we develop physically motivated weighting schemes for this fusion of the forward- and backward-time estimates. Of course, the time evolution of the ow between the two snapshots can also be estimated via linear interpolation, which does not require any underlying models or weighting schemes. In Section 2.3, we show that the reconstruction framework developed in this chapter yields velocity estimates that are signicantly more accurate than estimates obtained via direct interpolation. 2.2.2 Fusion of Forward and Backward Estimates The simplied linear equations obtained under RDT (2.3)-(2.4) and TH (2.5) are evolved forwards in time from the rst snapshot to yield a forward estimate for the velocity eld, ^ u f , and backwards in time from the second snapshot to yield a backward estimate, ^ u b . The reconstructed ow eld is given by a weighted average of these forward and backward 14 estimates: ^ u =k f ^ u f +k b ^ u b : (2.6) A simple way to fuse the forward and backward estimates is to use weights that vary linearly in time k f =k f (t) = 1 t T ; k b =k b (t) = t T : (2.7) Here t = 0 corresponds to the initial snapshot and t =T corresponds to the nal snapshot. This particular weighting scheme ensures that the forward-time estimate is weighted more heavily closer to the initial snapshot and the backward-time estimate is weighted more heavily towards the nal snapshot. This weighting scheme can be improved further by considering the mathematical nature of the equations emerging from RDT and Taylor's hypothesis. Since the hyperbolic advection term is expected to be dominant in wall-bounded turbulent ows, information is expected to propagate at a speed corresponding to the local mean velocity, U(y). This is illustrated in thext diagram shown in gure 2.2. The region of in uence (ROI) for the rst snapshot and the domain of dependence (DOD) for the second snapshot are determined by characteristics in thext plane that emanate from the upstream (x = 0) and downstream (x =L x ) edge of the snapshots and have slope dt=dx = 1=U(y). For advection-dominated ows, information propagation from the forward-time evolution is conned to the ROI of the rst snapshot, and information propagation from the backward-time evolution is conned to the DOD of the second snapshot. In other words, the forward-time estimate is expected to be accurate only in the ROI of the rst snapshot (i.e., green and blue regions in gure 2.2) while the backward-time estimate is expected to be accurate only in the DOD of the second snapshot (i.e., yellow and green regions in gure 2.2). Further, since the slope of the characteristics that dene the ROI and DOD is determined by the mean velocity, the size of the ROI and DOD also varies with y. To account for these eects, the weighting scheme in (2.7) can be modied as follows. 15 Snapshot 1 Snapshot 2 Common Region of Predictability Outside the DOD of S2 Outside the ROI of S1 l b l f x t L x T Figure 2.2: Schematic showing the combined Region of In uence (ROI) and Domain of Dependence (DOD) of two consecutive snapshots in the xt plane at a given wall-normal location. Characteristics with slope determined by the local mean velocity (dt=dx = 1=U(y)) are shown by the solid black lines at the yellow-green and green-blue interfaces. The green region shows where the ROI of Snapshot 1 (S1) and the DOD of Snapshot 2 (S2) coincide. In this region both the forward and backward estimates are used for fusion. The linear weighting scheme in (2.7) can be retained in the common region of predictability for both snapshots (green region in gure 2.2). The forward weight, k f , is set to 0 in the region outside the ROI of the rst snapshot (yellow region in gure 2.2) and the backward weight,k b , is set to 1. Similarly,k f = 1 in the region outside the DOD of the second snapshot (blue region in gure 2.2) and k b = 0. The resulting equations for the weights are: k f =k f (x;y;t) = 8 > > > > > > < > > > > > > : 0 0x<l f 1 t T l f xL x l b 1 L x l b <xL x ; (2.8) 16 and k b =k b (x;y;t) = 8 > > > > > > < > > > > > > : 1 0x<l f t T l f xL x l b 0 L x l b <xL x : (2.9) In the expressions above, l f = U(y)t, l b = U(y)(Tt), and L x is the streamwise extent of the PIV window. The upstream edge of the PIV window corresponds to x = 0 and the downstream edge corresponds to x = L x . Note that the weights shown in (2.8)-(2.9) are dependent on (x;y) as well as t, in contrast to temporal weighting scheme with k f = k f (t) andk b =k b (t) shown in (2.7). Moreover, notice that the weights are discontinuous at x =l f (yellow-green interface in gure 2.2) and at x =L x l b (green-blue interface in gure 2.2). Hence, this weighting scheme introduces spatial shocks in the reconstructed ow eld at those locations. 2.2.3 Numerical Evaluation and Error Quantication To test reconstruction accuracy for the forward, backward, and fused ow eld estimates, we use DNS data for turbulent channel ow at Re = u h= = h + = 1000 available from the JHTDB [39]. Here, h is the channel half-height, u is the friction velocity, and is kinematic viscosity. For consistency with basic PIV systems, we use 2D-2C velocity data in thexy plane that are sampled uniformly in time and space (with one exception where the data are sampled with logarithmic spacing in the y-direction, as discussed below). In other words, we use systematically sub-sampled DNS data as a surrogate for PIV measurements. The results presented below make use of three complementary DNS datasets (see Table 2.1). For each dataset, the PIV window extends across the entire height of the channel (2h). The streamwise extent is also set to L x = 2h. Dataset 1 in Table 2.1 is used to evaluate reconstruction accuracy for both RDT and TH for a benchmark test case. This case has a uniform spatial grid comprising 129 129 points across the PIV window including the walls at the top and bottom. A total of N = 17 Dataset Grid Resolution t + N Realizations Section 1 x + = y + 16 0.0649 256 48 2.3.1, 2.3.4 2 x + = y + 4 0.0649 512 1 2.3.2 3 x + 16 0.649 1024 12 2.3.3 Logarithmic in y Table 2.1: DNS datasets acquired from the JHTDB [39]. The time interval between indi- vidual snapshots in the dataset is t + and the the total number of snapshots is N. For datasets 1 and 3, multiple snapshot sequences are obtained for ensemble averaging purposes. Dataset 1 includes 48 dierent realizations obtained at 6 dierent spatial locations and over 8 dierent time windows. Dataset 3 includes 12 realizations obtained at dierent spatial locations. The last column lists the section(s) of the chapter in which results corresponding to each dataset appear. 256 snapshots are acquired at a sampling rate of t + = 0:0649. Only the rst and last snapshots in this dataset are used for reconstruction and so the prediction time horizon is T + =Nt + 16. The intervening snapshots from DNS are used to evaluate reconstruction accuracy. To evaluate statistical variations in reconstruction accuracy for this benchmark case, similar data are extracted at 6 dierent spatial locations and for 8 dierent time windows. In other words, this dataset comprises 48 independent realizations with identical spatiotemporal resolution and prediction time horizon. Reconstruction accuracy for this benchmark case is discussed in Section 2.3.1. Dataset 2 in Table 2.1 is acquired at higher spatial resolution with 513 513 uniformly sampled points across the 2h 2h PIV window. This dataset includes N = 512 snapshots obtained at intervals of t + = 0:0649, for a total time horizon of T + = Nt + 33. In Section 2.3.2, this dataset is sub-sampled systematically to evaluate the eect of spatiotem- poral resolution and prediction time horizon on reconstruction accuracy for both RDT and TH. Dataset 3 is used to compute turbulence statistics and spectra in Section 2.3.3. For this, we use logarithmic spacing with 129 grid points across the channel to better evaluate reconstruction accuracy in the near-wall region. In the streamwise direction, we use a uni- form grid spacing of x + 16, similar to the benchmark case. For improved statistical 18 convergence, a total of N = 1024 snapshots are acquired at intervals of t + = 0:649 over 12 dierent spatial locations of the DNS domain. Reconstruction is carried out using every 25th snapshot in this dataset, so that the prediction time horizon T + = 25t + 16 is comparable to the benchmark case. Reconstruction is only carried out using the fused TH model for this dataset. As a point of comparison, for a physical system with PIV analysis being carried out for 16 pixel by 16 pixel segments with 50% overlap (i.e., 8 pixels between data points), the benchmark case with 129 129 uniformly distributed points represents the use of a camera with approximately 1 Megapixel (1 MP) resolution. The high-resolution dataset with 513 513 uniformly sampled points across the PIV window represents the use of a camera with 17 MP resolution. Similarly, the time horizon for the benchmark case,T + 16, corresponds to water ow with friction velocityu = p T + =T 0:04 ms 1 for a PIV system capable of 100 Hz sampling rate (T = 0:01 s) andu 0:15 ms 1 for a system capable of 1000 Hz sampling rate (T = 0:001 s). For air ow, the corresponding friction velocity estimates are u 0:15 ms 1 for a 100 Hz system and u 0:5 ms 1 for a 1000 Hz system. These estimates assume a kinematic viscosity of 10 6 ms 2 for water and 1:5 10 5 ms 2 for air. For the reconstruction, the simplied equations (i.e., (2.3)-(2.4) for RDT and (2.5) for TH) are numerically integrated forwards in time from the rst snapshot and backwards in time from the last snapshots over the prediction horizon. A standard nite dierence scheme is used for this purpose. An explicit Euler method is used for time integration, a rst-order upwinding scheme is used for the advection terms, and a second-order central dierencing scheme is used for the diusion and coupling terms. Numerical evaluation is carried out at the spatial resolution of the snapshots, and a time step of t + = 0:0649. The Courant{Friedrichs{Lewy (CFL) condition is satised for all parameter combinations. Due to the linear nature of the governing equations and the relatively short prediction horizons (512 time steps at most), the results presented below are not particularly sensitive to the 19 choice of the numerical method. However, the nite dierence discretization does introduce additional articial viscosity. The magnitude of this viscosity increases with increasing grid spacing [82]. For TH, the eects of the articial viscosity introduced by the numerical discretization of (2.5) can be eliminated by using the method of characteristics to evolve the ow eld in time. Specically, the solution in the intervening period between the snapshots can be obtained by simply propagating the initial and nal ow elds along characteristics determined by the mean velocity. The forward-time evolution of the ow eld can be computed from the initial snapshot u (x; y; t = 0) as: ^ u f (x;y;t) = u (xU(y)t; y; 0): (2.10) Similarly, the backward-time evolution of the ow eld can be computed from the nal snapshot u (x; y; t =T ) as: ^ u b (x;y;t) = u (x +U(y)(Tt); y; T ): (2.11) Here, ^ u f and ^ u b are the forward- and backward-time estimates andU(y) is the mean velocity. By construction, the forward-time estimate is conned to the ROI of the rst snapshot (l f xL x in gure 2.2), and the backward-time estimate is conned to the DOD of the second snapshot (0 x L x l b in gure 2.2). The full reconstructed ow eld can be obtained by fusing these forward- and backward-time estimates using the spatio-temporal weights shown in (2.8)-(2.9). Table 2.2 summarizes the dierent forward-time, backward-time, and fused reconstruc- tion techniques used in this study. Reconstruction accuracy is quantied using the following 20 Technique Description Interpolation Linear interpolation in time RDT + Forward time integration of (2.3)-(2.4) RDT Backward time integration of (2.3)-(2.4) RDT t Forward and backward time integration of (2.3)-(2.4) Fused using temporal weights (2.7) RDT tx Forward and backward time integration of (2.3)-(2.4) Fused using spatio-temporal weights (2.8)-(2.9) TH + Forward time integration of (2.5) TH Backward time integration of (2.5) TH t Forward and backward time integration of (2.5) Fused using temporal weights (2.7) TH tx Forward and backward time integration of (2.5) Fused using spatio-temporal weights (2.8)-(2.9) TH tx Characteristics-based evolution of (2.5) Fused using spatio-temporal weights (2.8)-(2.9) Table 2.2: Description of the dierent reconstruction techniques used in this study. In general, superscripts (+,, or) represent time evolution while subscripts (t ortx) represent the fusion scheme. The () notation denotes characteristics-based evolution of the ow eld with TH. integrated error metrics. The time-varying global error is dened as (t) = R 2h x=0 R 2h y=0 (u ^ u) 2 + (v ^ v) 2 dxdy 1=2 R 2h x=0 R 2h y=0 (u 2 +v 2 )dxdy 1=2 ; (2.12) while the wall-normal variation in error over time is dened as (y;t) = R 2h x=0 (u ^ u) 2 + (v ^ v) 2 dx 1=2 R 2h x=0 (u 2 +v 2 )dx 1=2 : (2.13) 21 In the expressions above, ^ u and ^ v are the reconstructed velocity uctuations, u and v are the velocity uctuations from DNS `truth', subsampled to match the PIV spatial resolution. The lower wall of the channel is located at y = 0 and the upper wall is located at y = 2h. 2.3 Results and Discussion In this section, we evaluate reconstruction accuracy for all the techniques listed in Table 2.2. Reconstruction accuracy for the benchmark test case for the dierent models is evaluated in Section 2.3.1. Specically, we evaluate the time evolution of global error in Section 2.3.1, the wall-normal variation in error over time in Section 2.3.1, and compare reconstructed ow elds with DNS snapshots in Section 2.3.1. The eect of eld measurement spatial resolution (x + = y + ) and time horizon (T + ) on reconstruction accuracy is considered in Section 2.3.2. Reconstructed statistics and frequency spectra are shown in Section 2.3.3. This proof-of-concept study is based on eld data from DNS. In Section 2.3.4, we add Gaussian random noise to the rst and last DNS snapshots used as inputs in the RDT/TH models to evaluate the eect of noisy real-world measurements on reconstruction accuracy. 2.3.1 Reconstruction Accuracy for Benchmark Test Case Variation in Global Error over Time First, we compare reconstruction accuracy across dierent models for the benchmark case with spatial resolution x + = y + 16 and prediction horizon T + 16 using dataset 1 (see Table 2.1). The evolution in global error (2.12) for the simplest possible reconstruction technique | linear interpolation between the snapshots | is shown as the black line in gure 2.3. Note that the error is averaged over 48 dierent spatiotemporal realizations from the DNS database. The mean error across these 48 dierent cases is plotted as a solid line and the shading represents one standard deviation about the mean. As expected the error is zero at the beginning and the end of the prediction horizon, where snapshots of the ow 22 0 5 10 15 0.2 0.4 0.6 0.8 1 Figure 2.3: Time variation of integrated reconstruction error for the benchmark case with T + 16 and x + = y + 16. The error is computed using (2.12) and averaged over 48 dierent realizations. The solid, dashed, and dotted lines represent the mean and the color bands are spaced at 1 standard deviation from the mean line. Table 2.2 provides a description of the dierent reconstruction techniques used in this gure. eld are available. The error reaches a maximum of max = max() 0:8 at the middle of the time horizon. Next, we evaluate reconstruction accuracy for the forward time integration of the RDT equations with the snapshot at t + = 0 used as the initial condition (denoted RDT + ). With this technique, the global error is 0 initially and grows monotonically with time (solid blue line in gure 2.3). The reconstruction error exceeds the error from linear interpolation at t + 10. Similarly, the backward time RDT reconstruction (RDT ) yields 0 error at the end of the time horizon (i.e., att + =T + the initial condition for the backward time integration) and increases monotonically as time decreases. The reconstruction error fromRDT exceeds that from linear interpolation for times beforet + 6. In other words, for this case, the error dynamics are similar for both the forward and backward time RDT models. Note that linear interpolation outperforms the RDT-based models because it uses both the rst and last snapshots over the prediction horizon for reconstruction. In contrast, the forward RDT model uses only the rst snapshot for reconstruction while the backward RDT 23 model uses on the last snapshot. To improve reconstruction accuracy for the RDT-based models, we can fuse the forward and backward time estimates using the weighting functions shown in (2.7) and (2.8)-(2.9). The fused RDT model that uses the temporal weighting function (2.7) is denotedRDT t . The reconstruction error for this fused model (cyan line in gure 2.3) has a similar trend to the linear interpolation method, i.e., it has 0 error at the beginning and the end and reaches a maximum in the middle of the time domain. However, the maximum error is signicantly lower: max 0:5 for RDT t compared to max 0:8 for linear interpolation. Figure 2.3 shows that the error from the fused estimate (RDT t ) exceeds that from forward and backward reconstructions alone at early and late times. Ideally, any fusion would yield reconstructions that are as good as, or better than, the individual RDT + and RDT estimates over the entire time horizon. To a large extent, this can be achieved by taking into account the advection dominated nature of the ow under consideration. As illustrated schematically in gure 2.2, the fused reconstruction can be further rened using the spatiotemporal weighting scheme given in equations (2.8)-(2.9). This spatiotemporal reconstruction, termedRDT tx , has the lowest maximum error of all the techniques considered thus far, with max 0:3 (green line in gure 2.3). Moreover, theRDT tx reconstruction yields errors comparable to RDT + and RDT at early and late times, respectively. Next, we evaluate reconstruction accuracy for models grounded in Taylor's hypothesis (2.5), and compare these reconstructions against those obtained using RDT. Similar to the RDT models, equation (2.5) can also be discretized and integrated forwards or backwards in time. These forward and backward time reconstructions under TH are denoted TH + and TH , respectively. The forward and backward TH predictions can be combined using the temporal weights shown in (2.7) to yield the fused estimate TH t , or the spatiotemporal weights shown in (2.8)-(2.9) to yield the fused estimate TH tx . The fused reconstruction obtained using the method of characteristics (2.10)-(2.11) instead of numerical integration is denoted TH tx . Note that we only consider the fused TH tx reconstruction because the 24 characteristics-based evolution of (2.5) is not well dened outside the ROI of the rst snap- shot and the DOD of the second snapshot. Figure 2.3 shows the time evolution of error for the forward, backward, and fused re- constructions obtained using TH as dashed lines. When the discretized form of (2.5) is used to generate the forward- and backward-time estimates, reconstruction accuracy for the TH models is nearly identical to that for the corresponding RDT models. However, reconstruc- tion performance improves further when the method of characteristics is used to evolve the ow eld forwards and backwards in time as shown by the dotted magenta line. Maximum reconstruction error for TH tx is max < 0:2, compared to max 0:3 for TH tx and RDT tx . This observation suggests that the articial viscosity introduced by the numerical discretiza- tion of (2.5) leads to a deterioration in reconstruction performance for the relatively coarse spatial resolution used in the benchmark test case (x + = y + 16). We explore this issue further in Sections 2.3.2 and 2.3.3. Similar performance across the discretized RDT and TH models suggests that advection is the dominant physical mechanism retained in (2.3)-(2.4). The additional terms accounting for viscous eects and the interaction between wall-normal velocity uctuations and mean shear appear to be less important. However, keep in mind that the planar approximations to the RDT equations used here also neglect viscous eects due to out-of-plane gradients in the velocity, which are likely to be more important than the viscous eects arising due to streamwise gradients in velocity. In other words, the neglected (@ 2 u=@z 2 ) term is expected to be larger in magnitude than the (@ 2 u=@x 2 ) term. Unfortunately, accounting for out-of- plane gradients in streamwise and wall-normal velocity is not possible with access to planar PIV measurements alone. Moreover, as we show in Section 2.3.3, reconstructions based on discretized versions of the governing equations underestimate the intensity of the wall-normal uctuations in velocity. This may explain why inclusion of the additional v(@U=@y) term in (2.3) does not yield substantially dierent reconstructions relative to those based on TH. 25 Figure 2.4: Wall-normal variation in reconstruction error as a function of time for the (a) RDT tx (b) TH tx and (c) TH tx models. The error is calculated using (2.13) and averaged over 48 dierent realizations. Wall-normal Variation in Error over Time The time-evolution of error for the fused RDT tx , TH tx , and TH tx models is shown as a function of wall-normal distance (2.13) in gure 2.4. Consistent with the plots in gure 2.3, the error is 0 at eachy-location at the beginning and end of the time horizon and maximum in the middle for all models, and TH tx has the lowest reconstruction error at all wall-normal locations. The error is in general higher in the inner region of the ow (below y + 200 or y=h 0:2) where turbulence production and turbulent kinetic energy are higher, and maximum reconstruction error increases closer to the wall. However, keep in mind that the rst grid point is at y + 16 due to the linear distribution of grid points. This means that buer region of the ow is not resolved completely. Figure 2.4(c) also shows the presence of distinct temporal oscillations in reconstruction error for TH tx . The period for these oscillations decreases with increasing distance from the wall and closely matches the time scale x + =U + . This suggests that the oscillations arise from specic grid locations leaving the ROI of the initial snapshot or entering the DOD of the later snapshot as time advances (see gure 2.2). The increase in reconstruction error with decreasing y + could be attributed to the re- duction in turbulent timescales near the wall. Recall that linearization of the NSE is only accurate for predictions over short periods of time. For y + 10, the integral timescale 26 Figure 2.5: Snapshots of the velocity eld from DNS (a,b), the RDT tx reconstruction (c,d), and the TH tx reconstruction (e,f) in the middle of the time horizon, when error is maximum. Proles of streamwise and wall-normal velocity at y + 200 are plotted in panels (g) and (h). The bold solid lines show the reconstructed velocity eld for RDT tx , dotted lines are for TH tx , and the ne solid lines show the DNS velocity eld. The locations of the shocks in the weighting functions for the reconstruction are highlighted in red. 27 for streamwise velocity uctuations has been estimated to be T + u 20 while that for the wall-normal uctuations has been estimated to be T + v < 10 [78]. Assuming these timescales are representative of typical eddy turnover times, the applicability of RDT and TH is ques- tionable in the near-wall region for the benchmark case with prediction horizon T + 16. Integral timescales for the velocity uctuations are known to increase with distance from the wall and reachT + u 110 andT + v 50 aty + 200 [23, 56]. So, the assumptions underlying RDT and TH are better satised with increasing distance from the wall, which leads to an improvement in reconstruction accuracy. Comparison of Reconstructed Flow Fields with DNS To provide further insight into the reconstructed ow elds obtained from the RDT tx and TH tx models, gure 2.5 compares spatial snapshots of the uctuating velocity elds from DNS and the reconstructions in the middle of the reconstruction time horizon (t=T = 0:5), when the error is maximum. Panels (a), (c), and (e) show streamwise uctuations from DNS, the RDT tx reconstruction, and the TH tx reconstruction, respectively. Similarly panels (b), (d), and (f) show wall-normal uctuations from DNS, the RDT tx reconstruction, and the TH tx reconstruction, respectively. These spatial snapshots show that the reconstructed velocity elds obtained fromRDT tx qualitatively capture the large-scale structure. However, they do not reproduce small-scale features of the turbulent ow eld, particularly in the vicinity of the upper and lower walls. As shown in gures 2.5(g) and 2.5(h), the reconstructed ow elds for RDT tx appear to have gone through a spatial low-pass lter when compared to the DNS. In contrast, reconstructed ow elds obtained from the TH tx model match the DNS results much more closely. This is also evident in gures 2.5(g) and 2.5(h), which show that the TH tx proles retain nearly all the small-scale features present in the DNS results, and there is minimal attenuation of uctuation intensity. Reconstructed ow elds obtained from the discretized TH tx model (not shown in g- ure 2.5) closely resemble those obtained from theRDT tx model, i.e., they also appear low-pass 28 ltered relative to DNS. This observation conrms that the articial viscosity introduced by numerical discretization of the governing equations is responsible for the smoothing eect observed for the RDT tx andTH tx models. Recall that these reconstructions are carried out using a much coarser spatial grid (x + = y + 16; see Table 2.1) compared to DNS. Since the articial viscosity introduced by discretization is linearly proportional to grid resolu- tion, we anticipate signicant smoothing for the benchmark case. The characteristics-based evolution of (2.5) used for the TH tx reconstruction introduces no such articial viscosity. Finally, recall that the spatiotemporal weighting functions shown in (2.8)-(2.9) have spatial shocks. These spatial shocks are not visible in the reconstructed ow elds shown in gures 2.5(c)-(f). However, the streamwise proles of velocity at y + 200 shown in gures 2.5(g) and 2.5(h) do show the presence of minor discontinuities in the reconstructed velocity eld. The locations of the shocks in the spatiotemporal weighting functions (2.8)- (2.9) are highlighted in red for the reconstructed velocity proles. The streamwise gradient in velocity is discontinuous at both ends of the red regions, but smooth elsewhere. Together, gures 2.3-2.5 suggest that reconstruction accuracy is similar for the discretized RDT and TH models for the benchmark test case. However, reconstruction accuracy improves further when the method of characteristics is used to generate TH-based reconstructions. Next, we assess the eect of measurement spatio-temporal resolution on reconstruction performance. 2.3.2 Eect of Measurement Spatiotemporal Resolution In this section, we evaluate the accuracy of reconstruction as a function of spatial resolution and prediction time horizon using dataset 2 (see Table 2.1). For the fused RDT and TH models with spatiotemporal weights, the integrated error is 0 at the beginning and end of the prediction horizon, and reaches a maximum in the middle. This maximum error is plotted as a function of the grid resolution x + (= y + ) and time horizonT + in gures 2.6(a), (b), and (c) for the RDT tx , TH tx , and TH tx reconstructions, respectively. The grid resolution and time horizon for the benchmark case considered above are shown as a black cross. 29 10 1 10 2 5 15 25 10 1 10 2 5 15 25 10 1 10 2 5 15 25 0 0.2 0.4 0.6 0.8 1 Figure 2.6: Maximum reconstruction error ( max ) as a function of the grid resolution (x + ) and time horizon T + for (a) RDT tx , (b) TH tx , and (c) TH tx . The contour lines are shown at intervals of 0.05, and the benchmark case is shown as the black cross. Figure 2.7 shows the variation in integrated error as a function of time for the case marked with a black square (x + = 4;T + = 32). Consistent with the results shown in Section 2.3.1 for the benchmark case, reconstruc- tion accuracy is broadly similar for the discretized forms RDT tx and TH tx models. The characteristics-based reconstruction, TH tx , shows similar trends, though reconstruction er- rors are generally lower and less sensitive to spatial resolution. For a given time horizon T + , the maximum error increases gradually as a function of x + for the RDT tx and TH tx reconstructions. For the TH tx reconstruction, error is relatively insensitive to grid resolu- tion below x + 20 and only increases signicantly beyond this threshold value for x + . The gradual increase in reconstruction error as a function of x + for the RDT tx and TH tx models can be attributed to the articial viscosity introduced by numerical discretization. The magnitude of this articial viscosity is expected to increase as a function of x + [82]. The TH tx reconstruction does not introduce any articial viscosity. However, for higher x + , the initial and nal snapshots include less information from smaller scale turbulent ow features. The forward and backward reconstructions are therefore unable to resolve these smaller scales, and reconstruction accuracy deteriorates. For a given grid resolution x + , the maximum reconstruction error increases as a func- tion of time horizon for all the models. As an example, for grid resolutions corresponding to the benchmark case (x + = y + 16) the maximum error for RDT tx model increases 30 0 10 20 30 0 2.5 5 Figure 2.7: Variation in integrated error (2.12) as a function of time for the RDT, TH, and TH reconstructions at the grid resolution and time horizon marked with a black square (x + = 4;T + = 32) in gure 2.6. from max 0:3 for T + 16 to max 0:5 for T + 32. This observation is consistent with trends in gure 2.7, which show that the forward and backward reconstruction errors increase monotonically with increasing and decreasing time, respectively. Hence, as prediction time horizon increases, any fusion of the forward and backward estimates is also expected to yield larger errors. In general, as discussed earlier, any RDT- or TH-based reconstructions are expected to become less accurate as the prediction time horizon increases relative to typical turbulence timescales. Notably, gure 2.6(a) shows that the error for theRDT tx reconstruction grows dramati- cally with time horizon beyondT + = 20 at the lowest grid spacing considered here, x + 4. In contrast, gures 2.6(b) and 2.6(c) show no such increase in error for high grid resolutions and long prediction horizons for the TH tx and TH tx reconstructions. To provide further insight into this observation, gure 2.7 shows the time evolution of error for all the RDT and TH reconstructions for the case shown as a black square in gure 2.6, which corresponds to x + 4 and T + 33. Compared to the forward RDT predictions (solid blue line), error for the backward RDT predictions (solid red line) grows much more steeply with time. For 31 instance, at t + 16, the error associated with the forward RDT prediction is 1, while that for the backward RDT prediction is 5. This blowup can be attributed to the vis- cous terms in the RDT equations (2.3)-(2.4). As noted earlier, when integrating backwards in time using the transformation = Tt, this viscous diusion becomes negative. This steepens spatial gradients in velocity and eventually leads to instability (especially if addi- tional noise is introduced). In contrast, for the forward RDT integration, viscosity serves to smooth spatial gradients and damp external noise. The blowup in reconstruction error due to this negative diusion phenomenon is most prominent for lower grid spacing and longer prediction horizons. This can be attributed to two reasons. First, snapshots with ner grid resolutions are likely to include information from smaller-scale turbulent uctuations with larger spatial gradients. These large spatial gradients will be further amplied over time due to the negative eective viscosity in the backward RDT estimates. Second, the nite dierence discretization scheme used here contributes additional articial viscosity. The magnitude of this articial viscosity is proportional to the grid spacing x + [82]. Hence, for higher grid spacing the eect of negative diusion in the backward RDT estimates is partially oset by articial viscosity. However, for grid spacing x + 4, the eect of this articial viscosity is outweighed by the negative diusion over longer prediction horizons. Neglecting the viscous diusion terms in the RDT equations (2.3) and (2.4), along with the coupling term between horizontal and wall-normal velocity v(@U=@y) yields equation (2.5), corresponding to Taylor's hypothesis. This is an advection equation for the velocity uctuations with the mean velocity as the convection speed. This means that the veloc- ity uctuations are advected downstream with speed U(y) for the forward TH estimates and upstream with speed U(y) for the backward TH estimates. The dashed blue, red, and green lines in gure 2.7 show reconstruction accuracy as a function of time for the forward, backward, and fused TH models, respectively. When integrating forwards in time the perfor- mance of the TH model is comparable to that for the forward RDT model. However, when integrating backwards in time, the TH model far outperforms the RDT model. While the 32 backward RDT model led to 5 in the middle of the prediction horizon, the backward TH model yields 1. This is comparable to the reconstruction error for the forward TH model, suggesting that the forward and backward time dynamics are similar under Taylor's hypothesis (as expected from the governing equations). Thus, eliminating the viscous dif- fusion term alleviates the sharp increase in reconstruction error observed for the backward RDT estimates. This also means that the fused TH tx model leads to signicantly lower reconstruction error relative to the fusedRDT tx model for ne spatial resolutions and longer prediction horizons (see solid and dashed green lines in gure 2.7). Note that reconstruction accuracy for the TH tx model is similar to that for the TH tx at the grid resolution considered in gure 2.7. As mentioned above, the articial viscosity introduced by numerical discretization increases as a function of grid spacing [82]. Hence, at smaller grid resolutions both discretized and characteristics-based TH reconstructions yield similar results. Together, these observations suggest that the TH tx and TH tx models, which fuse the forward and backward TH estimates using the spatiotemporal weights shown in (2.8)-(2.9), yield more robust and accurate reconstructions relative to the other techniques tested in this chapter. Keep in mind that the TH and TH models use spatial information to infer time evolu- tion. As a result, reconstruction accuracy is limited by the spatial resolution of the snapshots. The accuracy of the temporal reconstruction is expected to improve only if the frequency corresponding to the spatial Nyquist limit, f + s = U + =(2x + ), is higher than the temporal Nyquist frequency of the acquisition system, f + t = 1=(2T + ). The expression for f + s above translates the spatial Nyquist limit (i.e., only structures longer than 2x + can be resolved) into a frequency using Taylor's hypothesis. The requirement that f + s be larger than f + t translates into the following condition for spatial resolution: x + k;2 >> k;m > 0 (3.8) and the orthonormality condition Z 1 1 ~ f k;l ~ f k;m dy = lm ; Z 1 1 ~ u k;l ~ u k;m dy = lm (3.9) A schematic of the resolvent analysis is shown in gure 3.1. As mentioned above the nonlinear term is treated as an external endogeneous forcing to the linear resolvent opera- tor. A fourier transform implies that this decomposition can be carried out at each k. A gain-based decomposition then identies the most amplied modes corresponding to that particular wavenumber and frequency combination k. These resolvent modes are then used as a basis onto which the velocity eld obtained from DNS / experiments is projected. The velocity eld can be expressed in terms of the weighted resolvent modes as u(x;y;z;t) = X k Re 0 B @ k ~ u k (y)e i( z }| { k x x +k z z!t) 1 C A (3.10) Note that k , ~ u k (y), and e i are complex. Further simplifying the above equation and 49 0 500 1000 0 50 100 -1 0 1 Figure 3.1: Schematic of the resolvent framework. The resolvent operator can be assembled by dierent k i.e., dierent wavenumber and frequency combinations. Resolvent modes then can be obtained via a gain-based decomposition of the resolvent operator. taking only the real part, we get: u(x;y;z;t) = X k 8 < : R k ~ u R k cos ~ u I k sin | {z } m + I k ~ u R k sin ~ u I k cos | {z } m+1 9 = ; (3.11) The superscripts R and I denote the real and imaginary parts respectively. The velocity eld is thus a linear superposition of the basis functions labeled `m' and `m+1'. Note that for each k, we get two basis functions. An eddy viscosity model is included in resolvent operator to incorporate the eect of Reynolds stresses [67]. Because turbulent channel ows are broad banded, a large database of resolvent modes was generated for turbulent channel ow at Re = 1000 over a span of 0:125j x j=h 16; 0:125 z =h 16; 0<c + p 22:63; (3.12) 50 A database of 53040 resolvent modes is compiled. To limit database size only rank 1 modes are used. Note that these resolvent modes are obtained directly from the governing equations with the only input as the mean velocity prole. Hence the grid resolution can be arbitrarily set. This is especially useful to reconstruct the unobserved quantities or ll in missing data. We now try to reconstruct the input velocity eld using a linear combination of these resolvent modes. Therefore the problem at hand is to identify dominant modes and calibrate their amplitudes and phase. 3.2.2 Problem Formulation Brie y the problem of resolvent-mode-based reconstruction can be stated as follows. Given a time series input data and a large database of resolvent modes computed at dierent k, we seek a sparse representation of the training data in terms of the resolvent modes in the database. In other words, a model is constructed for the input ow eld, approximating it as a linear combination of the dominant resolvent modes identied from the large database. The coecients of the identied modes are optimized to approximate both the input velocity eld and the statistics. Once the dominant modes are identied calibrated, we attempt estimating the velocity eld at other spanwise locations, where measurements are not made, using Eq. 3.10. Because the resolvent modes are directly computed from the governing equations, for each k, the mode shapes of each of the three velocity components and pressure is known apriori. This facilitates estimation of unobserved variables. For example, when the input velocity eld measurements are 2D2C planar streamwise and wall-normal velocity components, we attempt to estimate the spanwise component and pressure. The input data used can be 2D2C velocity eld snapshots from single xy plane or simultaneous measurements from multiple orthogonal planes with a limited eld of view. For the rest of this chapter we use x coordinate for streamwise, y coordinate for wall-normal, andz coordinate for spanwise directions. The schematic of the dierent types of input data used in this study is shown in Fig. 3.2 (a) - (c). As noted earlier we use DNS data for 51 turbulent channel ow at Re = 1000 available from the JHTDB [39] as the input or the surrogate for experimental PIV measurements. Although the DNS domain is much larger we use snapshots that are 2hh inxy plane, 2hh inyz plane, and 2h 2h inxz plane. Subplot (d) shows the full 3D velocity eld data we attempt to reconstruct. (a) (b) (c) (d) Figure 3.2: Subplots (a)-(c) show the dierent types of measurement planes we use for this study. (a) shows a single measurement plane inxy. (b) shows the measurement planes in two orthogonalxy andyz planes. (c) shows the measurement planes in two orthogonal xy and xz planes at y + = 200. Subplot (d) shows the target 3D velocity eld we aim to reconstruct. x, y, and z are the streamwise, wall-normal, and spanwise directions respectively. We assume that the total velocity eld can be represented as a linear combination of 52 the resolvent modes: ^ u(x;y;z;t) =< 0 @ X k k ~ u k (y) exp (i(k x x +k z z!t)) | {z } p k 1 A : (3.13) Here, ~ u k (y) is the resolvent mode computed for wavenumber-frequency combination k = (k x ;k z ;!). The complex coecients k determine the amplitude and phase of the mode. Note that each term in the summation comprises two basis functions and two coecients corresponding to the real and imaginary components of ~ u k and k , respectively. Thus,M=2 resolvent modes yield M basis functions and M coecients. The reconstruction problem essentially amounts to tuning the coecients k using mea- surement data. For this, the 3D ow elds for each resolvent mode, p k , are evaluated with the same spatial and temporal dimensions at the input data to form candidate basis func- tions d m in a large database D. For instance, for reconstruction from 2D-2C data in the xy plane, the streamwise and wall normal components of velocity (u;v) are extracted with the same spatial and temporal resolution as the PIV-like input data. The time series of K velocity measurements (or inputs) obtained at intervals of t, fz 1 ; z 2 ; ; z K g, is rst stacked into a single column vector of the form: z = 2 6 6 6 6 6 6 6 4 z 1 z 2 . . . z K 3 7 7 7 7 7 7 7 5 : (3.14) Note that each individual element of z represents a complete 2D2C spatial snapshot from all the measurement planes either single or multiple orthogonal planes. Next, we assume that the input data can be modeled as a linear combination of the 53 compatible resolvent mode sequences: z< X k c k p k ! : (3.15) Here,c k 2C are complex coecients that calibrate the amplitude and phase of each resolvent mode sequence to match the training data, p k are the resolvent basis functions given in equation 3.13, and<() represents the real component. Equation 3.15 can be rewritten as z X k (<(c k )<(p k )=(c k )=(p k )) (3.16) in which=() represents the imaginary component. A large database is now constructed by pooling all the resolvent modes computed for dierent k in the range of Eq. 3.12. To avoid computation with complex coecients, the real and imaginary parts of each compatible resolvent mode sequence, p k , are included as separate basis functions d m in this database. Thus, Eq. (3.16) is now reformulated as a linear regression problem as follows: 2 6 6 6 6 6 6 6 4 z 1 z 2 . . . z K 3 7 7 7 7 7 7 7 5 | {z } input data 2 6 6 6 6 6 6 6 6 6 6 4 d 1 d 2 d M 3 7 7 7 7 7 7 7 7 7 7 5 | {z } basis functions 2 6 6 6 6 6 6 6 4 b 1 b 2 . . . b M 3 7 7 7 7 7 7 7 5 | {z } coefficients (3.17) or equivalently, z Db: (3.18) Here, D represents the database ofM dierent resolvent mode sequences (or basis functions) and b is the column vector of coecients for each of these basis functions. We seek a sparse solution in which only a small number of coecients are non-zero. In other words, we seek 54 the dominant resolvent mode sequences that best approximate the input. 3.2.3 Identifying Dominant Resolvent Modes Using FROLS The vector of coecients b in Eqs. 3.17 and 3.18 can be computed using the standard least squares method. However, this will result in overtting due to the large number of modes used in the database. Therefore, we seek a sparse solution to this problem where the number of modes used to reconstruct the ow eld (M o ) is much less than the number used in the database (M). The reconstructed ow eld (^ z) can then be represented as: ^ z =b l 1 d l 1 +b l 2 d l 2 + +b l Mo d l Mo (3.19) in whichM o <<M andl m are the indices of the selected modes in the database. The above equation in matrix form is ^ z = 2 6 6 6 6 6 6 6 6 6 6 4 d 1 d 2 d M 3 7 7 7 7 7 7 7 7 7 7 5 | {z } basis functions 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 0 . . . b l 1 0 . . . b l Mo 0 . . . 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 | {z } sparse coeffs: (3.20) Such problems are typically solved using an iterative method that sequentially reduces the complexity of the model [see e.g., 84, 17]. These studies employ standard techniques such as lasso [95] or ridge [44] regression, which solve a least squares regression problem with regu- larization, sequentially eliminating coecients below a certain threshold until the required degree of sparsity is achieved. However, these algorithms initially start with a complex 55 (non-sparse) model and iteratively reduce complexity (or promote sparsity) until an optimal model is reached. If the database is huge, either due to the number or size of the basis functions, the computational cost incurred for iterating can be high. Instead, we use the Forward Regression with Orthogonal Least Squares (FROLS) algorithm, which starts with a simple model and increases complexity (or decreases sparsity) with the number of iterations. This algorithm has previously been used for nonlinear system identication problems [14]. We propose to use a similar approach to identify the dominant resolvent modes and also calibrate their amplitude and phase. The FROLS algorithm is an iterative process with two steps per iteration: 1) Ranking and 2) Gram-Schmidt orthogonalization. At each step a full search of the unselected modes is carried out to pick the best candidate from a large database. Recall that D is the full database of basis functions, D =fd 1 ; d 2 ;:::; d M g, which is composed of M candidate basis functions. We are looking for the subset of these basis functions, D Mo =fd l 1 ; d l 2 ;:::; d l Mo g, that best approximates the training snapshots denoted z so that M o << M. In the rst ranking step, the error reduction ratio, ERR, for each basis function is computed as ERR := z T d m d T m d m 2 d T m d m z T z = z T d m 2 jjzjj 2 jjd m jj 2 (3.21) which is equivalent to `cos 2 ()' between the input data and each of the basis functions. The basis function with the highestERR is selected; this is denoted o 1 , and is identical to d l 1 in D. This basis function is then added to the subset D Mo and a new set of orthogonal basis functions, O. In the next step, which is the Gram-Schmidt orthogonalization, the projection of each of the basis functions in D along d l 1 is subtracted to create a modied database ^ D 1 . The process is then repeated with ^ D 1 to identify the next orthogonal basis function o 2 . This continues until a satisfactory solution to the regression problem is achieved. The ranking and Gram-Schmidt orthogonalization steps are shown schematically in g- ure 3.3. Subplot (a) shows the training snapshot as the red vector and the (non-orthogonal) basis functions as blue vectors in three-dimensional space. The basis function with high- 56 Figure 3.3: Schematic showing successive steps in the FROLS algorithm. Subplot (a) shows the database of basis functions in blue. The training vector (input data) is shown in red and the closest match to the training vector is shown in green. Subplot (b) shows the remaining basis functions after the projection along the green vector in (a) is subtracted. est ERR is shown in green. In the next step, the Gram-Schmidt Orthogonalization, this component is subtracted from rest of the vectors. Now the basis functions span a two di- mensional space in contrast to three dimensions previously as shown in subplot (b). Again the basis function with highestERR is selected, shown as the green vector, and the process is repeated. Thus, we start with D = fd 1 ; d 2 ; d 3 ;:::; d M g and identify the orthogonal set O = fo 1 ; o 2 ; o 3 ;:::o Mo g with M o << M. Note that the orthogonal set O and D Mo are both subsets of D and span the same space. The predictor can be represented in terms of either O or D Mo as ^ z = Mo X i=1 f i o i = Mo X i=1 b l i d l i (3.22) in which the coecients for the orthogonal basis functions, f =ff 1 ;f 2 ; ;f Mo g, can be computed using the orthogonality condition: f i = z T o i o T i o i (3.23) Each of the basis functions in O can be represented as a linear combination of basis functions 57 in D Mo . In the next step, we map back the elements of O to D Mo . Note that the coecients b l =fb l 1 ;b l 2 ; ;b l Mo g can be computed from using the triangle equation Ab l = f; (3.24) where the upper triangular matrix A = 2 6 6 6 6 6 6 6 4 1 a 1;2 a 1;Mo 0 1 a 2;Mo . . . . . . . . . . . . 0 0 1 3 7 7 7 7 7 7 7 5 (3.25) is composed of elementsa r;s = o T r d ls =o T r o r anda s;s = 1. The convergence can be assessed in terms of the reconstruction accuracy for the input data = kz ^ zk kzk : (3.26) The process described above is the standard FROLS algorithm mentioned in Billings [14]. Recall that the real and imaginary parts of the resolvent modes are included in the database as separate basis functions (see Eqs. 3.16 and 3.17). We implement a minor change to the algorithm in the ranking and selection step to select a complete resolvent mode rather than just picking the real or imaginary part. This modication does not alter the underlying concept of the algorithm and is merely an adaptation for the problem at hand. This change is implemented as follows. After ranking each of the basis functions in the database, we select the basis function that has the highest ERR (Eq. 3.21) with the training data and add it to O. Then the rest of the basis functions in D are orthogonalized. Now, instead of ranking the modes in the database again, we pick the complex conjugate of the mode selected in the ranking step. For example, if the selected basis function is the real part of a particular resolvent mode, we select its imaginary part and add it to O. The ranking step is 58 then repeated after orthogonalizing the rest of the basis functions in D. This process ensures that both the real and imaginary parts of a particular resolvent mode are included in the set O. Note that inclusion of both the real and imaginary components of the resolvent modes ensures that arbitrary phase shifts can be captured in the reconstruction. Figure 3.4: (a) Overview of the sparse regression technique. The FROLS algorithm identies and calibrates the dominant resolvent modes from a large database that best represent the measurements. The weight vector is sparse. (b) Reconstructions with 200 resolvent modes in comparison with the DNS ground `truth'. To summarize, gure 3.4 shows a brief overview of the FROLS algorithm described. The measurements from the experiments are used for training. The database is compiled by pooling the resolvent modes computed at dierent k over a span of equation 3.12. The FROLS algorithm then identies the dominant modes that best represent the input training snapshots. 3.2.4 Enforcing Statistical Constraints Although the reconstructed ow eld (^ z = Db) is expected to represent the input data with reasonable accuracy, there is no guarantee that the resulting ow eld will be statistically sound. Thus, we solve an additional optimization problem to enforce statistical consistency. Specically, we minimize the following composite cost function: 59 C =C 1 + (1)C 2 (3.27) where is a hyper-parameter to adjust the t towards the velocity eld or the statistics, and C 1 =kz ^ zk 2 F =kzk 2 F ; C 2 =kz 2 b z 2 k 2 =kz 2 k 2 ; (3.28) represent cost functions for the velocity eld and second-order statistics. Note that the velocity components and statistical quantities used for training can be dierent depending on the data available. As an example, for reconstruction from 2D-2C velocity snapshots in thexy plane, we have z = (u;v) and z 2 = (u 2 ;v 2 ;uv). For this case, second-order statistics for the identied and calibrated resolvent modes are given by b z 2 = 2 6 6 6 6 4 b u 2 b v 2 c uv 3 7 7 7 7 5 = X k j k j 2 2 < 0 B B B B @ 2 6 6 6 6 4 ~ u k ~ u k ~ v k ~ v k ~ u k ~ v k 3 7 7 7 7 5 1 C C C C A (3.29) where ~ u k and ~ v k respectively represent the streamwise and wall-normal components of ~ u k , and a superscript denotes the complex conjugate. Note that second-order statistics for turbulent channel ow, z 2 , can be computed directly from the measurements or obtained from other sources (including models). We use the multistart optimization routine in MATLAB to nd the solution, i.e., the set of coecients =f k g) that minimize the composite cost function in Eq. 3.27: = arg min 0 C( 0 ). As a test of convergence, we systematically increased the number of starting points for the optimization algorithm from 5 to 30. This increase did not yield signicant changes in the error. As a result, for the remainder of this chapter, we present optimization results obtained using only 10 starting points. However, we 60 note that there is no guarantee that the multistart algorithm converges to a global optimum. 3.2.5 Numerical Evaluation The full DNS dataset acquired from JHTDB composes 256 snapshots obtained at intervals of t + = 0:0625, for both the xy and yz planes. Following standard notation, the superscript + denotes normalization with respect to u and. Thus, the total time interval between the rst and last snapshots is T + = 256t + 16. For the xy plane, we use a uniformly sampled grid with N x N y = 32 64 extending over a region of size 2hh. For theyz plane, we use aN y N z = 6432 grid extending over a region of sizeh2h. Finally for the xz plane we use a N x N z = 32 32 grid extending over a region of size 2h 2h Only the rst and last DNS snapshots are used to reconstruct the ow eld. The intervening snapshots are used to quantify reconstruction accuracy using separate error metrics for the observed (trained) variables and the unobserved (watched) variables. For reconstruction from 2D-2C data in the xy plane, the trained components are the velocity elds (u;v) and the statistical quantities (u 2 ;v 2 ;uv). The watched variables are the spanwise velocity component and pressure (w;p), and the statistical quantities (w 2 ;p 2 ;wp). Integrated error metrics for the trained and watched quantities are dened as z(train) = 1 2 ku ^ uk F kuk F + kv ^ vk F kvk F ; z(watch) = 1 2 kw ^ wk F kwk F + kp ^ pk F kpk F ; (3.30) and z 2 (train) = 1 3 ku 2 b u 2 k 2 ku 2 k 2 + kv 2 b v 2 k 2 kv 2 k 2 + kuv c uvk 2 kuvk 2 ! ; z 2 (watch) = 1 3 kw 2 c w 2 k 2 kw 2 k 2 + kp 2 b p 2 k 2 kp 2 k 2 + kup c upk 2 kupk 2 ! : (3.31) 61 Here the variables shown with a ( ^ ) are the reconstructed velocity or statistical components, while the unhatted variables are the DNS ground `truth'. Similarly, for reconstructions from data in both the xy and yz planes, the trained variables are the streamwise and wall- normal velocity components (u;v) for the xy plane, and the wall-normal and spanwise velocity components (v;w) for the yz plane. Statistical constraints are imposed on the quantities (u 2 ,v 2 ,w 2 ,uv). The watched components are (w;p) in thexy plane and (u;p) in the yz plane. The watched statistical quantities are (p 2 ;up). In this case, the error metrics in Eqs. 3.30-3.31 are modied to re ect the dierent trained and watched quantities. Table 3.1: Cases discussed in this study. In the second column, represents resolvent analysis with the molecular viscosity while T is for an eddy viscosity model. N is the number of snapshots used for training. The fourth column indicates the planes and velocity data available for training. is the hyperparameter used in constrained optimization, see Eq. 3.27 Case Visc. Model N Planes 1 2 x-y (u,v) 0.5 2 2 x-y (u,v), y-z (v,w) 0.5 3 T 2 x-y (u,v), y-z (v,w) 0.5 4 2 x-y (u,v), y-z (v,w) 0.9 5 8 x-y (u,v), y-z (v,w) 0.5 6 2 x-y (u,v), x-z (u,w) 0.5 3.2.6 Cases This chapter evaluates reconstruction for the 6 dierent data scenarios listed in Table 3.1. These cases systematically evaluate reconstruction accuracy for dierent measurement planes, viscosity models, training snapshot numbers, as well as the hyper-parameter () that sets the relative weighting for the velocity eld and statistics in the cost function used for the optimization algorithm (Eq. 3.27). 62 Case 1 represents the baseline case in which reconstruction is carried out using the molecular viscosity () resolvent modes from N = 2 snapshots obtained in the xy plane with a time interval of T + 16. Case 2 builds on Case 1 to consider reconstruction from simultaneous measurements in both the xy and yz planes. Case 3 makes use of the eddy viscosity ( T ) resolvent model for reconstruction. Case 4 considers a higher value for the hyper-parameter ( = 0:9 vs. = 0:5 for the remaining cases). This emphasizes reconstruction of the velocity eld. Case 5 evaluates reconstruction accuracy from N = 8 snapshots of the velocity eld. Finally, Case 6 uses simultaneous measurements from xy andxz plane aty + = 200 as shown in gure 3.2 (c). Note that this is in contrast to Cases 2 - 5 that use simultaneous measurements from xy and yz planes. Reconstruction accuracy is evaluated for both the recostructions from FROLS (which does not take into account statistical constraints) and from constrained optimization sep- arately. Since the watched variables and statistics are not used for either the FROLS or optimization, the reconstruction error for these quantities provides a sense for whether the resolvent-based reconstruction framework being developed here is able to generate volumetric reconstructions from limited planar data. 3.3 Results 3.3.1 Temporal Reconstructions Starting with Case 1, i.e., a singlexy plane with two snapshots in time. FROLS is used to identify 200 resolvent modes from a database of 53400 modes. The trianing snapshots can now be represented as linear combination ofD Mo . Figure 3.5 shows the training snapshots from DNS (subplots (a) and (b)) and the reconstruction using 10 modes (subplots (c) and (d)), 50 modes (subplots (e) and (f)), and 200 modes (subplots (g) and (h)). The streamwise velocities are shown on the left half of the gure and the wall-normal velocities are shown on the right half. As the number of modes is increased, the reconstructed snapshots start 63 to look more like the training snapshots. At 200 modes of reconstruction, the reconstructed velocity eld match the training data both qualitatively and quantitatively. This observation is consistent with the error of reconstruction computed using equation. 3.26 and is shown in gure 3.6. The error decreases monotonically with number of resolvent modes used and the rate of decrease of this error decreases with number of modes, In other words the error curve shows asymptotic nature and using any additional modes would not signicantly reduce this error. In the interest of time required for computation, we use only 200 resolvent modes. Although the training snapshots are reconstructed with sucient accuracy, it is not nec- essary that the reconstructed ow eld satises the statistical constraints in equation 3.29. Therefore the coecients of the identied modes from FROLS are now optimized to mini- mize the cost function in equation 3.27. Figure 3.7 shows the reconstructed statistics using equation 3.29. The blue line shows the statistics from FROLS reconstruction. The red dots correspond to the statistics after optimization. The DNS ground truth is shown in red dashed line. Notice that the plot has two axes. The left axis shown in blue is O(10 5 ), where as the right axis is shown in red is ofO(1). The statistics from FROLS reconstruction is several orders of magnitude larger than the DNS `truth'. After optimization, the recon- structed statistics qualitatively and quantitatively match the DNS data. The left hand side of Figure 3.8 show the error from the reconstructed velocity eld with time from FROLS reconstuction, after optimization, and linear interpolation between the snapshots as a base- line comparison. These errors are again computed using equation 3.26. The error from FROLS is almost constant at 0:35 with a slight waviness. The error after optimization increases to 0:42. The PIV snapshots are available at the beginning and end i.e., t + = 0 andt + 16. The intervening time is the reconstructed ow eld which is compared to DNS data to quantify the error using equation 3.26. Note that the resolvent modes are travelling wave solutions. When the ow eld is advanced forwards in time the modes are advected with their respective wave speed and results in constant error with time. The right hand side of this gure shows the reconstructed snapshots of streamwise and wall-normal velocity 64 Figure 3.5: Flow elds of stream wise and wall-normal uctuations of velocity. subplots (a) and (b) are for DNS and subplots (c)-(h) are for ow elds reconstructed by increasing the number of modes. in comparison with the DNS truth at the middle of the time horizon. The reconstructed snapshots qualitatively as well as quantitatively match the DNS ow eld. Although results here are shown for Case 1, the method of temporal reconstruction remains the same. The results are qualitatively similar for Cases 1 - 6. These trends are: error decreases with increasing number of modes for the training snapshots, reconstructed statistics from FROLS do not match the DNS statistics even qualitatively, however, after 65 0 50 100 150 200 No. of Modes 0.3 0.4 0.5 0.6 0.7 0.8 Figure 3.6: Error computed using equation 3.26 is plotted with increasing number of resolvent modes. The error decreases steeply in the beginning but later it asymptotes. 0 500 1000 0 1 2 3 4 10 5 0 2 4 6 8 10 FROLS Statistics enforced DNS 0 500 1000 0 0.5 1 1.5 2 2.5 10 4 0 0.5 1 1.5 0 500 1000 -10 -8 -6 -4 -2 0 2 10 4 -1.5 -1 -0.5 0 0.5 Figure 3.7: Reconstructed statistics from FROLS (blue line), after constrianed optimization (red dots), and the statistics from DNS (red dashed line) are shown for three statistical quantities for Case 1. The statistics from FROLS are several orders of magnitude larger and are shown on the left axis (blue). The statistics after constrained optimization match the DNS statistics closely. optimization the statistics match the DNS data, and the error of reconstructed ow eld after optimization is higher than that from FROLS. 3.3.2 Estimation of Unobserved Variables Table 3.2 shows reconstruction errors in the trained and watched variables for the training snapshots, along with the errors for the trained and watched statistical quantities. For each case, the top row shows the error from FROLS. The second row (in parenthesis) shows the error after constrained optimization. Starting with the baseline case, where a single 66 Figure 3.8: The error of reconstruction from FROLS, after constrained optimization, and linear interpolation between the snapshots is shown in gure on the left. The gure on the right shown the reconstructed ow eld for u and v components along with the DNS truth. The error after optimization is slightly higher than the FROLS reconstructions and the ow elds match qualitativley and quantitatively with the the DNS. 0 5 10 15 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Case 1 Case 2 Case 3 Case 4 Case 5 Case 6 0 5 10 15 0.5 1 1.5 2 2.5 a) b) Figure 3.9: Error evolution for trained variables in the xy plane (a) and yz plane (b). Solid lines represent reconstructions from FROLS and dashed lines from constrained optimization. xy plane is used, the errors for the training snapshots from FROLS are reasonably low ( z(train) < 0:4) and the ow elds, not shown here for brevity, are qualitatively similar to the DNS truth. However, the errors for the statistical quantities are several orders of magnitude higher,O(10 4 ). After optimization the error in the statistics decreases to O(10 1 ) with a slight increase from 0:36 to 0:46 in the error for the trained components of the velocity eld. 67 Table 3.2: Reconstruction errors in the velocity eld and statistical quantities for the trained and watched variables. For each case, the top row shows reconstruction error after the FROLS algorithm is used to identify and calibrate 200 resolvent mode amplitudes based on the training snapshots. The bottom row (in parentheses) shows reconstruction error after the optimization algorithm modies the amplitudes to enforce statistical constraints. Case z(train) z(watch) z 2 (train) z 2 (watch) 1 0:36 760:71 2:82 10 4 1:23 10 4 (0:46) (2:07) (0:15) (1:14) 2 0:51 1:15 12:50 224:16 (0:60) (1:25) (0:30) (1:54) 3 0:56 2:16 6:78 10 9 2:51 10 10 (0:86) (0:94) (0:18) (0:81) 4 0:51 1:15 12:50 224:16 (0:54) (1:14) (0:49) (1:78) 5 0:66 1:00 12:50 224:16 (0:75) (1:01) (0:31) (1:79) 6 0:52 1:10 2:18 6:61 (0:71) (1:56) (0:58) (2:60) It is interesting to note that the error for the watched variables (w;p) decreases substantially, from 760 to 2. Building on this case, we next use simultaneous measurements from both xy and yz planes. The error evolution for the reconstructions in the xy and yz plane are shown in Fig. 3.9(a) and (b), respectively. Reconstruction error for Case 2 in xy plane varies little over the time horizon, with a mean value ( z(train) 0:5) which is very close to the values reported for the training snapshots in Table 3.2. After optimization, the errors for both the trained and watched variables increases slightly, while the errors for the statistical quantities again decrease substantially. For instance, z 2 (train) decreases from 12:5 to 0:30. It is important to note that the use of two orthogonal planes results in a signicant decrease 68 in the error for watched quantities even before optimization. The error for the watched components of velocity eld changes from 760 for Case 1 to 1:15 for Case 2. Similarly, reconstruction error for the statistical quantities decreases by over two orders of magnitude for Case 2 relative to Case 1. When the T resolvent model is used (Case 3), reconstruction errors for the velocity eld from FROLS are comparable to those for the resolvent model (Case 2). However, errors for the statistical quantities (O(10 10 )) are several orders of magnitude larger than that for Case 2 ( O(10 4 )). After optimization, the errors for the statistical quantities are reduced to O(10 1 ). However, the error for the velocity eld rises signicantly when compared to Case 2. Figure 3.10 compares reconstructed statistics from FROLS for Case 2 and the constrained optimization results for Cases 2 and 3 with the DNS `truth'. As discussed before, the statistics from FROLS are several orders of magnitude higher, while the reconstructed statistics after optimization qualitatively match the DNS proles. Note that the statistics corresponding to Case 3 are smoother and more closely approximate the DNS proles compared to Case 2. This is consistent with the observations of Morra et al. [67], who showed that the T resolvent model yields better reduced representations for turbulence statistics. This is because the resolvent modes computed with an eddy viscosity model have greater support in the wall-normal direction. This results in smoother statistics. However, this also smooths out ne-scale features in the reconstructions which result in larger errors for the velocity elds. Figure 3.11 compares the reconstructed ow elds after optimization for Case 2 with the corresponding snapshots from DNS. Subplots (a)-(b) show the trained velocity compo- nents in the xy plane (u;v) for one of the 2 training snapshots. Subplots (c)-(d) show the watched variables in the xy plane, i.e., spanwise velocity and pressure. Subplots (e)- (f) show the velocity elds reconstructed for a xy plane that is oset by z=h = 0:125 in the spanwise direction from the training snapshots. In other words, panels (c)-(f) show reconstructions for quantities that are not available for training the model. The reconstruc- 69 0 500 1000 0 10 20 30 0 2 4 6 8 10 FROLS Opt. (Case 2) Opt. (Case 3) DNS 0 500 1000 0 10 20 30 0 0.5 1 1.5 2 0 500 1000 -10 -8 -6 -4 -2 0 -1 -0.5 0 0.5 0 500 1000 0 20 40 60 80 100 0 0.5 1 1.5 2 2.5 0 500 1000 0 20 40 60 80 0 5 10 15 0 500 1000 -10 -5 0 5 10 -2 0 2 4 6 Figure 3.10: Reconstructed statistics from FROLS for Case 2 and constrained optimization for Cases 2 and 3. The blue scale on the left corresponds to the FROLS reconstructions with the orange scale on the right corresponds to the optimization results and DNS proles. Note the substantial dierence in magnitude from the FROLS proles to the optimization and DNS proles. The statistics from eddy viscosity model are smoother and more closely approximate the DNS proles. tions are labelled with a ( ^ ). As expected, reconstructions for the training components (u;v) show good qualitative agreement with the DNS snapshots, though some ne-scale features are missing (see Fig. 3.11(a,b)). We do not see signicant similarity in structure for the untrained components, i.e., (w;p) in Fig. 3.11(c,d). However, the magnitudes are captured. Reconstructions in the spanwise oset plane (Fig. 3.11(e,f)) again do not match the DNS data but the magnitudes are within range. A comparison of Cases 1 and 2 shows that reconstruction error for the trained velocity components increases when using data from two orthogonal planes. However, errors for the watched variables and statistical quantities are reduced signicantly in comparison to those obtained using data from a single plane. This suggests that using simultaneous measurements in the xy and yz planes may help in estimating unobserved (untrained) variables and 70 potentially reconstructing volumetric ow elds from limited planar data. Figure 3.11: Reconstructed velocity and pressure elds after optimization for Case 2 are compared with DNS data. Subplots (a)-(b) show reconstructions for the observed (trained) components of velocity in the xy plane. Subplots (c)-(d) show reconstruction for the unobserved (watched) velocity and pressure elds. Subplots (e)-(f) show reconstructions at a dierent xy plane (z=h = 0:125) for the trained components (u;v) Next, we consider the in uence of the hyper-parameter () used in Eq. 3.27 on errors for the velocity eld and statistical quantities. Table 3.2 shows that, when compared to Case 2 ( = 0:5), reconstruction errors for the velocity eld for Case 4 ( = 0:9) are slightly better after optimization. However, errors for the statistical quantities are slightly worse. A similar trend can be seen for both the trained and watched variables. This is expected from Eq. 3.27. A higher value for places greater emphasis on the loss from the ow eld, which improves the velocity reconstructions but compromises the statistics. Case 5 considers the eect of using more training snapshots. We still use 2 snapshots for training the FROLS reconstructions. However, after the set of 200 resolvent modes is 71 identied, the data sequence is augmented with 6 additional snapshots obtained at the same spatio-temporal resolution. The coecients are then optimized to minimize the loss from these 8 snapshots and the statistics. Figure 3.9 shows that the error after optimization is signicantly increased for thexy plane when reconstruction is carried out from 8 snapshots. For the orthogonal (yz) plane, reconstruction error only increases slightly (see the change in dotted vs. dashed lines). Note that the error shown in Table 3.2 for Case 5 considers all 8 snapshots in the data sequence. However, the time series in gure. 3.9 only shows the reconstruction between the rst 2 snapshots. Finally, Case 6 considers the use of simultaneous measurements in xy and an xz plane at y + = 200. The results are qualitatively similar to Case 2. To summarize, the error for the train variables increases slightly after optimization, the error for the watch variables are ofO(1) even before optimization and after optimization they increase slightly. The error for the trained and watched statistical components are also of O(1) and again these errors decrease after optimization. In contrast to Case 2, we use statistical constraints only on u 2 ;v 2 ;uv, because only these statistical quantities can be computed over the channel height from the available measurements. The time evolution of errors is shown in gure 3.9. It is interesting to note that the error for the orthogonal plane is the lowest among the other cases. This is because the ow structures in the xz plane are advected in the streamwise direction with the wavespeed identied by FROLS. Whereas the other orthogonal plane (yz) is perpendicular to the streamwise direction. In general, errors for reconstructions from the xy plane are lower in comparison to errors for the yz plane. This is due to the fact that turbulent channel ow is advection- dominated. Since the ow advects in thexdirection, having training snapshots in thexy plane helps identify appropriate traveling wave resolvent modes. When interpolated in time using Eq. 3.13, these resolvent modes propagate downstream with the phase speedc =!=k x . While the snapshots in theyz plane provide some insight into spanwise structure, they may not contain the information necessary to characterize oblique wave motion. This indicates 72 that measurements made in the horizontal xz plane may yield improved out-of-plane reconstructions. Note that reconstruction errors generally increase with increasing input data (see Cases 1, 2, and 5). This may be attributed to the fact that all reconstructions are carried out using just 200 resolvent modes. A larger set of resolvent modes may be necessary to accurately represent the input data from multiple planes and snapshots. 3.4 Chapter conclusions The results presented in this chapter show that projection onto resolvent modes can yield useful reconstructions for wall-bounded turbulent ows when very limited measurements are available. The FROLS algorithm eciently identies resolvent modes that best represent the data from a large database. However, additional constraints are needed to ensure that the reconstructions are statistically sound. Errors from FROLS reconstructions are quite high for unobserved (watched) variables and statistical quantities when 2D-2C data from a single xy plane are used as input data. After optimization to minimize a composite loss function that accounts for the input data as well as second-order statistics, these errors decrease signicantly. The use of simultaneous measurements in multiple orthogonal planes also reduces this error. The use of an eddy viscosity resolvent model leads to an increase in reconstruction error for the input velocity snapshots but greatly improves the statistics. In general, reconstruction error for unobserved quantities (e.g., out of plane velocity, pres- sure) reduces as statistical constraints are introduced or simultaneous measurements from multiple orthogonal planes are used. Future work will consider reconstruction from input data consistent with stereoscopic and tomographic PIV (i.e., 2D-3C and 3D-3C) as well as wall-based measurements. 73 Chapter 4 Reconstruction and estimation for stratied wakes 4.1 Introduction 4.1.1 Motivation and Problem Statement Wakes behind blu bodies in natural environments frequently evolve in the presence of stable stratication in uid density. Such stratied wakes can be long lived and coherent [90]. The topology of these wakes can contain valuable information regarding the upstream blu body as well as the ow regime. Indeed, many aquatic organisms are thought to use the cues from wake signatures for navigation. For instance, it has been proposed that Brazilian green sea turtles make use of wake signatures from islands to reliably make a round trip of approximately 4400 km [77]. Fish, seals, and copepods also use such navigational cues for mate nding, predator evasion, and hunting [66, 81]. In general, the ability to identify the wake creator from limited measurements in the wake could aid object identication, navigation, and motion planning eorts in environments with, and without stratication. Stratied wakes are in uenced by inertial, viscous, and buoyancy forces. Hence, the relevant non-dimensional parameters governing these ows are the Reynolds number, Re, 74 which is the ratio between the inertial and viscous forces, and the Froude number,Fr, which is the ratio between the convective time scale and the buoyancy time scale. Perturbations to a stably stratied uid give rise to internal waves of characteristic frequency N, which is called the Brunt-V ais al a frequency. This characteristic frequency is dened asN = q g o @ @z , whereg is the acceleration due to gravity, is the uid density which decreases in the vertical (z) direction, and o is a mean density. Compared to non-stratied wakes, buoyancy eects in stratied wakes give rise to a rich variety of topological structures [55, 24]. Stable density stratication suppresses turbulent uctuations in the vertical direction [54, 89, 16, 80], while simultaneously enhancing the possibility of vertical momentum and energy transport through the internal waves excited either by the body itself, or by the wake. It has also been argued that late wakes from streamlined and sharp-edged bodies evolve in a self-similar manner when scaled appropriately, which suggests that information regarding the initial conditions is lost [62]. However, it is possible this loss of information is limited to statistical quantities and that the wake topology may yet retain information regarding its origins [90]. Several previous experimental and numerical studies have characterized the topology of stratied wakes as a function ofRe andFr [43, 55, 24, 74, 76, 25, 75]. Recent work by Ohh and Spedding [73] and Madison, Xiang, and Spedding [57] proposes that stratied sphere wakes can be classied into ve principal regimes across (Re;Fr) space. The Reynolds number is Re = U o D= and the Froude number is dened (based on radius) as Fr = U o =(NR). Here, U o is the freestream velocity, D is sphere diameter, R is the radius, and is the kinematic viscosity of the uid. As shown in gure 4.1, these ve regimes are labelled: vortex street, symmetric non-oscillation, asymmetric non-oscillation, planar oscillation, and spiral mode, describing the qualitatively dierent modes that emerge from the near-wake. Figure 4.1 includes the regime boundaries identied by Lin et al. [55] and Chomaz, Bonneton, and Hopnger [24] based on the body and near-wake ow features. Though the regimes delineated in these studies involved a dense coverage of parameter space that implied sharp transitions between regimes, the boundaries and regimes themselves are not identical and 75 Figure 4.1: Regime space in terms of Reynolds number (Re) and Froude number (Fr). Filled circles show the cases used in this study. Physical regimes identied in previous studies are denoted by the color of the circles [57, 73] and the background shading [55, 24]. Regimes are color-coded as: vortex street (purple), symmetric non-oscillation (red), asymmetric non- oscillation (yellow), planar oscillation (green), and spiral mode (blue). Overlapping shaded regions re ect dierences in regime boundaries identied in previous studies. we remain agnostic on whether or not transitions are sharp. Here the (Re;Fr) space is populated by experimental and numerical data at 24 discrete combinations reported in [73], and shown using lled circles in gure 4.1. Topological distinctions between the dierent regimes are illustrated in gure 4.2. For each regime, gure 4.2 shows contours of streamwise velocity in the horizontal (xy) and vertical (xz) planes obtained in DNS, along with cross-stream velocity components for the dominant DMD mode (see Sec. 4.2 for further detail on the simulations and DMD). The vortex street regime is illustrated in gure 4.2(a). Buoyancy forces due to stratication are the strongest in this regime (Fr = 0:5 is the lowest value) and excursions in the vertical plane are strongly suppressed. The dominant DMD mode in gure 4.2(a) shows the presence of a coherent vortex street in the horizontal plane while the slice in the vertical plane shows low amplitude traces of internal waves. WhenFr 1, the lee wave is steady in the body reference 76 20 25 30 -4 0 4 -4 0 4 0 5 10 15 20 25 30 -4 0 4 -4 0 4 20 25 30 -4 0 4 -4 0 4 0 5 10 15 20 25 30 -4 0 4 -4 0 4 20 25 30 -4 0 4 -4 0 4 0 5 10 15 20 25 30 -4 0 4 -4 0 4 -4 0 4 -4 0 4 20 25 30 -4 0 4 0 5 10 15 20 25 30 -4 0 4 -4 0 4 0 5 10 15 20 25 30 -4 0 4 -4 0 4 20 25 30 -4 0 4 (a) (b) (c) (d) (e) DMD mode Full-field Figure 4.2: Snapshots of normalized streamwise velocity (u) from DNS (left column) and cross-stream velocity components (v, w) for the dominant DMD mode (right column) in dierent regimes: (a) vortex street at (Re;Fr) = (1000; 0:5) , (b) symmetric non-oscillation at (Re;Fr) = (500; 1), (c) asymmetric non-oscillation at (Re;Fr) = (300; 4), (d) planar oscillation at (Re;Fr) = (500; 4), and (e) spiral mode at (Re;Fr) = (1000; 8). Shading for the DMD modes shows positive (red) and negative (blue) values for v andw, normalized by 0:75max(jvj;jwj). The Strouhal number, St =fD=U o , for the oscillation frequency f of the DMD modes is provided. 77 frame and dominates the near wake. Fluctuations are suppressed in both the horizontal and vertical planes, which results in a time-invariant ow eld that is classied as symmetric non-oscillation (gure 4.2(b)). With further increases in Fr, an asymmetry appears about the vertical plane for Re = 300. This regime, termed asymmetric non-oscillation, is shown in gure 4.2(c). Though only one instance of this regime is observed in the dataset being used here, at (Re;Fr) = (300; 4) it is a consistent feature in both laboratory and numerical experiment. As Fr increases further, the wake contains a mix of horizontal oscillations and vertical motions associated with internal waves (gure 4.2(d)). This regime is termed planar oscillation because most kinetic energy is associated with the horizontal oscillatory motion. These oscillations eventually appear in the vertical plane, replacing the internal wave modes, as Re and Fr rise, and the wake becomes fully three dimensional. This regime is termed the spiral mode (gure 4.2(e)) after the known shedding mode that appears at all higher Re Chomaz, Bonneton, and Hopnger [24]. The classication system described above is based on known and observed physical mechanisms, but it relies on topological distinctions identied and delineated by humans, who design the classication system. Here, we develop an expert-free and purely data-driven approach to identifying the dynamical regime (as gauged by an estimated Re and Fr pair) from sparse measurements in a stratied wake. This approach assumes the availability of signicant prior measurement or simulation data across the relevant (Re;Fr) space. As detailed below, regime identication then requires solution of a sparse regression problem that seeks to describe the wake structures present in the input (from an unknown regime) in terms of the wake structures present in the large database of prior measurements or sim- ulations (from known regimes). This approach also enables spatio-temporal reconstruction of the wake from limited data. 78 4.1.2 Previous Data-Driven Classication and Reconstruction Ef- forts There have been many prior classication and reconstruction eorts for unstratied wakes behind blu bodies. For instance, such wakes are often classied into dierent categories based on the number of alternating vortices shed per oscillation period [108]. For example, 2S refers to 2 single vortices shed per period of oscillation, 2P refers to 2 pairs of vortices of opposite signs shed period, S+P refers to an asymmetric version where a single vortex and a pair of vortices are shed per period, and so on. Similar to the problem being tackled in the present study, Colvert, Alsalman, and Kanso [27] classied the wake behind a apping airfoil into dierent regimes (2S, 2P+2S, etc.) from point measurements of vorticity using neural networks. In a continuing study, Alsalman, Colvert, and Kanso [5] evaluated four dierent kinds of ow sensors to classify the ow using neural networks and concluded that the use of a sensor that measures transverse velocity led to more robust classication. While pre-trained neural networks perform admirably for classication purposes, they have some important limitations. Specically, their utility is limited to the type of input data used for training purposes. For instance, a network trained using velocity snapshots cannot be used for point measurements and a network trained using point measurement data may not be used for velocity snapshots. In addition, such networks have limited ability to reconstruct velocity elds or generate predictions for input signals from regimes not included in the training data. The framework developed in the present study begins to address some of these limitations. For completeness we note that there have been many prior studies on super-resolution reconstruction, ow eld estimation from wall measurements, and state estimation from limited and non-intrusive measurements. Such studies have employed physics-based models [e.g. 36, 11, 47, 51, 96, 105, 104], data-driven modal representations [e.g. 18, 107, 99, 98], as well as machine learning tools [e.g. 35, 33, 34, 5, 40, 41]. A detailed review of these 79 eorts is not provided here for brevity. Recently, Gra et al. [38] used an approach similar to the present study to build dynamic models for Poiseuille ow using simulation data, and for a rotating n using PIV data. Specically, Gra et al. [38] extended the Least Angle Regression algorithm to select DMD modes from a library that best t the input data, and used this algorithm to develop reduced-complexity dynamic models. However, we note that these prior eorts on reconstruction, estimation, and the development of reduced-complexity dynamic models do not typically target regime or parameter identication. For stratied wakes, a number of experimental [55, 24] and computational [76, 75] stud- ies have made detailed descriptions of body-generated wake structures as a function of Re and Fr. In recent years, modal analysis techniques such as proper orthogonal decomposi- tion (POD) and dynamic mode decomposition (DMD) have been used to extract coherent structures and to build reduced order models [30, 110]. Snapshot POD yields spatially or- thogonal modes ordered based on energy, with temporal coecients that may contain a mix of frequencies [92]. To isolate modes of specic frequency, modal analysis methods that yield temporally orthogonal modes such as DMD [20, 87] or spatially and temporally orthogonal modes such as Spectral POD [97] can be used. Previous eorts have successfully imple- mented both snapshot POD and DMD to extract coherent structures in stratied wakes [73, 71]. Nidhan et al. [71] used DMD to analyze the wake structures from numerical datasets at two dierent Reynolds numbers with varying buoyancy eects (Fr). More recently, Ohh and Spedding [73] used DMD to classify wake structures at variousRe andFr from both numer- ical and experimental datasets. A custom algorithm informed by the physically-motivated classication system described in the previous section was developed to characterize wake topology. This algorithm interrogated the velocity eld, vorticity eld, and oscillation fre- quency associated with the most energetic DMD mode and made use of an expert-dened decision tree to classify wake topology into one of the ve dierent physical regimes shown in gure 4.1. To our knowledge, there are no previous works that directly estimate Re and Fr for stratied wakes from limited measurements. 80 4.1.3 Contribution and Outline In this study, we attempt regime identication and ow eld reconstruction from limited measurements in stratied wakes behind a sphere. We use a purely data-driven approach in which the input data (i.e., a small number of velocity snapshots or point measurements) are represented in terms of DMD modes obtained from numerical simulations across a wide range of Re and Fr. These DMD modes are pooled into a large database. We then seek a sparse representation of the input signal in terms of a linear combination of modes from this database. The forward regression with orthogonal least squares algorithm [FROLS 14] identies the dominant DMD modes that best represent the data. In other words, only a few dominant DMD modes are selected to represent the input data. The parameter regime for the input data is then estimated using a weighted average of Re and Fr values for the DMD modes selected. We also evaluate the predictive capability of the ow eld model built using this process and introduce a metric to evaluate condence in the predictions. There are a few important distinctions between the present study and previous eorts on ow reconstruction and classication. Flow reconstruction eorts typically use a library of POD or DMD modes from the same parameter regime as the measurements to reconstruct the velocity eld in space or time. Here, we use a library of DMD modes from various known regimes to reconstruct the velocity eld using limited input data from an unknown regime. In so doing, we can generate quantitative estimates for Re andFr. Unlike previous classication eorts, this approach also allows for interpolation in parameter values rather than binning into discrete regimes. We also note that, though the present eort focuses on stratied wakes, the regime identication framework can be generalized to other ows. The remainder of the chapter is organized as follows. In Sec. 4.2, we describe the simulations and experiments that generated the data used in this study together with the FROLS algorithm used for sparse regression. We also discuss the approach used to identify the parameter regime, and the metrics used to evaluate prediction accuracy and condence. Note that we attempt regime identication using a variety of input data from DNS and 81 laboratory experiments. This includes planar velocity snapshots (c.f., 2D measurements from Particle Image Velocimetry, PIV) as well as point measurements. In Sec. 4.3, we rst present the results for the baseline case in which planar velocity elds from DNS are used as the input data and regime identication is carried out using a comprehensive DMD mode database. We then use a sub-sampled database of DMD modes (i.e., a database that has only half the available Re and Fr cases) to evaluate predictive capability and test whether the regime identication framework is capable of interpolation. Next, we attempt regime identication with real-world experimental data as the input. We rst use snapshots of the velocity eld obtained from PIV for regime identication, and we evaluate predictive capability using dierent components of velocity. Next, we consider regime identication from point measurements also attempt to identify sensor location. In other words, the framework is used to not only identify Re and Fr based on the DMD modes selected, but also the physical placement of the sensors on the spatial DMD modes. We also evaluate the eect of database size for input data from 2D2C snapshots and point measurements. Finally, in Sec. 4.4 we present a brief summary of our ndings and conclu ding remarks. 4.2 Methods In this section, we rst describe the simulations used to generate the DMD database and evaluate the regime identication algorithm. We also describe the laboratory experiments used to generate additional datasets used for classication purposes. We subsequently present the problem formulation and the FROLS algorithm used to nd a sparse representation for the input data. Note that this eort makes use of the same experimental and numerical dataset as Ohh and Spedding [73], who developed a physics-informed algorithm to classify stratied wakes into the dierent physical regimes shown in gure 4.1. 82 ͳ ͸ ܦ ͳ Ͳ ܦ ͷ Ͳ ܦ 1/4 ܦ 1/8 ܦ 1/16 ܦ 1/32 ܦ ͳ ͸ ܦ ͳ Ͳ ܦ ͷ Ͳ ܦ 1/4 ܦ 1/8 ܦ 1/16 ܦ 1/32 ܦ 1/64 ܦ 1/128 ܦ (a) (b) Figure 4.3: Meshes used for numerical simulations (a) at Re =f200; 300; 500g and (b) at Re = 1000. The mesh used for the lower Re cases has four renement levels (2D [25] D). The mesh used for the higher Re case has six renement levels (2D [27] ). 4.2.1 Numerical Simulations and DMD Database Direct numerical simulations of stratied wakes behind spheres were carried out systemat- ically over a wide span of the (Re;Fr) parameter space. Four dierent Reynolds numbers, Re = U o D= =f200; 300; 500; 1000g, and six dierent Froude numbers, Fr = U o =(NR) = f0:5; 1; 2; 4; 8; 16g, were considered for a total of twenty four dierent cases (see lled circles in gure 4.1). Following Ohh and Spedding [73] and Madison, Xiang, and Spedding [57], stratied sphere wake ows in this parameter space can be classied into 5 broad categories, as denoted by the shaded regions in gure 4.1. The numerical simulations made use of the Boussinesq approximation to the Navier- Stokes equation together with a density transport equation. The dimensionless form of the governing equations is r u = 0; (4.1) @u @t + uru =rp 0 + 1 Re r 2 u + 1 Fr 2 0 g; (4.2) @ 0 @t + ur 0 =w d dz + 1 PrRe r 2 0 ; (4.3) where u is the velocity eld, t is time, g =ge z is the gravitational acceleration in the vertical (z) direction, w is the vertical velocity, and Pr is the Prandtl number. The mean background density eld is denoted and the uctuating density eld is 0 . Similarly, p 0 represents the pressure uctuations about the background hydrostatic pressure. The governing equations were solved using the PimpleFOAM solver in the OpenFOAM suite, 83 which uses a nite volume method with a pressure-velocity coupling algorithm. Second- order accurate temporal and spatial discretization schemes were used. To reduce internal wave re ection from the wall, the y-normal and z-normal walls of the domain were set to zero-gradient conditions. To allow for possible breaks in the axisymmetry of the wake, the sphere was oscillated once along each of the three axis as an initial condition. The extent of the computational domain was [10D; 50D] in the streamwise (x) direction, [8D; 8D] in the lateral (y) direction, and [8D; 8D] in the vertical direction. The positionx = 0,y = 0, z = 0 corresponds to the center of the sphere. Multiple grid renement levels were used in an unstructured mesh, where the cell size was halved at each level closer to the sphere, as shown in gure 4.3. The resolution of the coarse mesh used for Re =f200; 300; 500g was roughly 3 million cells while the ne mesh used for Re = 1000 consisted of roughly 17 million cells. Velocity elds were sampled at a convective time increment of t DNS =R=U o for Re =f200; 300; 500g and t DNS = 0:4R=U o for Re = 1000. The number of snapshots saved ranged between 100-1250 depending on the complexity of the ow. However, in all cases snapshots were saved over a minimum of ve dominant oscillation periods. For further computational details, the reader is referred to Ohh and Spedding [73]. For each of then r = 24 cases in the (Re;Fr) parameter space, dynamic mode decomposi- tion was carried out using data obtained over the regionx2 [20:25R; 29:75R],y2 [4R; 4R], andz2 [4R; 4R] with a uniform grid spacing of (1=8)R in each direction. This region was chosen as the observation window since it is still within the near-wake but does not include developments immediately behind the sphere (i.e., the recirculation zone). DMD was car- ried out using the method described in Chen, Tu, and Rowley [20] from the time-resolved sequence of S + 1 velocity eld snapshots from DNS,fu 0 ; u 1 ; ; u S g, which generated a set of DMD spatial modes (or Ritz vectors)fv 1 ; v 2 ; ; v S g and the corresponding tempo- ral coecients (or Ritz values)f 1 ; 2 ; ; S g. Since DMD is an established technique, the details of this procedure are not repeated here. For the lower Reynolds number cases (Re =f200; 300; 500g), DMD was carried out using a sequence of S + 1 = 40 snapshots 84 from DNS obtained at time intervals of t DNS = R=U o . Thus, the total duration of the sequence was 39(R=U o ), which is roughly 4 ow through times for the streamwise extent of the snapshots used for DMD (9:5R). The use of 40 snapshots generatedS = 39 DMD modes for each case. The DNS time step forRe = 1000 cases was lower, and the velocity snapshots were available at t DNS = 0:4R=U o . Therefore we used a total ofS + 1 = 100 snapshots for DMD to keep the total duration constant at roughly 4 ow through times. This resulted in 99 DMD modes for each ofRe = 1000 cases. It is important to note that t DNS is the time dierence between consecutive snapshots in the DNS database and not the time step used in the DNS. For eachReFr combination, the rst 35 modes sorted in decreasing order of energy were retained and used to construct a library (or database) of DMD modes for sparse regression. We used an equal number of modes for each case to eliminate any regression bias due to inclusion of more modes from a particular parameter combination. Note that the time-evolution of the snapshots can be expressed in terms of the complex Ritz vectors and values as: u k = S X j=1 k j v j ; k = 0; ;S 1: (4.4a) u S = S X j=1 S j v j + r; r? spanfu 1 ; u 2 ; ; u S1 g (4.4b) where r is the resulting residual for the last snapshot. Since the Ritz values determine the time evolution from one snapshot to the next, the Strouhal number for the j th DMD mode is given by St j := ! j R 2U o = arg( j )R 2t DNS U o (4.5) where t DNS is the time interval between two consecutive snapshots in the DNS database and! j represents the oscillation frequency. Note that the snapshots used here are the three- dimensional (3D) ow elds from DNS and so the spatial modes obtained via DMD are also 3D ow elds with the same spatial extent. These modes must therefore be sub-sampled to 85 50D ͳ ʹ ܦ ʹ ʹ ܦ U Screw driven translation stage ͳ Ͳ ܦ ͳ Ͳ ܦ ͳ Ͳ ܦ PIV FOV ݔ ݕ ݖ U ݐ ௡ A A A A A A B (a) (b) Figure 4.4: (a) Schematic showing experimental setup. The annotated scale is based on a typical sphere diameter of D = 4 cm. (b) Schematic showing temporal data reconguration where the colored box A is moving in reference to the sphere. yield 2D snapshots or point measurements when the DMD mode database is compiled for regression and regime identication. Further details on how the database is compiled for dierent input data are provided below in Sec.4.2.3. 4.2.2 Laboratory Experiments Experiments were carried out at values ofRe andFr that match the simulations (gure 4.1). A sphere was towed horizontally in a tank of size 1m1m2:5m, lled with stably stratied uid as illustrated in gure 4.4(a). To achieve a linear density gradient while avoiding optical distortions, a salt and ethanol water mixture was used to match refractive index [111]. The sphere consisted of a 3D-printed polylactic acid (PLA) shell around a steel ball. This sphere was submerged to the center of the tank in the lateral and vertical directions with three shing wires (d = 0:5mm). To adjust Re and Fr, the buoyancy frequency, sphere radius and tow speed were varied systematically over the following ranges: N2 [0:13; 1:0] rad/s, R2 [0:72; 5:5] cm, and U o = [0:37 43] cm/s. Planar velocity measurements were obtained in the wake of the sphere using Particle Image Velocimetry (PIV). The tank was seeded with titanium dioxide particles of aver- age density 4:23 g/cm 3 and diameter 15 m. The vertical and horizontal center planes 86 were illuminated with a pulsed laser (Nd:YAG, LaVision NANO L100-50PIV). Two cam- eras (LaVision-Imager sCMOS), each having a resolution of 2560 2160 pixels, were aligned side-by-side with 5% overlap to obtain an extended eld of view (FOV) in the streamwise direction. For measurements in the vertical (xz) plane, the FOV for each camera was approximately 350 mm 300 mm. Given the varying sphere sizes, this translates into a scaled FOV that varies from 6:5R 5:5R to 49R 42R. For measurements in the horizontal (xy) plane, the FOV for each camera was roughly 230 mm 190 mm, which translates into a scaled eld of view ranging from 4:2R 3:5R to 32R 26R. Images were taken at a single-pulsed time interval ranging between t PIV = 0:050:2 s, which was chosen to ensure that the particle pixel displacements between images produced a well-resolved velocity eld. To ensure consistency with the numerical simulations, the experimental data were trans- lated to a sphere-xed frame. A xed sphere reference frame was achieved by stacking half the streamwise span of the FOV, markedA in gure 4.4(b), into a temporal sequence. Note that, box A moved with the sphere as it is towed across the FOV. Once box A reaches the upstream edge of the FOV, the next half of the FOV (marked B in gure 4.4(b)) becomes available. The duration of the time series available in the xed-sphere reference frame was restricted to the period in which box A remained within the FOV. For the experiments, this ranged between 2:2 21:6 s and the number of snapshots available ranged between m = x FOV =(2U o t PIV ) = 60 300 frames. Where possible, velocity measurements were extracted over a streamwise extent corresponding roughly to the DMD domain from these re-sequenced datasets. However, it was not always possible to ensure an exact overlap with the DMD domain of x2 [20:25R; 29:75R]. Indeed, for several low Froude number cases, the streamwise window available from PIV measurements was signicantly smaller (approx- imately 3:5R) than the DMD domain. 87 22 24 26 28 -4 -3 -2 -1 0 1 2 3 4 0.8 0.85 0.9 0.95 1 0.7 0.8 0.9 1 -0.1 0 0.1 0 10 20 30 -0.04 0 0.04 Figure 4.5: (a) Streamwise velocity contours normalized by freestream velocity U o . Three point sensors are shown in blue, green, and red. These sensors are displaced from a common central location, (x;y;z) = (24R;0:25R;0:25R), by a distance of 1:25R in the streamwise, lateral, and vertical directions, respectively (n.b., the green sensor is not in the same xz plane as the blue and red sensors). (b) Streamwise, vertical, and spanwise velocity signals are plotted against time for all the three sensors. The signals are color coded by sensor. 4.2.3 Problem Formulation Brie y, the regime identication problem can be stated as follows. Given a time series of input data from an unknown regime and a large database of DMD modes from known regimes, we seek a sparse representation of the training data in terms of the DMD modes in the database. In other words, a model is constructed for the input ow eld, approximating it as a linear combination of the dominant DMD modes identied from the large database. The parameter regime is then estimated by taking a projection-weighted average of the known Re andFr values for the selected DMD modes. We consider regime identication using two dierent types of input data: (a) planar (two-dimensional, two-component; 2D2C) velocity elds or (b) point measurements of all three-components of velocity from a limited set of sensors. The problem formulation is similar in both scenarios, though the DMD modes must be sub-sampled dierently depending on the input data prior to forming the database. 88 The 2D2C snapshots used for regime identication are obtained in the central xz (streamwise-vertical) plane of the sphere, i.e., corresponding to y = 0. If 2D2C snapshots from DNS are used, these include data over the same streamwise and vertical extents as the DMD modes. If 2D2C snapshots from experiments are used, the streamwise extent varies slightly from case to case but is comparable to that for the DMD modes. The point measurements used for regime identication are extracted from the DNS database. In this case, we assume the availability of 3 sensors, each capable of measuring all components of velocity. These 3 sensors are oset by a distance of 1:25R in streamwise, lateral, and vertical directions, respectively, from a common central point, which we term the sensor location. This setup is shown in gure 4.5. The time series of K velocity measurements (or inputs) obtained at intervals of t, fz 1 ; z 2 ; ; z K g, is rst stacked into a single column vector of the form: z = 2 6 6 6 6 6 6 6 4 z 1 z 2 . . . z K 3 7 7 7 7 7 7 7 5 : (4.6) Note that each individual element of z represents a complete 2D2C spatial snapshot or all the velocity readings from the point sensors. Since these input measurements are obtained at dierent spatial and temporal resolutions relative to the DMD dataset, the spatial DMD modes must be sub-sampled and the time evolution must be modied before the database of DMD modes is created. To ensure compatibility with the input data, a time sequence from 89 DMD mode j in regime i is created as follows: p j;i;h = 2 6 6 6 6 6 6 6 4 0 j;i S h v j;i t j;i S h v j;i . . . (K1)t j;i S h v j;i 3 7 7 7 7 7 7 7 5 (4.7) where t is the ratio of the sampling rates of training data and the DNS database used for DMD (i.e., t = t=t DNS ), j;i and v j;i are the complex Ritz values and vectors for DMD mode j from regime i, and S h is the h th sub-sampling matrix that is spatially compatible to the input data. For the 2D2C input data, we simply extract the central plane from the DMD spatial modes and so only a single sub-sampling matrix is used. For the point sensor data, we extract all compatible 3-point combinations from the DMD spatial modes, i.e., we consider all potential values for the central sensor location. In this case, the regime identication algorithm also provides an estimate for the sensor location. Note that factor of t in the exponential for the Ritz values ensures time evolution in the DMD modes at a rate consistent with the input data. Next, we assume that the input data can be modeled as a linear combination of the compatible DMD mode sequences: z< X j;i;h c j;i;h p j;i;h ! : (4.8) Here, c j;i;h 2 C are complex coecients that calibrate the amplitude and phase of each DMD mode sequence to match the training data, and<() represents the real component. Equation 4.8 can be rewritten as z X j;i;h (<(c j;i;h )<(p j;i;h )=(c j;i;h )=(p j;i;h )) (4.9) in which=() represents the imaginary component. 90 A large database is now constructed by pooling all the DMD mode sequences from each regime. To avoid computation with complex coecients, the real and imaginary parts of each compatible DMD mode sequence, p j;i;h , are included as separate basis functions d m in this database. Thus, equation (4.9) is now reformulated as a linear regression problem as follows: 2 6 6 6 6 6 6 6 4 z 1 z 2 . . . z K 3 7 7 7 7 7 7 7 5 | {z } input data 2 6 6 6 6 6 6 6 6 6 6 4 d 1 d 2 d M 3 7 7 7 7 7 7 7 7 7 7 5 | {z } basis functions 2 6 6 6 6 6 6 6 4 b 1 b 2 . . . b M 3 7 7 7 7 7 7 7 5 | {z } coefficients (4.10) or equivalently, z Db: (4.11) Here, D represents the database of M dierent DMD mode sequences (or basis functions) and b is the column vector of coecients for each of these basis functions. We seek a sparse solution in which only a small number of coecients are non-zero. In other words, we seek the dominant DMD mode sequences that best approximate the input. Recall that a total of n r = 24 dierent (Re;Fr) combinations were considered in DNS corresponding to 4 dierent Re values and 6 dierent Fr values. For the simulations at Re =f200; 300; 500g, 40 snapshots were used for DMD, yielding 39 DMD modes. For the simulations atRe = 1000, 100 snapshots were used for DMD, yield 99 DMD modes. However, only the rst 35 modes from each case were used to construct the database. Thus, there were a total of 35 24 = 840 dierent DMD modes available for the regression problem. For the 2D2C snapshots of velocity, each DMD mode yielded 2 dierent basis functions (real and imaginary part of p j;i;h ), for a total ofM = 1680 basis functions in the database D. For the point measurement input data, there were 252 dierent possible sensor locations (i.e., sub- sampling matrices S h ) spaced at a distance of 1:25R in thex, y, andz directions across the 91 3D spatial modes obtained from DMD. As a result, each DMD mode yielded 252 2 = 504 potential basis functions, which led to M = 840 504 = 423360 basis functions in the database. For regime identication with 2D2C snapshots from DNS, we obtain the streamwise and vertical components of velocity (u, w) over a grid of size n x n z = 77 65 in the streamwise and vertical directions. In general, we use 10 consecutive snapshots of the velocity eld obtained at a sampling rate of t = 2R=U o , such that t = (t=t DNS ) = 2 for Re =f200; 300; 500g and t = 5 forRe = 1000. Note that these data do not coincide with the DNS time series used for DMD. For regime identication from experimental data, we use a comparable spatial resolution though the sampling rate, t, is dierent for each case. However, for each case, we ensure that the DMD modes in the database are compiled at the sampling rate corresponding to the input data. The input data are split into training, validation and testing sets in the ratio of 60% : 20% : 20%, i.e., 6 snapshots for training, 2 for validation, and 2 for testing. A model is constructed for the training data, the so- called hyper-parameters (in this case, the number of modes retained in the regression) are ne-tuned based on the validation data, and the prediction condence is evaluated using the testing data. For regime identication using point measurements, 3 dierent sensor locations in the wake were considered. The rst location was directly in the center of the wake. The second and third locations were oset in the vertical and spanwise directions, respectively. Specics regarding sensor location are provided in Sec. 4.3.3. We assume each sensor measures all 3 velocity components at the sampling rate used for 2D2C snapshots, t = 2R=U o . A time series of 40 measurements is used as the input for regime identication. As before, these data are further split into training, validation, and testing sets in the ratio of 60% : 20% : 20%, i.e., 24 point measurements for training, 8 for validation and 8 for testing. 92 4.2.4 Identifying Dominant DMD Modes Using FROLS The vector of coecients b in Eqs. 4.10 and 4.11 can be computed using the standard least squares method. However, this will result in overtting due to the large number of modes used in the database. Therefore, we seek a sparse solution to this problem where the number of modes used to reconstruct the ow eld (M o ) is much less than the number used in the database (M). The reconstructed ow eld (^ z) can then be represented as: ^ z =b l 1 d l 1 +b l 2 d l 2 + +b l Mo d l Mo (4.12) in whichM o <<M andl m are the indices of the selected modes in the database. The above equation in matrix form is ^ z = 2 6 6 6 6 6 6 6 6 6 6 4 d 1 d 2 d M 3 7 7 7 7 7 7 7 7 7 7 5 | {z } basis functions 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 0 . . . b l 1 0 . . . b l Mo 0 . . . 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 | {z } sparse coeffs: (4.13) Notice that this problem is similar to identifying dominant resolvent modes from a large database in section 3.2.3 of chapter 3. We use the same sparse regression algorithm, FROLS, with the same modication. A brief recap of the steps of FROLS with modication is reproduced here for completeness. Recall that the real and imaginary parts of the DMD modes are included in the database as separate basis functions (see Eqs. 4.9 and 4.10). We implement a minor change to the algorithm in the ranking and selection step to select a 93 measurements a) DMD mode database weights b) a) measurements Figure 4.6: (a) Overview of the sparse regression technique. The FROLS algorithm identies and calibrates the dominant DMD modes from a large database that best represent the input measurements. The weight vector is sparse. (b) Reconstructed snapshots in comparison with the the input ow elds. complete DMD mode rather than just picking the real or imaginary part. This modication does not alter the underlying concept of the algorithm and is merely an adaptation for the problem at hand. This change is implemented as follows. After ranking each of the basis functions in the database, we select the basis function that has the highest ERR (Eq. 3.21) with the training data and add it to O. Then the rest of the basis functions in D are orthogonalized. Now, instead of ranking the modes in the database again, we pick the complex conjugate of the mode selected in the ranking step. For example, if the selected basis function is the real part of a particular DMD mode, we select its imaginary part and add it to O. The ranking step is then repeated after orthogonalizing the rest of the basis functions in D. This process ensures that both the real and imaginary parts of a particular DMD mode are included in the set O. Note that inclusion of both the real and imaginary components of the DMD modes ensures that arbitrary phase shifts can be captured in the reconstruction. This subtlety does not apply for DMD modes corresponding to the mean ow eld. Figure 4.6 summarizes this sparse regression algorithm. The FROLS algorith identies the dominant DMD modes from a large database to represent the input ow eld. The weight vector is now sparse where only a few elements are non-zero. Subgure (b) shows 94 the reconstructed snapshots in comparisons with the input training snapshots. 4.2.5 Estimating Regime and Quantifying Accuracy and Con- dence Recall that the input data is split into three blocks: a training set, a validation set, and a testing set. The training data is used to construct the model, the validation set is used to ne-tune the hyperparameters (i.e., the number of FROLS iterations, M o ), and the testing set is used to evaluate the condence in the predicted regime. To evaluate the optimum number of iterations, the error is tracked for the training and validation sets using Eq. 3.26. For the training set, we expect the error to decrease monotonically with the number of iterations. For the validation set, we expect the error to decrease initially but then increase once the model begins to overt to the training data. The iteration that yields the lowest error for the validation set is therefore deemed optimal and the calibration coecients are evaluated using Eq. 3.24 for this value of M o . From the identied and calibrated modes in the database, the (Re;Fr) regime of the input signal (together with location for the point measurements) is then estimated by computing a projection-weighted average Re w := P Mo i=1 W i Re l i P Mo i=1 W i ; Fr w := 10 P Mo i=1 W i log 10 Fr l i P Mo i=1 W i ! ; (4.14) in which Re l i and Fr l i are the Reynolds number and Froude number for the selected basis function l i . Note that the weighted average for Fr is evaluated in logarithmic space rather than linear space to re ect the available DNS data, which spanned Fr = 0:5 to Fr = 16. The scalar projection of the calibrated mode on to the training data is used as the weight: W i = (b l i d l i ) z kzk 2 (4.15) 95 For regime identication using point measurements, the sensor location is also estimated us- ing a projection weighted average of the sensor coordinates for the identied basis functions. To quantify prediction accuracy over the parameter range covered by the DNS and experiments (see gure 4.1), we transform theRe andFr space into a unit square using the following transformation: Re 0 = Re 200 1000 200 ; Fr 0 = log 10 Fr log 10 0:5 log 10 16 log 10 0:5 : (4.16) We then compute the distance between the prediction and the actual value and normalize it with the diagonal of the unit square. Thus, accuracy is dened as: Accuracy: = 1 2 6 4 Re 0 Re 0 w Fr 0 Fr 0 w 3 7 5 = p 2 (4.17) Since we use a weighted average, the worst possible prediction cannot exceed the diag- onal of the unit square. This choice of normalization is to ensure that the error lies between 0 and 1. An error of 0 implies perfect prediction and an error of 1 implies the worst possible prediction. Condence in the constructed model, and hence the estimated regime, is evaluated using the testing data. The condence is dened as the alignment between the predictor (^ z) and the input data (z) over the testing data: Condence: = z test ^ z test kz test kk^ z test k ; (4.18) where the subscript is used to indicate that this metric is computed only over the testing data. Further, we only compute condence using the vertical velocity component (w) for recon- structions based on 2D2C velocity snapshots, and using both vertical and lateral components (v andw) for reconstructions based on point measurements. The streamwise component (u) 96 1 8 16 24 24 16 8 1 0.99 0.992 0.994 0.996 0.998 1 1 8 16 24 24 16 8 1 0 0.2 0.4 0.6 0.8 1 Figure 4.7: Component-wise alignment between mean ow elds is evaluated between each of the (Re;Fr) cases in gure 4.1. There are a total of 24 cases, yielding a 24 24 grid. (a) Alignment between streamwise mean velocity elds,ju i u j j=(ku i kku j k), and (b) Alignment between vertical mean velocity elds,jw i w j j=(kw i kkw j k). Here, the () notation implies an element-wise product and sum. is not included in the calculation because there is a very high degree of alignment between the mean streamwise elds (u) of all the regimes, which dominates the condence metric. Figure 4.7 shows the alignment in mean ow elds between each of the dierent (Re;Fr) cases for which we have DNS data. There are a total of 24 cases and so we have a 2424 grid showing alignment between cases. As expected, the plot is symmetric and all the diagonal cases (i = j) show perfect alignment. Figure 4.7(a) shows that there is a very high degree of alignment between the streamwise mean velocity elds . In contrast, gure 4.7(b) shows that, for the mean vertical velocity, there is very little alignment between cases. Therefore, if the streamwise velocity is included in the computation of the condence metric in Eq. 4.18, then the reconstructions show high condence levels irrespective of whether the model has the correct modes. Thus, usingw (andv) to evaluate condence is more appropriate for the reconstructions. Note that we also pursue reconstructions in which the mean ow eld is removed entirely from the input data. This is discussed in greater detail below. We recognize that prediction accuracy and condence can be evaluated using metrics other than the ones dened in Eq. 4.17-4.18. Nevertheless, we believe that these metrics provide useful insight into performance. 97 4.3 Results and Discussion In this section we discuss the results for various types of input data. Section 4.3.1 deals with cases in which 2D2C Velocity elds from numerical simulations are used as the input. We rst evaluate the accuracy of the method described in the previous section using input data from a known regime and quantify the accuracy and condence over the available (Re;Fr) space. We then attempt regime identication with a sub-sampled database. Specically, we include DMD modes from half the available cases in the database and again evaluate the regime predictions to test if the method is able to interpolate between regimes. In Section 4.3.2, we use 2D2C PIV snapshots from experiments and the same database as before, i.e., DMD modes obtained via DNS. The experimental data is inherently noisy. In addition, the spatial windows over which velocity elds are available do not always overlap with the spatial windows used for the DMD modes in the database. The accuracy and condence is again quantied for the experimental data. We see the results obtained in this section as an indicator of regime identication performance for real-world measurements. Finally, in Section 4.3.3 we make use of three point measurements of velocity as the input and evaluate predictions using both the full and the sub-sampled DMD databases. In contrast to the database used for reconstruction of 2D2C velocity elds, the database for this scenario contains many possible triplets of velocity data from each DMD mode, corresponding to dierent sensor locations, as separate basis functions. Therefore, in addition to identifying the regime, we also attempt to identify the sensor location. Regime estimation and ow eld reconstruction is rst carried out for (Re;Fr) = (500; 4). We focus on this particular case because it falls at the intersection between three dierent regimes (see gure 4.1) which makes reconstruction challenging. Following this, regime estimation, accuracy, and condence are evaluated for all 24 (Re;Fr) combinations in gure 4.1. Regime estimation for each case is performed using the full ow eld as well as the mean-subtracted ow eld. Note that mean subtraction is performed only for the input 98 data used for regime identication, i.e., the 2D2C snapshots or point measurements that constitute z in Eq. (4.11). The DMD modes used to create the database of candidate basis functions, D in Eq. (4.11), are computed using the full ow eld from DNS. A summary of results for the dierent scenarios considered is provided in Table 4.1. Full ow eld Mean-subtracted Input data Source Database (%) (%) (Fr) (%) (%) (%) (Fr) (%) 2D-2C DNS Full 92.8 87.3 97 83.8 68.8 91.2 Snapshots (10) Sub-sampled 88.1 84.4 93 77.4 50.5 87.1 2D-2C DNS Full - - - - - - Snapshots (30) Sub-sampled 88.2 80.4 92 76.1 35 82.8 2D-2C EXP Full 71.8 66.7 78.6 68.3 38.1 79 Snapshots (10) Sub-sampled - - - - - - 2D-2C EXP Full 73.2 73.6 79.7 70.1 34.5 76.5 Snapshots (30) Sub-sampled - - - - - - 2D-1C EXP Full 74.4 70.3 86.2 69.2 38.6 81.2 Snapshots (10) Sub-sampled - - - - - - 2D-1C EXP Full 73.9 79.4 83.8 67.3 36.5 76.7 Snapshots (30) Sub-sampled - - - - - - Point Meas. DNS Full 95.5 83.4 95.7 84.2 49.3 91.1 Location 1 Sub-sampled 75.9 71.8 80 71.7 43.5 80.7 Point Meas. DNS Full 67.9 80.6 80.6 76.8 63.5 79.9 Location 2 Sub-sampled 66.6 74.4 78.5 70.2 51.7 76.5 Point Meas. DNS Full 70.8 78.3 77.1 76.8 66.6 88.4 Location 3 Sub-sampled 71.8 68.0 78.6 68.7 56.2 80.6 Table 4.1: The accuracy (, Eq. 4.17) and condence (, Eq. 4.18) metrics, averaged over all 24 (Re;Fr) cases, using the full ow eld and the mean-subtracted ow eld. Results are shown for dierent input data and database sizes. The numbers in parentheses (10 or 30) indicate the number of 2D snapshots used in the input dataset. For the cases in which 10 snapshots are used, the accuracy and condence metrics are averaged over 3 dierent realizations. 4.3.1 2D2C PIV Snapshots from DNS We rst begin with input data from DNS for one of the 24 cases, (Re;Fr) = (500; 4), and quantify the accuracy of prediction and also assess the quality of the reconstructed ow eld. We use this specic parameter combination since it lies at the intersection of several 99 22 24 26 28 -4 -2 0 2 4 22 24 26 28 -4 -2 0 2 4 0.8 0.85 0.9 0.95 1 22 24 26 28 -4 -2 0 2 4 22 24 26 28 -4 -2 0 2 4 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 1 50 100 10 -3 10 -2 20 40 60 80 200 300 500 1000 0.5 1 2 4 8 16 Figure 4.8: Comparisons between streamwise (a, b) and vertical velocities (c, d) over the testing data for (Re;Fr) = (500; 4). (e) Evolution of error with number of iterations for training and validation data sets. The minimum point for the validation error (i.e., the optimal iteration) is highlighted with a red circle. (f) The evolution of the predictions Re w and Fr w is shown until the optimal iteration. 100 dierent dynamic regimes (see gure 4.1). For these reconstructions, we use a full database, i.e., with DMD modes from all 24 cases. The FROLS algorithm is then used to identify the DMD modes from this database that best represent the training data, and to calibrate their amplitude and phase. After calibration, the weighted DMD modes are used to reconstruct the ow eld. Importantly, note that the 10 snapshots used as the input data (6 for training, 2 for testing, 2 for validation) are dierent to the snapshots used for DMD. As expected, the method yields perfect predictions for Re andFr if the snapshots used as input data are from the set used for DMD; these predictions are not shown here for brevity. Figure 4.8 shows comparisons between the DNS ow eld and the reconstructed ow eld for one of the test data snapshots. Panels (a) and (b) show comparisons for the stream- wise velocity, while panels (c) and (d) show comparisons for the vertical velocity. The reconstructions are labelled with a ^ (). The reconstructed velocity eld shows very good qual- itative agreement with the velocity elds from DNS. More quantitative results are shown in panel (e). The reconstruction error for the training data (Eq. 3.26) decreases monotonically whereas the error for the validation data reaches a minimum at the 86 th iteration, which suggests that the regression algorithm begins to overt to the training data beyond this point. Hence, we reconstruct the velocity eld and identify the parameter regime using only the modes selected by the algorithm up to iteration 86. The corresponding regime predic- tion using the projection weighted average (Eq. 4.14) is shown in panel (f). Note that the algorithm yields a very close prediction, Re w = 497 andFr w = 3:9. The condence in these predictions, evaluated using Eq. 4.18, is 99%. Recall that the condence is a measure of how similar the vertical velocity uctuations in the reconstructed ow eld are to the testing data. In this case, the vertical velocity elds for the testing data and the reconstruction look very similar. As a result, we have a high level of condence in the trained model and hence the predicted regime. Next, we generate predictions for each of the 24 dierent (Re;Fr) cases in the database. Figure 4.9 shows the accuracy and the condence of prediction computed using Eqs. 4.17 101 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 Figure 4.9: Plots showing prediction accuracy and condence averaged over three dierent sets of input data for all (Re;Fr) cases in the database. Regime predictions made using full ow elds as input data are shown in (a,c) and predictions generated using mean-subtracted ow elds are shown in (b,d). Large boxes represent prediction accuracy (Eq. 4.17) while smaller inset boxes represent condence (Eq. 4.18). Panels (a,b) show the accuracy computed using both Re and Fr, while panels (c,d) show the accuracy computed using Fr alone. and 4.18. The color of the larger square box represents accuracy while the color of the small inset box is an indicator for condence. Here, we also compare regime predictions made using the full velocity eld with predictions made using the mean-subtracted velocity eld, i.e., using only the uctuating velocity eld. The accuracy and condence metrics shown in gure 4.9 are averaged over three dierent sets of input data. Overall, the prediction accuracy is higher when the mean ow is included in the training data, suggesting that 102 the mean ow eld contains useful information when it comes to distinguishing between (Re;Fr) values. Specically, the average accuracy across all 24 cases is approximately 93% for regime predictions made using the full velocity eld, compared to 84% for the predictions made using the uctuating velocity eld alone (see Table 4.1). A comparison of the results presented in gure 4.9(a) and gure 4.9(b) shows that, when the mean-subtracted velocity eld is used, prediction accuracy deteriorates signicantly for the cases that are classied as symmetric or asymmetric non-oscillation (red and yellow regions in gure 4.1). This is expected since the dominant DMD modes for these regimes have St = 0 (gure 4.2). The FROLS algorithm is unlikely to select these St = 0 modes without mean ow information in the input data. There is a qualitative correlation between prediction accuracy and the so-called condence metric used in this study. The condence metric is higher for cases in which the regime prediction is accurate, suggesting that the trained model is a useful future predictor of velocity, i.e., it yields accurate reconstructions for the testing data. Panels (a) and (b) in gure 4.9 indicate that regime identication is less successful for cases with high Reynolds and Froude numbers. Recall that the accuracy metric plotted in panels (a) and (b) represents the normalized distance between the actual (Re;Fr) and the predicted values (Re w ;Fr w ) over the parameter space considered in this study (Eq. 4.17). Panels (c) and (d) in gure 4.9 show accuracy computed using only the normalized distance between actual and predicted Froude number values, Fr and Fr w . A comparison between panels (a,b) and (c,d) shows that the framework developed here generates much better predictions for Froude number compared to Reynolds number. For example, the average prediction accuracy for the high Reynolds and Froude number cases (Re = 1000 and Fr 2) is 67% when considering the normalized distance in (Re;Fr) space compared to 88% when only considering Fr. This indicates that the dominant source of error lies in the Re prediction. One potential explanation for this is the fact that the high Re and Fr cases fall into a similar 3D dynamic regime (see gure 4.1), and so the 2D2C snapshots used for regime identication and ow reconstruction may not have sucient information to 103 distinguish between dierent Reynolds numbers. Another possible explanation is that the input sequence comprising 10 snapshots may not be sucient to capture the wider range of timescales expected for the higher Re cases. This is explored further below. As a further test, we evaluate prediction accuracy using a more limited database of DMD modes. In other words, we test whether the method is able to generate reasonable (Re;Fr) predictions even if DMD modes from that specic parameter combination are not available in advance. Specically, we make use of a database with DMD modes from only half the available (Re;Fr) cases; every other case in gure 4.1 is eliminated. We assess prediction accuracy for cases that are not in the database as well as the ones that are included. Figure 4.10 again shows predictions for input data from the case with (Re;Fr) = (500; 4). However, these predictions are made using the more limited database, which did not contain DMD modes for this specic parameter combination. Panels (a)-(d) show that, even though the reconstructions are composed purely of DMD modes corresponding to other regimes, the reconstructed velocity elds qualitatively match the velocity elds from DNS for both the streamwise and vertical components. Interestingly, panel (f) shows that the FROLS algorithm picks a DMD mode corresponding to (Re;Fr) = (1000; 4) in the rst iteration for this case. However, as the algorithm selects more modes the prediction evolves and eventually converges to a value of (Re w ;Fr w ) = (664; 3:9), which is reasonably close to (Re;Fr) = (500; 4). The nal condence level predicted using Eq. 4.18 is approximately 95%. Figure 4.11 shows the accuracy and condence metrics for all 24 (Re;Fr) cases based on predictions made using the sub-sampled DMD database. Unsurprisingly, prediction accuracy and condence are higher for input data from (Re;Fr) combinations included in the database (marked with a). For predictions made using the full ow eld (gure 4.11(a)), the average prediction accuracy is 95% for (Re;Fr) cases included in the database and 81% for cases excluded from the database. As before, predictions made using mean-subtracted input data are less accurate than regime predictions made using the full velocity eld (average accuracy 104 22 24 26 28 -4 -2 0 2 4 22 24 26 28 -4 -2 0 2 4 0.8 0.85 0.9 0.95 1 22 24 26 28 -4 -2 0 2 4 22 24 26 28 -4 -2 0 2 4 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 1 50 100 10 -2 2 10 -2 5 10 15 200 300 500 1000 0.5 1 2 4 8 16 Figure 4.10: Comparison of reconstructions between streamwise (a, b) and vertical velocities (c, d) over the testing data for (Re;Fr) = (500; 4) using a sub-sampled DMD database. (e) Evolution of error with number of iterations for training and validation data sets. The minimum point for the validation error (i.e., the optimal iteration) is highlighted with a red circle. (f) The evolution in Re w and Fr w until the optimum iteration. The nal prediction is (Re w ;Fr w ) = (664; 3:9) with a condence of 94%. 105 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 Figure 4.11: Plots showing prediction accuracy and condence averaged over three dierent sets of input data. These predictions are made using a database only consisting of DMD modes from the cases marked with a cross (). Regime predictions made using full ow elds as input data are shown in (a) and predictions generated using mean-subtracted ow elds are shown in (b). Large boxes represent prediction accuracy (Eq. 4.17) while smaller inset boxes represent condence (Eq. 4.18). is roughly 88% for results shown in panel (a) and 77% for results shown in panel (b)). Importantly, with the exception of a few outliers, prediction accuracy and condence are reasonably high even for the (Re;Fr) cases not included in the database. This conrms that the method developed here is able to reconstruct ow elds from regimes not seen a priori, as long as the database contains basis functions that can reasonably approximate the input data. Similar to the results shown in gure 4.9, the prediction accuracy is signicantly higher for Fr alone (93% on average for the full velocity eld) compared to the accuracy for the combined (Re;Fr) metric (88% on average), i.e., the Re predictions are less accurate than the Fr predictions. For brevity, these predictions are not shown here explicitly. Further, prediction accuracy is again lower for the cases with Re = 1000. As noted earlier, one possible explanation for this is that the short input sequence comprising 10 snapshots may not be sucient to capture the more broadband temporal signature associated with wakes at higher Re. To evaluate whether a longer input sequence leads to an improvement in performance, we also attempted regime identication with 30 snapshots from DNS and the 106 sub-sampled database of DMD modes. However, as shown in Table 4.1, this did not yield signicant changes in prediction accuracy or condence. Table 4.1 provides a summary of prediction accuracy and condence for the dierent scenarios tested. 4.3.2 2D2C PIV Snapshots from Experiments Next, we evaluate whether the method developed here is able to identify parameter regimes using 2D2C PIV snapshots from experiments as input data. As before, we make use of 10 snapshots split into training, testing, and validation data. Compared to regime identication using 2D2C snapshots from DNS, there are three major challenges associated with the use of input data from experiments. First, the experimental data tend to be noisier. Second, though the experiments were carried out over the same approximate (Re;Fr) range as the DNS, the Re and Fr values do not match exactly with the DNS cases used to generate the DMD mode database. Third, the spatial extent of the PIV windows in experiment does not coincide with the streamwise extent of the DNS data used for DMD. To remove small-scale noise from the experimental data, we pursue snapshot POD for the experimental data [see e.g., 92] and reconstruct the velocity eld by retaining the 7 leading POD modes. In all cases, these leading modes capture more than 99% of the total kinetic energy. We do not take any explicit actions to address the second and third challenges and so the results presented below are indicative of real-world reconstruction performance. Figure 4.12 shows reconstruction and regime identication performance using input data from experiments conducted at (Re;Fr) = (533; 3:7), similar to the (Re;Fr) = (500; 4) case considered in gures 4.8 and 4.10. The streamwise extent of the wake captured in the experimental data, x=R2 [26:5; 36], overlaps partially with the streamwise extent for the DMD modes, x=R2 [20:25; 29:75]. Despite the dierences in (Re;Fr) values and spatial extent, the method developed here yields reasonable reconstructions for the experimental velocity elds using DMD modes obtained from DNS (see panels (a)-(d) in gure 4.12). The condence metric, which is a measure of reconstruction accuracy for the test data, is 107 22 24 26 28 -4 -2 0 2 4 28 30 32 34 36 -4 -2 0 2 4 0.75 0.8 0.85 0.9 0.95 1 22 24 26 28 -4 -2 0 2 4 28 30 32 34 36 -4 -2 0 2 4 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 0.04 0.05 1 50 100 0.01 0.02 0.03 0.04 5 10 15 20 200 300 500 1000 0.5 1 2 4 8 16 Figure 4.12: Comparisons between measurements and reconstructions for the streamwise (a, b) and vertical velocity (c, d) over the testing data for the experimental dataset at (Re;Fr) = (533; 3:7). Note the diering spatial extent of the wakes for the experimental and DNS data. (e) Evolution in error for the training and validation data sets. The minimum point for the validation error (i.e., the optimal iteration) is highlighted with a red circle. (f) The evolution of Re w and Fr w until the optimal iteration. The nal prediction is (Re w ;Fr w ) = (349; 7:5) with a condence of 83% per Eq. 4.18. 108 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 Figure 4.13: Prediction accuracy (large boxes, Eq. 4.17) and condence (smaller inset boxes, Eq. 4.18) averaged over three snapshot sequences from experiments, using the full DMD mode database. As before, panel (a) shows results for the full ow eld while panel (b) shows results obtained using the mean-subtracted ow eld. Experimental data for the Fr = 16 cases are unavailable. However, the corresponding DMD modes from DNS have been included in the database. 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 Figure 4.14: Prediction accuracy (large boxes, Eq. 4.17) and condence(smaller inset boxes, Eq. 4.18) using experimental measurements of only the vertical componentw, averaged over 3 dierent sets of input data. As before, (a) shows results obtained using the full ow elds while (b) shows results obtained using the mean-subtracted ow elds. estimated to be approximately 83%. Figure 4.12(f) shows that the regime predictions con- verge to (Re w ;Fr w ) (349; 7:5). At rst glance, there is a signicant discrepancy between the estimated Froude number and the actual value, Fr = 3:7. However, it is important to 109 keep in mind that the DMD mode database spans nearly 2 orders of magnitude in Fr. In this context, the algorithm does a reasonable job of localizing parameter values. Moreover, per the classication in gure 4.1 the estimated ReFr combination falls along the same interface between physical regimes as the actual parameter combination. We next compute the accuracy and condence for all the available experimental datasets. Results averaged over three runs are shown in Figure 4.13. Note that the PIV FOV for the experimental training data and, hence, the overlap with the DMD modes varies for each case. The wake window is smallest for (Re;Fr) (1000; 0:5) with a streamwise extent of 3:125R and larger for higher Fr cases. For cases in which experimental data are available over a larger spatial window, the streamwise extent is capped at 9:5R, which is the window used for DMD. Reconstruction and regime identication is carried out with the smaller of the two windows (i.e., experimental and DMD) both in streamwise extent as well as in the vertical direction. For some cases there is signicant overlap in the streamwise direction (e.g., approximately 8R for (Re;Fr) = (500; 2)) and for others there is little or no overlap (e.g., for (Re;Fr) = (200; 0:5)). In general, the regime identication accuracy and condence are lower with the experimental input data compared to training data from DNS. As shown in Table 4.1, regime identication accuracy is roughly 72% on average with the full ow eld from experiments compared to 93% with the DNS snapshots. The lower accuracy and condence could be attributed to the noise inherent in the experimental data and dierences in streamwise extent between the experimental data and the DMD modes. Nevertheless, these results serve as proof-of-concept for regime identication with real-world experimental data. As before, prediction accuracy is lowest for the high Re cases, and a comparison of the results in gure 4.13(a) and (b) shows that predictions improve when the full velocity eld from experiments is used as the input rather than just the uctuations (72% vs. 68% per Table 4.1). It is worth noting here that reconstruction carried out using just the vertical component of velocity (w) yields similar results to reconstruction based on both the components of 110 velocity. This is illustrated by a comparison between gure 4.13, which shows results obtained using both components of velocity, and gure 4.14, which shows results obtained using only the vertical component. As shown in Table 4.1, overall regime identication accuracy improves by roughly 3% with the use of w alone, and there is a nearly 8% improvement in predictions for Fr. As was the case with the snapshots from DNS, prediction accuracy is higher for Fr alone compared to the composite metric dened in Eq. 4.17. Similarly, per Table 4.1, using a sequence of 30 snapshots from the experiments as the input does not lead to any noticeable improvement in prediction accuracy. 4.3.3 Point Measurements As a nal test of the framework developed here, we now attempt regime identication using point measurements (extracted from DNS data). We present results corresponding to three nominal sensor locations: (x;y;z) 1 = (24;0:25;0:25), (x;y;z) 2 = (25:875;0:25; 1:25), and (x;y;z) 3 = (25:875; 1:25;0:25). Recall that these nominal sensor locations represent the central point around which 3 separate point sensors are arranged; the point sensors are respectively placed at a distance of 1:25R in the x, y, and z directions from the nominal sensor location. For example, the coordinates of the three sensors corresponding to sensor location 1 are (25:25;0:25;0:25), (24;0:25; 1), and (24; 1;0:25). Sensor location 1 corresponds to measurements made directly downstream of the object. Sensor locations 2 and 3 are oset in the vertical and spanwise directions, respectively, from the center of the wake. We assume that each point sensor is capable of measuring all three components of velocity. The input data consists of a time series of 40 measurements from each sensor, with 24 measurements dedicated for training, 8 for testing, and 8 for validation. Though we do not pursue an extensive survey of sensor locations, the results obtained in this section provide useful insight into sensor placement for regime identication. We rst start with regime identication based on measurements directly in the wake of the object. Figure 4.15 shows comparisons between the point measurements extracted 111 0.7 0.8 0.9 1 -0.1 0 0.1 0 10 20 -0.04 0 0.04 0.7 0.8 0.9 1 -0.1 0 0.1 0 5 -0.04 0 0.04 1 50 100 10 -3 10 -2 5 10 -2 5 10 15 20 200 300 500 1000 0.5 1 2 4 8 16 Figure 4.15: Comparisons between the predictor (^ z) and the DNS (z) over the training (a) and testing (b) datasets. These reconstructions correspond to sensor location 1. The evolution of reconstruction error with the number of iterations for the training and validation sets is shown in (c). The corresponding variation of Re w and Fr w is shown in (d). The input data corresponds to (Re;Fr) = (500; 4) and the regime prediction corresponds to (Re w ;Fr w ) = (500; 4) with a condence of 97%. from DNS (z) and the reconstructions (^ z) over the training (a) and testing (b) data for the representative case corresponding to (Re;Fr) = (500; 4). The full DMD database is used for reconstruction. Panel (c) shows the evolution in reconstruction error with the number of iterations. There is a clear minimum in reconstruction error for the validation data after 20 112 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 Figure 4.16: Prediction accuracy (large boxes, Eq. 4.17) and condence (smaller inset boxes, Eq. 4.18) with point measurements as input data using the full DMD mode database. Sub- plots (a), (b), and (c) correspond to sensor locations 1, 2, and 3 respectively. 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 200 300 500 1000 0.5 1 2 4 8 16 0 0.25 0.5 0.75 1 Figure 4.17: Error in predicted sensor location normalized by R. Subplots (a), (b), and (c) correspond to sensor locations 1, 2, and 3 respectively. Note that the error is capped at 1R and red shading indicates lower error. iterations. Panel (d) shows that the predictions converge toRe w = 500 andFr w = 4 after 20 iterations, yielding a perfect match with the input data. Recall that the method developed here is also able to estimate the sensor location based on the identied basis functions. For this case, the identied sensor location matches perfectly with the actual sensor location (see gure 4.17(a)). The close match between the input and reconstructed ow elds for the testing data also yields a high level of condence (97%) per the metric dened in Eq. 4.18. 113 Prediction accuracy and condence for all 24 (Re;Fr) cases are shown in gure 4.16 for the three dierent sensor locations. As shown in Figure 4.16(a), for sensor location 1 the method yields high-accuracy and high-condence predictions for nearly all cases. Indeed, the average prediction accuracy is over 95% for this scenario (see Table 4.1). However, the method does not yield accurate predictions for the case with (Re;Fr) = (1000; 16). Figure 4.17 shows the normalized error in predicted sensor location. The normalized error is dened as max(=R; 1), in which is the distance between the actual and estimated sensor location. Since the wake develops slowly in the streamwise direction, the distance is computed using the vertical and spanwise coordinates only; the streamwise coordinate is not included. Consistent with the regime predictions, gure 4.17(a) shows that the method is able to identify the correct sensor location in nearly all the cases based on measurements at sensor location 1. Recall that this location is directly downstream of the body, roughly centered in the wake. Figures 4.16(b,c) show that prediction accuracy and condence deteriorate for input data from sensor locations 2 and 3. The average prediction accuracy is roughly 68% for sensor location 2 and 71% for sensor location 3 (c.f., 95% for location 1; see Table 4.1). Figures 4.17(b,c) show a similar deterioration in sensor location estimates. This deterioration of performance is perhaps not surprising given that sensor locations 2 and 3 are oset from the center of the wake. Interestingly, the use of mean-subtracted ow elds alleviates some of the performance degradation in regime identication from measurements at sensor location 2. Table 4.1 shows that the average prediction accuracy improves from 68% to 77% if only the velocity uctuations from sensor location 2 are used for regime identication. The improvement of performance with the use of mean-subtracted ow elds is also evident for sensor location 3; average prediction accuracy improves from roughly 71% for the full ow eld to 77% for uctuations alone. In general, the results presented in Table 4.1 show a smaller variation in prediction accuracy across the dierent sensor locations when mean- subtracted ow elds are used as input data. This suggests that the mean velocity eld is 114 useful for regime identication from point measurements, but primarily for measurements made close to the wake centerline. For measurements made away from the center, the mean component may have a confounding in uence. For completeness, we note that we also pursued regime identication from point mea- surements using the sub-sampled DMD mode database. Similar to the trends observed for 2D2C snapshots, the use of the smaller database led to a reduction in overall prediction accuracy. Further, prediction accuracy was generally better for just Fr compared to the composite ReFr metric dened in Eq. 4.17. These results are summarized in Table 4.1. 4.4 Chapter Conclusions In this chapter, we developed a data-driven framework for identifying the dynamic regime (i.e., ReFr values) for stratied wakes from limited velocity snapshots or point measure- ments. This approach makes use of a library of DMD modes available from prior DNS. The FROLS algorithm is used to identify DMD modes that best represent the input data, and ReFr are then estimated via a projection-weighted average of parameter values for the selected DMD modes. We demonstrated the utility of this framework using input data from both numerical simulations and laboratory experiments. Regime identication accuracy was very high (> 90%; see Section 4.3.1) when velocity snapshots from numerical simulations were used as input data. Moreover, prediction accuracy deteriorated only slightly when the library (or database) of DMD modes was sub-sampled. This conrms that the method is ca- pable of interpolating between parameter values when the specicReFr combination is not included in the library. The average regime identication accuracy from velocity snapshots obtained in laboratory experiments was lower (> 70%; see Section 4.3.2). This deteriora- tion in performance was expected given that the laboratory measurements were inherently noisier and did not have the same spatial extent as the DMD modes obtained from numer- ical simulations. Given these challenges, the reasonable prediction accuracy suggests that this technique is robust and could be used for real-world regime identication. The results 115 presented in Section 4.3.3 show that this framework is capable of regime identication from point measurements. When using point measurements, the method also yields reasonable estimates for sensor location. We recognize that the parameter combinations considered here correspond to low Reynolds number ows that are relatively narrow-banded. Regime identication is expected to become more challenging for higherRe ows that exhibit a wider range of spatiotemporal scales. In other words, the limited datasets typically obtained from laboratory and eld measurements may not contain enough information to allow for meaningful delineation between regimes at higher Re. Even in the present study, a deterioration in performance is evident at Re = 1000 for inputs from both numerical simulations and experimental measurements. Yet another challenge for high Re applications is the limited availability of high-delity simulation and experimental data. This would restrict the library of candidate basis functions that can be used for regression. However, these challenges associated with data scarcity can be overcome with the development of better measurement systems, the use of emerging techniques for optimal sensor placement and data fusion [58, 105], and continued advances in our ability to simulate higher Reynolds number ows. In many ways, this study is complementary to the neural network based classiers [e.g., 27, 5] and expert-dened decision trees [e.g., 72, 73] developed in recent work on wake ows. However, the approach developed here has some advantages. First, the library-based framework developed here is exible and can be used with dierent forms of input data (e.g., snapshots or point measurements with varying spatiotemporal delity). Feed-forward neural networks must be developed and trained independently for dierent input data and similarly decision trees must be modied to account for diering data types. Second, in addition to classifying the data into dierent regimes, the framework developed here can also reconstruct the velocity eld. The condence in the constructed model can therefore be evaluated based on how well the testing data are predicted. 116 Note that the framework developed here is not limited to the use of DMD modes. The FROLS algorithm simply identies the best basis functions that represent the data. These basis functions can be Fourier modes, DMD modes, POD modes, or indeed, any set of candidate basis functions that is suciently rich to characterize the input data. Previous work shows that, in addition to DMD, Spectral POD (SPOD; [97]) can also be used to characterize dierent wakes [71, 73]. Additional work is needed to see if the use of SPOD modes improves regime identication and reconstruction accuracy for stratied wakes. Of course, any such evaluation must also account for the increased data requirements for SPOD, which relies on ensemble averaging. Future work also includes a systematic evaluation of reconstruction accuracy as a function of spatial and/or temporal resolution of the input data as well as the duration of the available time series. Nevertheless, this manuscript shows that the library-based sparse regression formulation developed here holds promise for the development of automated data-driven uid pattern classiers. 117 Chapter 5 Conclusions The objective of this thesis is to develop techniques that can be used to reconstruct the velocity eld or estimate dynamic parameters from limited measurements for shear ows. These estimates may be in time or space, or just the regime of the ow. In chapter 2, we used models grounded in Rapid Distortion Theory and Taylor's hypothesis, both of which are derived from the governing Navier-Stokes equations with modeling assumptions. The results from chapter 2 show that both RDT and Taylor's hypothesis are useful dynamic models to make short time predictions for wall-bounded turbulent ows from non-time- resolved PIV snapshots. The only required input for these reconstructions is the mean velocity prole. Both these models reconstruct the ow eld with reasonable accuracy and are computationally inexpensive. When using the TH model with characteristics-based evolution of ow eld, the eects of articial viscosity are completely eliminated. This results in the preservation of the small scales or the high frequency components of the ow eld, which are otherwise smoothed out when integrating the discretized form of the equations. This results in higher accuracy of reconstructions and the reconstructed statistics and spectra more closely resemble the DNS ground `truth'. It is worth noting here that the premultiplied spectra for frequencies that are above the Nyquist limit of the acquisition rate of the PIV- like data are also resolved. However, because we are inferring temporal reconstructions from spatial information, the frequencies that can be resolved is dependent on the PIV resolution. 118 Next, in chapter 3, the velocity eld measurements were projected onto resolvent modes, which serve as the basis functions. The resolvent modes are obtained via a gain-based de- composition of the governing equations for turbulent channel ow at Re = 1000. A large database of resolvent modes computed at dierent wavenumbers and wavespeeds are com- piled into a large database. We now employ a two step process to identify and calibrate these modes. First, the dominant modes are identied using a sparse regression algorithm, FROLS. The coecients of these basis functions are optimized to minimize a composite cost function. This cost function is composed of cost from the input training velocity eld snapshots and another cost from the statistical constraints. This way, the reconstructed ow eld approximates the input training snapshots as well as the statistics. To estimate the ow eld between the two consecutive measurements, the identied modes are advected forward in time with their wavespeeds. Because the resolvent modes are computed directly from the governing equations, each mode contains information about the wall-normal variation of all the three velocity components and pressure. Once the modes have been identied and calibrated via limited eld measurements, it should be possible to estimate the unmeasured or unobserved variables such as the third component of velocity and pressure. Because the resolvent modes are three dimensional travelling wave ow structures, the training process identies a three-dimensional mode from two-dimensional data. Therefore volumetric recon- structions should also be possible. However, the data requirement to enable such estimation is an open question. In this study, we methodically increase and vary the type of data on which we train the model. We start with just 2 temporal snapshots in time. We then use simultaneous measurements from 2 orthogonal planes. These include streamwise - wall- normal, and streamwise - spanwise planes. We also test the reconstructions with varying number of training snapshots, molecular viscosity, and eddy viscosity models. The results show that the reconstruction error from unobserved variables decreases signicantly when the statistical constraints are imposed. Although the estimated unobserved variables do not match qualitatively with the DNS data, the magnitudes are within range. One reason for 119 this shortcoming could be due to the use of a limited set of resolvent modes. This includes using only rank 1 modes, using a limited number of wavenumber and wavespeed combina- tions. However, using an extensive set of resolvent modes increases the computational cost incurred and is beyond the scope of the present study. Finally, in chapter 4 we introduced a purely data-driven method solve the inverse prob- lem of nding the regime from the ow structures. For this, we use DMD modes from prior DNS to reconstruct the input velocity eld measurements and in-turn estimate the regime; both Re and Fr. First, we compile a large database by poling the DMD modes obtained from prior DNS. These simulations are carried out at 24 dierent locations inReFr space. A sparse regression algorithm, FROLS, is used to identify the dominant modes that best represent the input training data. Based on how well the reconstructed ow eld matches the set of input testing snapshots, we evaluate the condence in these predictions. Proof- of-concept for this technique is tested with a time series of 2D2C velocity eld snapshots or a time series of three point measurements. This is in contrast to previous numerical and experimental studies that used measurements in two orthogonal planes or had access to full 3D ow eld to make the classication. Importantly, this study directly estimates the parameter values instead of a class-based categorization. 5.1 Implications and Future Work The dynamical models developed in this thesis can be fused with real-world measurements to reconstruct the uid ows. Specically, we focused on developing dierent physics-based dynamical model to estimate ow elds. Each of these models has its own advantages and disadvantages. For example, in chapter 2 we develop a reduced-order-model (ROM) derived directly from the Navier-Stokes equations that can estimate the measured velocity elds at higher temporal resolution. Both RDT and TH are very inexpensive computationally. Be- cause of the simplicity of TH model, the method of charateristics can be used to estimate the solution in time thus retaining the high frequency uctuations. However, it cannot be used 120 to estimate these variables at locations other than the sensor locations. For this, we project the velocity eld on to resolvent modes. Once the dominant resolvent modes have been iden- tied and calibrated, they can be used to estimate unmeasured or unobserved variables and also estimate velocity elds at other spanwise locations where meaurements are not made. This interesting feature could also be exploited for other estimation problems. For example, in of compressible ows, the three velocity components, density, and temperature may be chosen as the state variables for the resolvent modes [8]. Therefore for these ows, density and temperature could perhaps be estimated from velocity measurements alone. Similarly, for incompressible ows when working with free or forced convection problems, including the energy equation in the resolvent mode computation, and temperature in the state vari- able can be used to estimate temperature from velocity eld measurements. Although, the data requirements for these estimations may be more. Recently, there have been extensive studies that use physics-informed neural networks (PINN) for ow reconstruction[79, 48]. In the context of ow reconstructions, a PINN is dened as a deep neural network that takes the coordinates of the collocation points as the input and outputs the two velocity eld components. A composite loss function is dened that incorporates the deviation from input snapshots (or the boundary conditions) and another loss from the PDE constraint. The PINN is trained to minimize this loss. This way, the reconstructed ow eld satises the boundary conditions as well as the physics, which has been incorporated via the PDE constraint. Using conventional PIV systems, it is very hard to get velocity elds in the near- wall region. Preliminary results show that a PINN can be trained with the PDE constraints discussed in chapter 2 to predict the velocity elds at higher temporal resolutions as well as the near-wall region. However, training a PINN can be computationally more expensive than evolving the TH model forwards in time. The models developed in this thesis are complementary to each other. Their predictions can be combined using data assimilation to generate better estimates for higher resolution ow elds. Another important limitation of high Re reconstructions is the limited availability of 121 both high-delity simulations as well as experimental data. However, these high-dimensional ow elds have a low-dimensional representation that enables sparse sensor placement to ex- tract the most information [58]. Previous studies have successfully reconstructed turbulent channel ows at Re = 1000 using a multi-sensor fusion approach[105]. They used mea- surements from a few time-resolved point measurements, non-time-resolved PIV snapshots, estimates from a dynamical model (RDT), all fused using a Kalman Filter. Thus any future research should also address the problem of sensor placement in addition to the underly- ing dynamical model to make data assimilation of high Re turbulent ows practical and cost-eective. Regime identication and ow reconstruction is expected to be much more challenging at higher Reynolds numbers. Although investigating higher Reynolds number ows is beyond the scope of this study, the framework developed here is likely to remain applicable as relevant measurement and simulation data become available. Note that, as the Reynolds number increases, the ow becomes more broad-banded. This means a wider range of ow features and a larger number of modes are required for reconstruction. Generating such features or modes requires more DNS simulations in the regime space (see gure 4.1). In other words, we need a richer database for reconstruction and estimation. However, DNS at high Re becomes increasingly dicult and time consuming as Re increases. To alleviate this problem, at least partially, the database may be composed of resolvent modes which can be directly obtained from the governing equations. Increasing Re does not appreciably in- crease computational cost for the resolvent modes. However, it is important to note that this approach poses problems of its own: resolvent mode computations are inexpensive for simple geometries and get increasingly more computationally intensive for complex geometries. The ideas developed in this thesis will potentially be useful to reconstruct higher res- olution ow elds, estimate unmeasured variables, or estimate the regime, all from limited measurements. The author hopes that the reliance on use of limited measurements makes the techniques developed in this thesis relevant for use in practical applications. 122 Bibliography [1] R. J. Adrian. \Hairpin vortex organization in wall turbulence". In: Phys. Fluids 19.4 (2007), p. 041301. [2] R. J. Adrian. \Particle-imaging techniques for experimental uid mechanics". In: Annu. Rev. Fluid Mech. 23.1 (1991), pp. 261{304. [3] R. J. Adrian. \Stochastic estimation of conditional structure: a review". In: Appl. Sci. Res. 53.3-4 (1994), pp. 291{303. [4] R. J. Adrian and P. Moin. \Stochastic estimation of organized turbulent structure: homogeneous shear ow". In: J. Fluid Mech. 190 (1988), pp. 531{559. [5] M. Alsalman, B. Colvert, and E. Kanso. \Training bioinspired sensors to classify ows". In: Bioinspir. Biomim. 14.1 (2018), p. 016009. [6] Filipe R Amaral et al. \Resolvent-based estimation of turbulent channel ow using wall measurements". In: Journal of Fluid Mechanics 927 (2021). [7] Nadine Aubry et al. \The dynamics of coherent structures in the wall region of a turbulent boundary layer". In: Journal of uid Mechanics 192 (1988), pp. 115{173. [8] H Jane Bae, Scott TM Dawson, and Beverley J McKeon. \Resolvent-based study of compressibility eects on supersonic turbulent boundary layers". In: Journal of Fluid Mechanics 883 (2020). [9] Z. Bai et al. \Low-dimensional approach for reconstruction of airfoil data via com- pressive sensing". In: AIAA J. 53.4 (2014), pp. 920{933. [10] G. K. Batchelor and I. Proudman. \The eect of rapid distortion of a uid in turbulent motion". In: Q. J. Mech. Appl. Math. 7.1 (1954), pp. 83{103. [11] S. Beneddine et al. \Unsteady ow dynamics reconstruction from mean ow and point sensors: an experimental study". In: J. Fluid Mech. 824 (2017), pp. 174{201. 123 [12] G. Berkooz, P. Holmes, and J. L. Lumley. \The proper orthogonal decomposition in the analysis of turbulent ows". In: Annu. Rev. Fluid Mech. 25.1 (1993), pp. 539{575. [13] M. G. Berry et al. \Low-dimensional and data fusion techniques applied to a super- sonic multistream single expansion ramp nozzle". In: Phys. Rev. Fluids 2.10 (2017), p. 100504. [14] S. A. Billings. Nonlinear system identication: NARMAX methods in the time, fre- quency, and spatio-temporal domains. John Wiley & Sons, 2013. [15] J. P. Bonnet et al. \Stochastic estimation and proper orthogonal decomposition: com- plementary techniques for identifying structure". In: Exp. Fluids 17.5 (1994), pp. 307{ 314. [16] K. A. Brucker and S. Sarkar. \A comparative study of self-propelled and towed wakes in a stratied uid". In: J. Fluid Mech. 652 (2010), pp. 373{404. [17] S. L. Brunton, J. L. Proctor, and J. N. Kutz. \Discovering governing equations from data by sparse identication of nonlinear dynamical systems". In: Proc. Nat. Acad. Sci. 113.15 (2016), pp. 3932{3937. [18] T. Bui-Thanh, M. Damodaran, and K. Willcox. \Aerodynamic data reconstruction and inverse design using proper orthogonal decomposition". In: AIAA J. 42.8 (2004), pp. 1505{1516. [19] E. J. Candes. \The restricted isometry property and its implications for compressed sensing". In: C R Math. 346.9-10 (2008), pp. 589{592. [20] K. K. Chen, J. H. Tu, and C. W. Rowley. \Variants of dynamic mode decomposi- tion: boundary condition, Koopman, and Fourier analyses". In: J. Nonlinear Sci. 22.6 (2012), pp. 887{915. [21] Vamsi Krishna Chinta and Mitul Luhar. \Statistically consistent resolvent-based re- construction of turbulent channel ows from limited measurements". In: Twelfth In- ternational Symposium on Turbulence and Shear Flow Phenomena, TSFP. 2022. [22] Vamsi Krishna Chinta et al. \Regime identication for stratied wakes from limited measurements: A library-based sparse regression formulation". In: Physical Review Fluids 7.3 (2022), p. 033803. [23] J. Choi, K. Yeo, and C. Lee. \Lagrangian statistics in turbulent channel ow". In: Phys. Fluids 16.3 (2004), pp. 779{793. 124 [24] J. M. Chomaz, P. Bonneton, and E. J. Hopnger. \The structure of the near wake of a sphere moving horizontally in a stratied uid". In: J. Fluid Mech. 254 (1993), pp. 1{21. [25] K. Chongsiripinyo, A. Pal, and S. Sarkar. \On the vortex dynamics of ow past a sphere at Re= 3700 in a uniformly stratied uid". In: Phys. Fluids 29.2 (2017), p. 020704. [26] D. R. Cole, M. N. Glauser, and Y. G. Guezennec. \An application of the stochastic estimation to the jet mixing layer". In: Phys. Fluids A 4.1 (1992), pp. 192{194. [27] B. Colvert, M. Alsalman, and E. Kanso. \Classifying vortex wakes using neural net- works". In: Bioinspir. Biomim. 13.2 (2018), p. 025003. [28] J. C. Del Alamo and J. Jim enez. \Estimation of turbulent convection velocities and corrections to Taylor's approximation". In: J. Fluid Mech. 640 (2009), pp. 5{26. [29] D. J. C. Dennis and T. B. Nickels. \On the limitations of Taylor's hypothesis in constructing long structures in a turbulent boundary layer". In: J. Fluid Mech. 614 (2008), pp. 197{206. [30] P. J. Diamessis, R. Gurka, and A. Liberzon. \Spatial characterization of vortical structures and internal waves in a stratied turbulent wake using proper orthogonal decomposition". In: Phys. Fluids 22 (2010), p. 086601. [31] S. Discetti, M. Raiola, and A. Ianiro. \Estimation of time-resolved turbulent elds through correlation of non-time-resolved eld measurements and time-resolved point measurements". In: Exp. Therm. Fluid Sci. 93 (2018), pp. 119{130. [32] M. P. Encinar and J. Jim enez. \Logarithmic-layer turbulence: A view from the wall". In: Phys. Rev. Fluids 4.11 (2019), p. 114603. [33] N. B. Erichson et al. \Shallow neural networks for uid ow reconstruction with limited sensors". In: Proc. R. Soc. Lond. A 476.2238 (2020), p. 20200097. [34] K. Fukami, K. Fukagata, and K. Taira. \Machine-learning-based spatio-temporal su- per resolution reconstruction of turbulent ows". In: J. Fluid Mech. 909 (2021). [35] K. Fukami, K. Fukagata, and K. Taira. \Super-resolution reconstruction of turbulent ows with machine learning". In: J. Fluid Mech. 870 (2019), pp. 106{120. 125 [36] F. G omez, A. S. Sharma, and H. M. Blackburn. \Estimation of unsteady aerodynamic forces using pointwise velocity data". In: J. Fluid Mech. 804 (2016). [37] F. G omez et al. \A reduced-order model of three-dimensional unsteady ow in a cavity based on the resolvent operator". In: J. Fluid Mech. 798 (2016). [38] J. Gra et al. \Reduced-order modeling for dynamic mode decomposition without an arbitrary sparsity parameter". In: AIAA J. 58.9 (2020), pp. 3919{3931. [39] J. Graham et al. \A web services accessible database of turbulent channel ow and its use for testing a new integral wall model for LES". In: J. Turbul. 17.2 (2016), pp. 181{215. [40] L. Guastoni et al. \Convolutional-network models to predict wall-bounded turbulence from wall quantities". In: J. Fluid Mech. 928 (2021). [41] A. Guemes et al. \From coarse wall measurements to turbulent velocity elds through deep learning". In: Phys. Fluids 33.7 (2021). [42] Y. G. Guezennec. \Stochastic estimation of coherent structures in turbulent boundary layers". In: Phys. Fluids A 1.6 (1989), pp. 1054{1060. [43] H. Hanazaki. \A numerical study of three-dimensional stratied ow past a sphere". In: J. Fluid Mech. 192 (1988), pp. 393{419. [44] A. E. Hoerl and R. W. Kennard. \Ridge regression: Biased estimation for nonorthog- onal problems". In: Technometrics 12.1 (1970), pp. 55{67. [45] Philip Holmes et al. Turbulence, coherent structures, dynamical systems and symme- try. Cambridge university press, 2012. [46] J. C. R. Hunt and D. J. Carruthers. \Rapid distortion theory and the `problems' of turbulence". In: J. Fluid Mech. 212 (1990), pp. 497{532. [47] S. J. Illingworth, J. P. Monty, and I. Marusic. \Estimating large-scale structures in wall turbulence using linear models". In: J. Fluid Mech. 842 (2018), p. 146. [48] Chiyu" Max" Jiang et al. \MeshfreeFlowNet: a physics-constrained deep continuous space-time super-resolution framework". In: Proceedings of the International Con- ference for High Performance Computing, Networking, Storage and Analysis. 2020, pp. 1{15. 126 [49] J. Jim enez and P. Moin. \The minimal ow unit in near-wall turbulence". In: J. Fluid Mech. 225 (1991), pp. 213{240. [50] J. Kim and F. Hussain. \Propagation velocity of perturbations in turbulent channel ow". In: Phys. Fluids A 5.3 (1993), pp. 695{706. [51] C. V. Krishna et al. \Reconstructing the time evolution of wall-bounded turbulent ows from non-time-resolved PIV measurements". In: Phys. Rev. Fluids 5.5 (2020), p. 054604. [52] P- A. Krogstad, J. H. Kaspersen, and S. Rimestad. \Convection velocities in a turbu- lent boundary layer". In: Phys. Fluids 10.4 (1998), pp. 949{957. [53] L. Lesshat et al. \Resolvent-based modeling of coherent wave packets in a turbulent jet". In: Phys. Rev. Fluids 4.6 (2019), p. 063901. [54] J. T. Lin and Y. H. Pao. \Wakes in Stratied Fluids". In: Annu. Rev. Fluid Mech. 11 (1979), pp. 317{338. [55] Q. Lin et al. \Stratied ow past a sphere". In: J. Fluid Mech. 240 (1992), pp. 315{ 354. [56] J. Luo et al. \Lagrangian dispersion in turbulent channel ow and its relationship to Eulerian statistics". In: Int. J. Heat and Fluid Flow 28.5 (2007), pp. 871{881. [57] T.J. Madison, X. Xiang, and G.R. Spedding. \Laboratory and numerical experiments on the near wake of a sphere in a stably stratied ambient". In: J. Fluid Mech. 933 (2022). [58] K. Manohar et al. \Data-driven sparse sensor placement for reconstruction: Demon- strating the benets of exploiting known patterns". In: IEEE Control Syst. Mag. 38.3 (2018), pp. 63{86. [59] B. J. McKeon. \Applications of resolvent analysis in uid mechanics". In: (2019). [60] B. J. McKeon. \The engine behind (wall) turbulence: perspectives on scale interac- tions". In: J. Fluid Mech. 817 (2017). [61] B. J. McKeon and A. S. Sharma. \A critical-layer framework for turbulent pipe ow". In: J. Fluid Mech. 658 (2010), pp. 336{382. [62] P. Meunier and G. R. Spedding. \A loss of memory in stratied momentum wakes". In: Phys. Fluids 16.2 (2004), pp. 298{305. 127 [63] R. Moarref et al. \Model-based scaling of the streamwise energy density in high- Reynolds-number turbulent channels". In: J. Fluid Mech. 734 (2013), pp. 275{316. [64] P. Moin. \Revisiting Taylor's hypothesis". In: J. Fluid Mech. 640 (2009), pp. 1{4. [65] Parviz Moin and Krishnan Mahesh. \Direct numerical simulation: a tool in turbulence research". In: Annual review of uid mechanics 30.1 (1998), pp. 539{578. [66] J. C. Montgomery, C. F. Baker, and A. G. Carton. \The lateral line can mediate rheotaxis in sh". In: Nature 389.6654 (1997), pp. 960{963. [67] Pierluigi Morra et al. \On the relevance of Reynolds stresses in resolvent analyses of turbulent wall-bounded ows". In: Journal of Fluid Mechanics 867 (2019), pp. 969{ 984. [68] N. E. Murray and L. S. Ukeiley. \Estimation of the oweld from surface pressure measurements in an open cavity". In: AIAA J. 41.5 (2003), pp. 969{972. [69] N. E. Murray and L. S. Ukeiley. \Modied quadratic stochastic estimation of resonat- ing subsonic cavity ow". In: J. Turbul. 8 (2007), N53. [70] A. M. Naguib, C. E. Wark, and O. Juckenh ofel. \Stochastic estimation and ow sources associated with surface pressure events in a turbulent boundary layer". In: Phys. Fluids 13.9 (2001), pp. 2611{2626. [71] S. Nidhan et al. \Dynamic mode decomposition of stratied wakes". In: AIAA Avia- tion 2019 Forum. 2019, p. 3330. [72] C. Y. Ohh and G. Spedding. \Automated stratied wake classication using Dynamic Mode Decomposition". In: Bulletin of the American Physical Society (2020). [73] Chan-Ye Ohh and Georey R Spedding. \Wake identication of stratied ows using dynamic mode decomposition". In: Physical Review Fluids 7.2 (2022), p. 024801. [74] T. S. Orr et al. \Numerical simulations of the near wake of a sphere moving in a steady, horizontal motion through a linearly stratied uid at Re= 1000". In: Phys. Fluids 27.3 (2015), p. 035113. [75] A. Pal et al. \Direct numerical simulation of stratied ow past a sphere at a sub- critical Reynolds number of 3700 and moderate Froude number". In: J. Fluid Mech. 826 (2017), pp. 5{31. 128 [76] A. Pal et al. \Regeneration of turbulent uctuations in low-Froude-number ow over a sphere at a Reynolds number of 3700". In: J. Fluid Mech. 804 (2016). [77] F. Papi and P. Luschi. \Pinpointing'Isla Meta': the case of sea turtles and albatrosses". In: J. Exp. Biol. 199.1 (1996), pp. 65{71. [78] M. Quadrio and P. Luchini. \Integral space{time scales in turbulent wall ows". In: Phys. Fluids 15.8 (2003), pp. 2219{2227. [79] Maziar Raissi, Alireza Yazdani, and George Em Karniadakis. \Hidden uid me- chanics: Learning velocity and pressure elds from ow visualizations". In: Science 367.6481 (2020), pp. 1026{1030. [80] J. A. Redford, T. S. Lund, and G. N. Coleman. \A numerical study of a weakly stratied turbulent wake". In: J. Fluid Mech. 776 (2015), pp. 568{609. [81] L. Ristroph, J. C. Liao, and J. Zhang. \Lateral line layout correlates with the dier- ential hydrodynamic pressure on swimming sh". In: Phys. Rev. Lett. 114.1 (2015), p. 018102. [82] Patrick J Roache. \On articial viscosity". In: Journal of Computational Physics 10.2 (1972), pp. 169{184. [83] Clarence W Rowley, Tim Colonius, and Richard M Murray. \Model reduction for compressible ows using POD and Galerkin projection". In: Physica D: Nonlinear Phenomena 189.1-2 (2004), pp. 115{129. [84] S. H. Rudy et al. \Data-driven discovery of partial dierential equations". In: Sci. Adv. 3.4 (2017), e1602614. [85] K. Sasaki et al. \Transfer functions for ow predictions in wall-bounded turbulence". In: J. Fluid Mech. 864 (2019), pp. 708{745. [86] A. M. Savill. \Recent developments in rapid-distortion theory". In: Annu. Rev. Fluid Mech. 19.1 (1987), pp. 531{573. [87] P. J. Schmid. \Dynamic mode decomposition of numerical and experimental data". In: J. Fluid Mech. 656 (2010), pp. 5{28. [88] A. J. Smits, B. J. McKeon, and I. Marusic. \High{Reynolds number wall turbulence". In: Annu. Rev. Fluid Mech. 43 (2011). 129 [89] G. R. Spedding. \Anisotropy in turbulence proles of stratied wakes". In: Phys. Fluids 13.8 (2001), pp. 2361{2372. [90] G. R. Spedding. \Wake signature detection". In: Annu. Rev. Fluid Mech. 46 (2014), pp. 273{302. [91] S. Symon, D. Sipp, and B. J. McKeon. \A tale of two airfoils: resolvent-based mod- elling of an oscillator versus an amplier from an experimental mean". In: J. Fluid Mech. 881 (2019), pp. 51{83. [92] K. Taira et al. \Modal analysis of uid ows: An overview". In: AIAA J. 55.12 (2017), pp. 4013{4041. [93] G. I. Taylor. \The spectrum of turbulence". In: Proc. Royal Soc. London. 164.919 (1938), pp. 476{490. [94] J. A. Taylor and M. N. Glauser. \Towards practical ow sensing and control via POD and LSE based low-dimensional tools". In: J. Fluid Eng-T ASME 126.3 (2004), pp. 337{345. [95] Robert Tibshirani. \Regression shrinkage and selection via the lasso". In: J. R. Stat. Soc. Series B. Stat. Methodol. 58.1 (1996), pp. 267{288. [96] A. Towne, A. Lozano-Dur an, and X. Yang. \Resolvent-based estimation of space{time ow statistics". In: J. Fluid Mech. 883 (2020). [97] A. Towne, O. T. Schmidt, and T. Colonius. \Spectral proper orthogonal decomposi- tion and its relationship to dynamic mode decomposition and resolvent analysis". In: J. Fluid Mech. 847 (2018), pp. 821{867. [98] J. H. Tu et al. \Integration of non-time-resolved PIV and time-resolved velocity point sensors for dynamic estimation of velocity elds". In: Exp. Fluids 54.2 (2013), pp. 1{ 20. [99] J. H. Tu et al. \On dynamic mode decomposition: Theory and applications". In: arXiv preprint arXiv:1312.0041 (2013). [100] J. H. Tu et al. \Spectral analysis of uid ows using sub-Nyquist-rate PIV data". In: Exp. Fluids 55.9 (2014), p. 1805. [101] L. Ukeiley et al. \Dynamic surface pressure based estimation for ow control". In: IUTAM Symposium on Flow Control and MEMS. Springer. 2008, pp. 183{189. 130 [102] L Ukeiley et al. \Examination of large-scale structures in a turbulent plane mixing layer. Part 2. Dynamical systems model". In: Journal of Fluid Mechanics 441 (2001), pp. 67{108. [103] S. B. Vadarevu et al. \Coherent structures in the linearized impulse response of tur- bulent channel ow". In: J. Fluid Mech. 863 (2019), pp. 1190{1203. [104] M. Wang and T. A. Zaki. \State estimation in turbulent channel ow from limited observations". In: J. Fluid Mech. 917 (2021). [105] M. Wang et al. \Model-based multi-sensor fusion for reconstructing wall-bounded turbulence". In: Theor. Comput. Fluid Dyn. 35 (2021), pp. 683{701. [106] J. Westerweel, G. E. Elsinga, and R. J. Adrian. \Particle image velocimetry for com- plex and turbulent ows". In: Annu. Rev. Fluid Mech. 45 (2013), pp. 409{436. [107] K. Willcox. \Unsteady ow sensing and estimation via the gappy proper orthogonal decomposition". In: Comput. Fluids 35.2 (2006), pp. 208{226. [108] C. H. K. Williamson and A. Roshko. \Vortex formation in the wake of an oscillating cylinder". In: J. Fluids Struct. 2.4 (1988), pp. 355{381. [109] J. C. Wyngaard and S. F. Cliord. \Taylor's hypothesis and high{frequency turbu- lence spectra". In: J. Atmos. Sci. 34.6 (1977), pp. 922{929. [110] X. Xiang, K. K. Chen, and G. R. Spedding. \Dynamic mode decomposition for esti- mating vortices and lee waves in a stratied wake". In: Exp. Fluids 58.5 (May 2017), p. 56. [111] X. Xiang et al. \The turbulent wake of a towed grid in a stratied uid". In: J. Fluid Mech. 775 (2015), pp. 149{177. [112] X. I. A. Yang and M. F. Howland. \Implication of Taylor's hypothesis on measuring ow modulation". In: J. Fluid Mech 836 (2018), pp. 222{237. 131

Abstract (if available)

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Pattern generation in stratified wakes

PDF

Model based design of porous and patterned surfaces for passive turbulence control

PDF

Near wake characteristics of towed bodies in a stably stratified fluid

PDF

Physics-based and data-driven models for bio-inspired flow sensing and motion planning

PDF

Modeling and analysis of parallel and spatially-evolving wall-bounded shear flows

PDF

Aerodynamics at low Re: separation, reattachment, and control

PDF

Flow and thermal transport at porous interfaces

PDF

Data-driven multi-fidelity modeling for physical systems

PDF

Physics-based data-driven inference

PDF

An experimental study of the elastic theory of granular flows

PDF

Supersonic flow study of a capsule geometry using large-eddy simulation

PDF

A Boltzmann model for tracking aerosols in turbulent flows

PDF

Re-assessing local structures of turbulent flames via vortex-flame interaction

PDF

Passive flight in density-stratified fluid environments

PDF

On the simulation of stratified turbulent flows

PDF

Large eddy simulations of turbulent flows without use of the eddy viscosity concept

PDF

Wake modes of rotationally oscillating circular cylinder in cross-flow and its relationship with heat transfer

PDF

Improved brain dynamic contrast enhanced MRI using model-based reconstruction

PDF

Large eddy simulations of laminar separation bubble flows

PDF

RANS simulations for flow-control study synthetic jet cavity on NACA0012 and NACA65(1)412 airfoils.

Asset Metadata

Creator Chinta, Vamsikrishna (author)

Core Title Reconstruction and estimation of shear flows using physics-based and data-driven models

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Mechanical Engineering

Degree Conferral Date 2022-12

Publication Date 10/18/2022

Defense Date 08/10/2022

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag data assimilation,flow reconstruction,machine learning,OAI-PMH Harvest,regime estimation,stratified geophysical flows,turbulence,wall-bounded flows

Format theses (aat)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Luhar, Mitul (committee chair), Oberai, Assad (committee member), Savla, Ketan (committee member), Spedding, Geoffrey (committee member)

Creator Email vamsi482me@gmail.com,vchinta@usc.edu

Permanent Link (DOI) https://doi.org/10.25549/usctheses-oUC112124636

Unique identifier UC112124636

Identifier etd-ChintaVams-11274.pdf (filename)

Legacy Identifier etd-ChintaVams-11274

Document Type Dissertation

Format theses (aat)

Rights Chinta, Vamsikrishna

Internet Media Type application/pdf

Type texts

Source 20221019-usctheses-batch-987 (), University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the author, as the original true and official version of the work, but does not grant the reader permission to use the work if the desired use is covered by copyright. It is the author, as rights holder, who must provide use permission if such use is covered by copyright. The original signature page accompanying the original submission of the work to the USC Libraries is retained by the USC Libraries and a copy of it may be obtained by authorized requesters contacting the repository e-mail address given.

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA

Repository Email cisadmin@lib.usc.edu