Prediction of the Optimal Vortex in Synthetic Jets

Abstract

This article presents three different low-order models to predict the main flow patterns in synthetic jets. The first model provides a simple theoretical approach based on experimental solutions explaining how to artificially generate the optimal vortex, which maximizes the production of thrust and system efficiency. The second model is a data-driven method that uses higher-order dynamic mode decomposition (HODMD). To construct this model, (i) Navier–Stokes equations are solved for a very short period of time providing a transient solution, (ii) a group of spatio-temporal data are collected containing the information of the transitory of the numerical simulations, and finally (iii) HODMD decomposes the solution as a Fourier-like expansion of modes that are extrapolated in time, providing accurate predictions of the large size structures describing the general flow dynamics, with a speed-up factor of 8.3 in the numerical solver. The third model is an extension of the second model, which combines HODMD with a low-rank approximation of the spatial domain, which is based on singular value decomposition (SVD). This novel approach reduces the memory requirements by 70% and reduces the computational time to generate the low-order model by 3, maintaining the speed-up factor to 8.3. This technique is suitable to predict the temporal flow patterns in a synthetic jet, showing that the general dynamics is driven by small amplitude variations along the streamwise direction. This new and efficient tool could also be potentially used for data forecasting or flow pattern identification in any type of big database.

1. Introduction

Environmental problems such as the greenhouse effect or the global warming combined with the running out of fossil fuels, encourage researchers to look for alternative energy sources and transportation devices minimizing the environmental impact. During recent years, the continuous search for developing alternative propulsion systems inspired by, or even mimicking, animal motion, has become a research topic of high interest [1,2], especially in the field of marine locomotion [3,4,5]. Propulsion systems found in nature are mainly driven by a vortex ring producing thrust [5]. For example, the key physical feature of animal swimming motion lies in the presence of an optimal vortex behind their bodies. The formation and evolution of such vortex ring distinguish the efficiency of the propulsion system and the amount of impulse generated by the fluid producing thrust.

Synthetic jets are devices able to artificially create vortex rings modelling the formation and transition of the optimal vortex driving the animal swimming motion [6]. A synthetic jet is a fluid stream formed by a group of vortex rings ejected periodically from a cavity. The cavity contains a membrane or a piston oscillating and forcing the fluid to leave and go into the same cavity through an orifice. The net momentum flux that passes through the jet orifice is non-zero, although the net mass flux is zero. Synthetic jets produce non-zero mean streamwise momentum without the need for additional mass injection [7], modelling the swimming motion of some animals such as squids, jellyfish, or salps [8], but also this fact is very attractive for industry using synthetic jet devices for fluid mixing [9], heat transfer enhancement [10], active flow control of boundary layer separation [11,12], or for plasma actuators [13].

Synthetic jets distinguish two flow topologies depending on the oscillation phase, generally periodic, of the device (membrane or piston) found in the cavity upstream the jet nozzle. When the flow is ejected from the cavity, the flow separates and rolls up forming a vortex ring moving downstream due to self-induced velocity, and producing thrust. This is known as the injection phase. On the contrary, when the flow enters the cavity again, a saddle point separates the flow re-entering into the cavity from the one that continues moving downstream. This is known as the suction phase [14]. Depending on two non-dimensional parameters, the Reynolds number, which compares convective and viscous terms, and the Strouhal number, measuring oscillation frequencies, Carter & Soria [15] distinguished four different flow regimes, each one related to the generation of a different type of vortex. These regimes include the generation of (i) laminar jets, (ii) vortex rings, (iii) transitional jets and (iv) turbulent jets. The vortex ring presented in regimes (i) and (ii) is generally driving the swimming motion of jellyfish or squids [16]. Smith & Glezer shows that the vortex leaves the cavity, travels downstream, and finally breaks down, suggesting that the mechanism triggering such break down could be an instability occurring along the azimuthal spatial component. The location of the vortex breakdown is modified by the flow conditions, thus in regimes (i) and (ii) the vortex travels downstream and breaks down in the far field, while in regimes (iii) and (iv) the vortex ring breaks down in the near field of the jet, being this vortex almost undistinguishable in regime (iv). It is remarkable that the Reynolds number and the Strouhal number are crucial to categorize the type of vortex formation, when the shape of the cavity is rounded [17] (the length and diameter of the cavity does not play any fundamental role on these flow regimes, although the shape of the cavity influences the type of vortex [18]). Holman et al. [19] identified the critical value of these parameters to guarantee the jet formation, for values upper or below this criterion, the flow continuously leaves and enter the cavity, blocking the vortex travelling upstream.

Due to their close connections with the animal swimming motion, researchers have studied in detail the generation of vortices described in regimes (i) and (ii) [20]. Dabiri et al. [16] studied carefully the vortex producing thrust in the swimming motion of seven different types of jellyfish. The authors found that the jellyfish producing the largest thrust were swimming in regime (i), while the highest efficient jellyfish were swimming in regime (ii). The key goal for optimal propulsion is finding how to generate artificially the optimal vortex that maximize the efficiency of the system producing the largest thrust. Gharib et al. [6] studied in detail the optimal time for vortex formation. The authors found that the generation of an optimal vortex, defined as the upper limit to acquire the maximum circulation by the vortex itself, was dependent on the bulk length ($L_b$) and jet nozzle diameter (D). In particular, they found in a wide range of cases that this optimal vortex was formed for values of $L_b/D\simeq 4$ at Reynolds number $\simeq {10}^4$. For simplicity, in this paper the ratio $L_b/D$ is known as the vortex formation number. Afterwards, Krueger [21] determined experimentally that there was a specific time (non-dimensional) for the vortex formation were the system was producing the highest thrust. Such time was in good agreement with the definition of the generation of the optimal vortex, formed for values of $L_b/D\simeq 4$ (at Reynolds number $\simeq {10}^4$). In other words, he found that it was possible to optimize the efficiency of momentum transport by maximizing the size of the vortex ring.

Based on the knowledge of the existence of an optimal vortex for optimal propulsion, it is necessary to find the way to easily generate this vortex artificially using synthetic jets. Bandyopadhyay & Beal [22] related the optimal vortex formation with the Strouhal number. However, these authors demonstrated that in some cases the frequency is insufficient to predict the optimal propulsion, making more difficult to relate the results presented by Carter & Soria (regimes (i) and (ii)) and the definition of the optimal vortex presented by Gharib et al. [6].

The main goal of this work is to present a simple way to predict the optimal vortex and a new method for data forecasting to reduce the computational cost in numerical simulations describing synthetic jet flows. The article presents three methodologies for data forecasting with the aim at identifying the main global patterns and flow structures characteristic of this type of flow. On the one hand, the article introduces a new (to the author knowledge) theoretical approach to predict the optimal set-up of the device (St and Re) to generate artificially the optimal vortex (maximum thrust and efficiency), which is validated with results obtained with numerical simulations. On the other hand, the article will also introduce a data-driven method to create a reduced order model (ROM) to predict the main flow dynamics in a synthetic jet, with the aim at reducing the computational cost in numerical simulations. The model is based on dynamic mode decomposition (DMD) [23] and higher-order DMD (HODMD) [24] (an extended version of the method), for simplicity known as DMD-based ROM. The good performance of this method has already been tested for data forecasting in the three-dimensional wake of a cylinder [25], in the turbulent wake of a wind turbine [26], and in the compressible flow around a profile with a shock wave subject to buffeting condition [27]. This DMD-based ROM will be combined with a low-rank approximation of the spatial computational domain, reducing in 70% the memory size of the database analyzed and reducing in a factor of 3 the computational time. This new method, for simplicity known as DMD low rank-ROM (DMDlr-ROM), presented for the first time in this article (to the author knowledge), could be potentially used for the analysis of turbulent flow databases [28], generally formed by large size memory files, requiring the latest technology computer facilities for processing the data. In this article, both the DMD-based ROM and the DMDlr-ROM have been applied for data forecasting in numerical simulations, obtaining the main global patterns describing a synthetic jet flow with relative errors ~3%, presenting a speed-up factor of the numerical simulations of ~8.

The article is organized as follows. Section 2 introduces the model. Section 3 and Section 4 presents the algorithm for the DMD-based ROM and DMDlr-ROM. Section 5 describes a theoretical approach to predict the artificial generation of the optimal vortex. Finally, the main results and conclusions are presented in Section 6 and Section 7.

2. Model Description

Synthetic jets are characterized by the jet velocity scales (U), the time scales (T) and the length scales, which are determined by the bulk length ($L_b$) and the diameter of the jet nozzle (D). Such flow scales are determined by two parameters, the Reynolds (Re) and Strouhal (St) numbers, defined as $Re=UD/\nu$ and $St=fD/U$, respectively, where D is the diameter of the jet orifice, $\nu$ the kinematic viscosity, f is the piston oscillation frequency and U is the characteristic velocity scale represented by the mean velocity at the jet orifice. This parameter is based on the mean momentum velocity, defined as

$$U={\left[\frac{1}{AT}{\int }_0^T{\int }_0^{D/2}{\int }_0^{2\pi }2\pi rv(r,\theta ,t)d\theta drdt\right]}^{1/2}.$$

(1)

where T is the piston oscillation period ($T=\frac{1}{f}$), $A=\pi D^2/4$ is the area of the jet orifice, and $v(r,\theta ,t)$ is the streamwise velocity at the exit of the orifice, $x=0$ [15]. It is remarkable that in this article the jet velocity is defined in a cylindrical coordinate system, being r, x, and $\theta$ the radial, streamwise, and azimuthal coordinates, respectively. The solution of the incompressible Navier–Stokes equations provides the value of $v(r,\theta ,t)$.

The average axial momentum flow across the jet orifice is given by

$$AU=\frac{A_p{\widehat{U}}_p}{2},$$

(2)

where $A_p=\pi D_p^2/4$ is the area of the piston, with diameter $D_p$, and ${\widehat{U}}_p$ is the peak amplitude of the piston velocity, defined as $v_p\left(t\right)={\widehat{U}}_{p}sin\left(\frac{2\pi t}{T}\right)$. This assumption provides a simple way to define the mean velocity at the jet orifice, as

$$U=\frac{D_{p}D}{\sqrt{2}\widehat{U_p}},$$

(3)

see [14,29] for more details.

Numerical simulations have been carried out to model the axi-symmetric flow generated by a synthetic jet with a cylindrical cavity and circular jet nozzle at Re = 1000 and St = 0.02, 0.022 and 0.03. The database generated numerically in [30] at St 0.03 is used in this article and two new data bases have been generated to study the cases at St 0.02 and 0.022 using the same numerical code and working conditions, briefly presented below. Based on the diagram presented by Carter & Soria [15] and the results presented in [14,30,31] the flow is expected to be axi-symmetric at this flow conditions (regimes (i) or (ii) of the vortex generation in [15]). Non-linear Navier–Stokes equations have been solved using the open source numerical solver Nek5000 [32]. This solver uses spectral elements as spatial discretization discretizing each macro-element conforming the computational domain with Gauss-Lobatto-Legendre points of polynomial order $\Pi$. The temporal discretization uses an explicit second order extrapolation scheme for the non-linear terms and an implicit second order backwards differentiation scheme for the viscous terms (see [33]).

Figure 1 shows the computational domain of the numerical simulations, where axi-symmetric condition is imposed in the bottom side of the domain. The remaining boundary conditions and dimensions of the computational domain are set as in [30,34]. Regarding the dimensions of the domain: the diameters of the jet nozzle and the piston are set as $D=1$ and $D_p=5D$, respectively; the cavity length, jet stroke, cavity radius and jet orifice radius are defined as $L_c=5D$, $L_o=0.2D$, $H_c=D_p/2=2.5D$ and $R=0.5D$, respectively; finally, to avoid boundary reflections the length and radius of the computational domain are set as $L=80D$ and $H=40D$, respectively. Regarding the boundary conditions: the inlet boundary (piston) is defined as $U_x=v_p\left(t\right)={\widehat{U}}_p·sin\left(2\pi ft\right)$ and $U_r=0$ for the streamwise and radial velocity components, respectively, and Neumann boundary conditions for the pressure; zero-stress and Dirichlet boundary conditions are imposed for the velocity ($\nabla \mathbf{u}·\mathbf{n}$ = 0, with $\mathbf{n}=$ unit normal) and pressure ($p=0$), respectively, in the top and outflow surfaces of the domain, while non-slip boundary conditions ($U_x=U_r=0$) are imposed in the wall of the jet and left surfaces of the domain.

A grid independence study and the validation of the numerical simulations with experiments have been carried out in [30]. The simulations have been performed varying the polynomial order as $\Pi =4,6,8,10$. The dominant frequency calculated in several points of the computational domain is St $=0.03$ (similar to the experimental results) in the cases with $\Pi =8$ and 10, varying in the 6th decimal digit, thus the polynomial order is set as $\Pi =8$.

3. The DMD-Based ROM: A Reduced Order Model Based on HODMD

Higher-order dynamic mode decomposition (HODMD) [24] is an extension of DMD [23] recently introduced for the analysis of noisy experimental data, non-linear dynamical systems and complex flows [35,36,37]. As with DMD, HODMD decomposes spatio-temporal data collected at time $t_k$, ${\mathit{v}}_k=\mathit{v}(x,r,t_k)$, as an expansion of M Fourier-like modes as

$${\mathit{v}}_k\simeq {\displaystyle \sum _{m=1}^M}{a}_m{\mathbf{u}}_m(x,r)e^{({\delta }_m+i{\omega }_m)t_k},k=1,\dots ,K.$$

(4)

where ${\mathit{u}}_m(x,r)$ are the DMD modes, weighted with the amplitudes $a_m$, oscillating and growing/decaying in time with frequency ${\omega }_m$ and growth rate ${\delta }_m$. HODMD is a data-driven method, applied to analyze some data collected in the time interval $t\in [t_1,t_K]$ to construct the previous DMD expansion. Nevertheless, for accurate calculations of such expansion, it is possible to either reconstruct the original data sequence (interpolation) or to predict temporal events (extrapolation), by means of adjusting the temporal term $t_k=t_{\xi }\gg t_K$. This is the basic idea behind using HODMD as a ROM for data forecasting, with the aim at reducing the computational time in numerical simulations [25]. In other words, (i) a set of data are collected in the transitory of a numerical simulation, (ii) HODMD is applied to analyze the data and to construct the DMD expansion (4), (iii) to predict the attractor, only the permanent modes are retained to construct the previous expansion, these are the modes whose growth rate is close to zero (in the transitory all the modes either grow or decay in time, but the real dynamics is represented by a group of modes with relatively small growth rate, see [25]), (iv) the temporal term $t_k$ is adjusted to predict the attractor. A good way to select the permanent modes composing (4) is to establish a criterion based on a tolerance $\gamma$ for the growth rates, thus the permanent modes are those whose growth rate is defined as $|{\delta }_m|<\gamma$. Once the DMD expansion (4) is constructed, the growth rate of such modes could be forced to be exactly equal to zero, ${\delta }_m=0$, although for small values of $\gamma$, it is possible to maintain the real growth rate of such modes [26,27].

Algorithm

HODMD algorithm is briefly introduced in this section. A more detailed description can be found in Ref. [24]. For simplicity, before proceeding with HODMD algorithm, a group of K spatio-temporal snapshots are collected in the following snapshot matrix as

$${\mathit{V}}_1^K=[{\mathit{v}}_1,{\mathit{v}}_2,\dots ,{\mathit{v}}_k,{\mathit{v}}_{k+1},\dots ,{\mathit{v}}_{K-1},{\mathit{v}}_K].$$

(5)

The dimension of the matrix is $J\times K$, with $J=N_x\times N_r$, being $N_x$ and $N_r$ the number of grid points defined along the streamwise and radial spatial components of the domain. The data are equidistant in time, with time interval $\Delta t$. In this problem, each snapshot represents the velocity vectors ${\mathit{v}}_x$ and ${\mathit{v}}_r$ corresponding to the streamwise and radial spatial directions, organized in columns.

HODMD algorithm can be encompassed in three main steps.

• Step 1: SVD. With the aim at removing spatial redundancies or noise, SVD is applied to the snapshot matrix (5), as

  $${\mathit{V}}_1^K\simeq \mathit{U}\Sigma {\mathit{T}}^{\top }=\mathit{U}{\widehat{\mathit{V}}}_1^K,$$

  (6)

  where ${\mathit{U}}^T\mathit{U}={\mathit{T}}^T\mathit{T}=$ the $N\times N$ unit matrix, the diagonal of matrix $\Sigma$ contains the singular values ${\sigma }_1,\cdots ,{\sigma }_N$, and ${\widehat{\mathit{V}}}_1^K$ is the reduced snapshot matrix, of dimension $N\times K$. Based on a tolerance $\epsilon$ (set by the user), it is possible to determine the number of N singular values retained in the previous equation, as

  $$\frac{{\sigma }_{N+1}}{{\sigma }_1}\le \epsilon .$$

  (7)

  At this step, the spatial complexity of the data J is reduced to a set of N linearly independent vectors, determining the spatial dimension of the system, generally $N\le J$.
  In the analysis of complex data, SVD is replaced by a high-order singular value decomposition (HOSVD) [38], leading to the called multi-resolution HODMD algorithm. This decomposition organizes the original snapshot matrix $V_1^K$ in tensor form as $V(l,x_i,r_j,t_k)=V_{lijk}$, with $l=1,2$ being the number of velocity components (${\mathit{v}}_x$ and ${\mathit{v}}_r$), for $i=1,\cdots ,N_x$; $j=1,\cdots ,N_r$ and $k=1,\cdots ,K$. HOSVD algorithm applies an SVD to the three matrices whose columns are formed by each one of the 3 data variables (similar to the fibers of a tensor). Thus, this method provides the following decomposition for each velocity component l (for simplicity the solution is particularized for a single velocity component)

  $$V_{ijk}\simeq {\displaystyle \sum _{p_1=1}^{P_1}}{\displaystyle \sum _{p_2=1}^{P_2}}{\displaystyle \sum _{n=1}^N}{S}_{p_1{p}_{2}n}{W}_{i{p}_1}^{\left(x\right)}{W}_{j{p}_2}^{\left(r\right)}{T}_{kn},$$

  (8)

  where ${\mathit{S}}_{p_1{p}_2}$ is a third-order tensor (called the core tensor) and the columns of the matrices ${\mathit{W}}^{\left(x\right)}$, ${\mathit{W}}^{\left(r\right)}$ and $\mathit{T}$ are known as the modes of the decomposition (two spatial and one temporal modes, respectively). The rank reduction in Equation (7) is then applied to each one of these modes. The temporal modes $\mathit{T}$ are used to construct the reduced snapshot matrix. Due to the complexity of the data analyzed in this article, this multi-resolution algorithm is used to construct the DMD-based ROM (see more details about the algorithm in [14]).
• Step 2: the DMD-d algorithm. DMD algorithm [23] is combined with the Takens’ delay embedding theorem [39], providing the following high order Koopman assumption that is applied to the reduced snapshot matrix as

  $${\widehat{\mathit{V}}}_{d+1}^K\simeq {\widehat{\mathit{R}}}_1{\widehat{\mathit{V}}}_1^{K-d}+{\widehat{\mathit{R}}}_2{\widehat{\mathit{V}}}_2^{K-d+1}+\dots +{\widehat{\mathit{R}}}_d{\widehat{\mathit{V}}}_d^{K-1}.$$

  (9)
  The dynamics of the system is contained in the reduced Koopman operators ${\widehat{\mathit{R}}}_1$, ${\widehat{\mathit{R}}}_2$, ⋯, ${\widehat{\mathit{R}}}_d$. These operators are encompassed into a single matrix. The eigenvalues of such matrix provide the frequencies ${\omega }_m$ and growth rates ${\delta }_m$ of the DMD expansion (4), and the eigenvectors are used to construct the DMD modes ${\mathit{u}}_m$. It is remarkable that for values of $d=1$, this algorithm is similar to standard DMD [23], providing similar solutions. Nevertheless, in the analysis of complex dynamical system, composed by large number of frequencies, transient modes, or even noisy data, values of $d>1$ will provide accurate approximations of the expansion (4), required to construct a ROM and to predict the attractor (note that the errors grow exponentially with time).
• Step 3: Amplitude calculations and RRMS error. The amplitude $a_m$ weighting the DMD modes are calculated in this step. For proper reconstruction of the original data (interpolation), a good way to calculate the amplitudes is by least squares fitting of the expansion (4). This amplitude will determine the number of M DMD modes retained in the expansion as function of a second tolerance ${\epsilon }_1$, as

  $$a_{M+1}/a_1<{\epsilon }_1.$$

  (10)
  Once the DMD expansion (4) is approximated, the relative root mean square error (RRMSE) is used to quantify the difference between the original data and the approximated solution (4), as

  $$RRMSE\sim \frac{\Vert {\mathit{V}}_{original}-{\mathit{V}}_{DMD}{\Vert }_F}{\Vert {\mathit{V}}_{original}{\Vert }_F},$$

  (11)

  where $\Vert ·{\Vert }_F$ represents the Frobenius norm.
  To minimize this error, the method is applied iteratively. In other words, the algorithm is applied over the solution reconstructed in (4), obtaining new reconstructions, until the number of singular values N is maintained similar in two consecutive iterations. It is possible to download a MATLAB [40] version of the HODMD algorithms in Ref. [41].

4. The DMDlr-ROM: A Low-Rank Spatial Approximation for the Data-Driven ROM

Generally, the storage and analysis of databases describing complex or turbulent flows implies using sophisticated and powerful computers. The size of each snapshot composing the snapshot matrix (5) varies from a few megabytes to hundred gigabytes in data generated in numerical simulations describing turbulent flows [28], making prohibitively expensive using HODMD algorithm for the analysis of flow structures. This section presents a low-rank approximation of the snapshot matrix (5) to construct the data-driven ROM based on HODMD previously introduced, which will reduce the computational cost. This new algorithm, presented for the first time in this article to the author knowledge, is known as DMD low-rank ROM (DMDlr-ROM). The methodology can be encompassed in three main Stages.

• Stage 1: rank reduction and generation of the new snapshot matrix. SVD is applied to reduce the spatial dimension of the original snapshots. In contrast to Step 1 from HODMD algorithm, SVD is applied separately to each one of the snapshots ${\mathit{v}}_k$ composing the snapshot matrix (5). Thus, the information contained in each snapshot ${\mathit{v}}_k={\mathit{v}}_k(x,r)$ is organized in tensor form as ${\mathit{v}}_k(x_i,r_j)={\mathit{v}}_{kij}$ for $i=1,\cdots ,N_x$ and $j=1,\cdots ,N_r$. For each snapshot k the SVD leads as

  $${\mathit{v}}_{kij}={\mathit{U}}_k{\Sigma }_k{\mathit{T}}_k^{\top },$$

  (12)

  where the dimension of ${\mathit{U}}_k$ is $N_x\times N$, the dimension of ${\mathit{T}}_k^{\top }$ is $N\times N_r$ and the dimension of ${\Sigma }_k$ is $N\times N$. This decomposition has been carried out using the function svd in MATLAB [40]. The spatial matrices ${\mathit{U}}_k$ and ${\mathit{T}}_k$ contain information related to the streamwise and radial spatial components, respectively, and ${\Sigma }_k$ is the matrix containing the singular values ${\sigma }_{k1},{\sigma }_{k2},\cdots ,{\sigma }_{kN}$. The spatial matrices are weighted with the square root of the singular values as

  $${\overline{\mathit{U}}}_k={\mathit{U}}_k\overline{\Sigma },$$

  (13)

  and

  $${\overline{\mathit{T}}}_k=\overline{\Sigma }{\mathit{T}}_k^{\top },$$

  (14)

  being $\overline{\Sigma }$ the matrix composed by the square root of the singular values $\sqrt{{\sigma }_{k1}},\sqrt{{\sigma }_{k2}},\cdots ,\sqrt{{\sigma }_{kN}}$.
  The same process is repeated for each one of the K snapshots of the initial snapshot matrix (5) for each velocity component ${\mathit{v}}_x$ and ${\mathit{v}}_y$. Two new snapshot matrices are constructed with the results previously obtained, the called X-snapshots defined as

  $${\mathit{X}}_1^K=[{\overline{\mathit{U}}}_1,{\overline{\mathit{U}}}_2,\cdots ,{\overline{\mathit{U}}}_K],$$

  (15)

  and the called Y-snapshots, defined as

  $${\mathit{Y}}_1^K=[{\overline{\mathit{T}}}_1,{\overline{\mathit{T}}}_2,\cdots ,{\overline{\mathit{T}}}_K].$$

  (16)
  It is remarkable that matrices (15) and (16) are composed by sub-matrices of dimension $N_x\times N$ and $N\times N_r$, respectively, thus the most efficient way to organize this information is in tensor form as $\mathit{X}(x_i,s_n,t_k)={\mathit{X}}_{ink}$ and $\mathit{Y}(s_n,r_j,t_k)={\mathit{Y}}_{njk}$, being $s_n$ the variable defining the number of singular values, for $n=1,\cdots ,N$. For each velocity component l the previous expressions lead as ${\mathit{X}}_{link}$ and ${\mathit{Y}}_{lnjk}$.
  The reduction of memory size in the previous snapshot matrices is dependent on the number of singular values N retained in each one of such matrices. To ensure the proper performance of the model proposed, N can be determined according to a certain tolerance, as in Equation (7), or simply by identifying changes in the tendency of the singular values obtained in the SVD analysis carried out for each snapshot.
• Stage 2: modal decomposition using HODMD. The multi-resolution HODMD algorithm is applied to each one of the tensors obtained in the previous step, ${\mathit{X}}_{link}$ and ${\mathit{Y}}_{lnjk}$. The following DMD expansion is obtained for the X-snapshots (for simplicity particularized of one velocity component)

  $${\mathit{X}}_{ink}={\displaystyle \sum _{m=1}^{M1}}{a}_m^x{\mathit{u}}_m^x(x_i,s_n)e^{({\delta }_m+i{\omega }_m)t_k},$$

  (17)

  and for the Y-snapshots (for one velocity component)

  $${\mathit{Y}}_{njk}={\displaystyle \sum _{m=1}^{M2}}{a}_m^y{\mathit{u}}_m^y(s_n,r_j)e^{({\delta }_m+i{\omega }_m)t_k},$$

  (18)

  where $M1$ and $M2$ are the DMD modes retained for each expansion, $a_m^x$ and $a_m^y$ are the amplitudes and ${\mathit{u}}_m^x$ and ${\mathit{u}}_m^y$ are the DMD modes. Please note that the number of modes composing each expansion can be different, meaning that the amplitude of each mode describing the dynamics of the system varies in each analysis.
• Stage 3: temporal predictions and data reconstruction. The attractor is predicted by simply adjusting the time as $t_{\xi }\gg t_k$ in (17) and (18). Then, the new data predicted by the previous models are reconstructed as (for each velocity component)

  $$\mathit{V}(x,y,t_{\xi })={\mathit{V}}_{\xi }={\displaystyle \sum _{i=1}^{N_x}}{\displaystyle \sum _{n=1}^N}{\displaystyle \sum _{j=1}^{N_r}}{\mathit{X}}_{in\xi }{\mathit{Y}}_{nj\xi },$$

  (19)

  in matrix form defined as

  $$\mathit{V}(x,y,t_{\xi })={\mathit{X}}_{\xi }{\mathit{Y}}_{\xi },$$

  (20)

  where the dimension of ${\mathit{X}}_{\xi }$ is $N_x\times N$ and the dimension of ${\mathit{Y}}_{\xi }$ is $N\times N_r$.
  As previously mentioned, the number of DMD modes defining each spatial matrix can be different, since the amplitude of the most relevant flow dynamics varies in the spatial domain. Thus, for some groups of tolerances ${\epsilon }_1$, it is possible obtaining $M1$ DMD modes defining expansion (17) and $M2=0$ modes in expansion (18). Similarly, for some tolerances it is possible that $M1=0$ and $M2>0$. Depending on the complexity of the data analyzed, and the influence of the global dynamics related to each spatial direction, a good alternative to reconstruct the solution predicted in time is using the following expansion

  $$\mathit{V}(x,y,t_{\xi })={\mathit{V}}_{\xi }={\displaystyle \sum _{i=1}^{N_x}}{\displaystyle \sum _{n=1}^N}{\displaystyle \sum _{j=1}^{N_r}}{\mathit{X}}_{ink}{\mathit{Y}}_{nj\xi },$$

  (21)

  where ${\mathit{X}}_{ink}$ is the tensor obtained in Stage 1 of this algorithm. Similarly, the predicted solution could be reconstructed as

  $$\mathit{V}(x,y,t_{\xi })={\mathit{V}}_{\xi }={\displaystyle \sum _{i=1}^{N_x}}{\displaystyle \sum _{n=1}^N}{\displaystyle \sum _{j=1}^{N_r}}{\mathit{X}}_{in\xi }{\mathit{Y}}_{njk},$$

  (22)

  where ${\mathit{Y}}_{njk}$ is the tensor obtained in Stage 1.

5. Generating Artificially the Optimal Vortex: A Way to Predict the Perfect Set-Up in Synthetic Jets

This section presents a new theoretical and simple approach suitable to predict the artificial generation of the optimal vortex using synthetic jets. Gharib et al. [6] identified the formation number to generate the optimal vortex as $L_b/D\simeq 4$ for values of Re $\simeq {10}^4$. These authors also relate such number with the running mean of the piston velocity $U_p$ and the piston oscillation period $T=\frac{1}{f}$ as

$$\frac{L_b}{D}=\frac{U_{p}T}{D}\simeq 4.$$

(23)

Thus, it is possible to relate $U_p$ with the oscillation frequency and the formation number as $U_p=4fD$. The oscillatory movement of the piston can be approximated by a simple sinusoidal function, encouraging the flow to continuously leave and re-enter into the cavity, as presented in Section 2. Since this function is even (zero net mass flow), the proper way to define the piston oscillation velocity is considering half of the oscillation period ($T/2=2f$), as

$${\overline{U}}_p=8fD.$$

(24)

It is remarkable that these authors identified the piston velocity as the velocity scale defining the flow dynamics, nevertheless the present article considers the mean flow velocity leaving the jet nozzle as the characteristic velocity scale of the flow (see more details in Refs. [14,15,30]).

Gharib et al. [6] introduced the definition of the running mean of the piston velocity as

$$U_p=\frac{1}{T}{\int }_0^T{v}_p\left(t\right)dt.$$

(25)

Following the work presented in [14,15,30], as introduced in Section 2, the piston velocity varying with time is defined as $v_p\left(t\right)={\widehat{U}}_{p}sin\left(2\pi ft\right)$. Introducing this equation into (25), leads to the following integral, defining the piston velocity in half period of the oscillatory movement as

$${\overline{U}}_p=\frac{1}{T/2}{\int }_0^{T/2}{\widehat{U}}_{p}sin\left(2\pi ft\right)dt=\frac{4}{T}\frac{T}{2\pi }{\left[-cos\left(\frac{2\pi }{T}t\right)\right]}_0^{T/2}=\frac{2{\widehat{U}}_p}{\pi }.$$

(26)

Comparing terms in Equations (24) and (26), it is possible to obtain the following expression for the frequency

$$f=\frac{{\overline{U}}_p}{8D}=\frac{2{\widehat{U}}_p/\pi }{8D}=\frac{{\widehat{U}}_p}{4\pi D}.$$

(27)

This value defines the optimal oscillation frequency of the piston required to generate the optimal vortex, which depends on the jet nozzle diameter and the amplitude of the piston oscillations.

To relate this solution to the bifurcation diagram presented in Ref. [15], it is necessary to calculate the frequency in terms of the non-dimensional Strouhal number, defined in Section 2 as $St=\frac{fD}{U}$. The characteristic velocity scale is defined by the mean flow velocity leaving the jet nozzle, which was related with the piston velocity and dimensions in Equation (3) as $U=\frac{D_{p}D}{\sqrt{2}{\widehat{U}}_p}$. Taking into account the relationship between the diameter of the piston and the jet nozzle, defined as $D_p=5D$, it is possible to approximate the mean flow velocity as

$$U=\frac{D_{p}D}{\sqrt{2}{\widehat{U}}_p}=\frac{5D^2}{\sqrt{2}{\widehat{U}}_p}.$$

(28)

Using Equations (27) and (28) it is possible to approximate the Strouhal number as

$$St\simeq \frac{\sqrt{2}{\widehat{U}}_p^3}{20\pi D^2}.$$

(29)

This expression defines a way to calculate the (non-dimensional) frequency required to generate artificially an optimal vortex using synthetic jets with rounded cavity and nozzle, depending on the amplitude of the piston velocity and the diameter of the jet nozzle. Generally, for non-dimensional values, the jet diameter is set as $D=1$, thus the previous expression only depends on the piston velocity amplitude.

For values of ${\widehat{U}}_p>1$ the frequency is St $\ge 0.03$ (regime(i) in Ref. [15]). On the contrary, if ${\widehat{U}}_p<1$ then the frequency is St $\le 0.02$ (regime (ii) in Ref. [15]). As seen, the efficiency of the system is defined by the amplitude of the piston oscillation, which is in good agreement with the quantity of energy necessary to produce such amplitude, using for instance an engine. The artificial generation of the optimal vortex assumes amplitude values as ${\widehat{U}}_p=1$ (minimizing the quantity of external energy), meaning that the oscillation frequency required to generate such optimal vortex is St ≃ 0.022. This result is in good agreement with the bifurcation diagram presented by Carter & Soria [15] and the experiments carried out by Daribi et al. [16]. On the one hand, Figure 5 in Ref. [15] shows that for values of Re $\le {10}^4$ and St $=0.03$ the flow regime is described by a laminar jet (regime (i)), while for values of St $\le 0.018$ the flow regime is conformed by vortex rings (regime (ii)). From these experimental results it is possible to guess that the transition point between these two regimes is defined by Strouhal values defined in the interval St $\in ]0.018,0.03[$. On the other hand, when Daribi et al. [16] studied in detail the swimming motion of 7 types of jellyfish, the authors showed that the most efficient jellyfish were producing a jet stream defined in regime (ii) (St $\le 0.018$), while the highest impulse was generated by the jellyfish producing a jet stream defined in regime (i) (St $=0.03$). This result connects system efficiency with the energy introduced to the system, by means of controlling the piston oscillation amplitude. As seen, the generation of the optimal vortex lies in the transition between regions (i) and (ii), meaning that this vortex will produce the highest impulse using the minimum amount of energy possible. Increasing the frequency to values larger than $0.022$ will lead to energetic rises coming from the higher amplitude of the piston oscillation required to generate the vortex, and consequently diminishing the system efficiency.

Finally, it is remarkable that the formation number given by Gharib et al. [6] was calculated for values of Re $\simeq {10}^4$ (these authors defined the Reynolds number as function of the vortex circulation $\Gamma$ as $R{e}_{\Gamma }=\Gamma /\nu =2800$, which is equivalent to Re $\simeq {10}^4$ using the jet nozzle diameter and the jet nozzle mean flow velocity as characteristic scales). As explained by Bandyopadhyay & Beal [22], the frequency is not the only parameter that guarantees the presence of the optimal vortex. Thus, based on practical experience, to guarantee the generation of the optimal vortex it is necessary to follow the diagram by Carter & Soria [15] and to fix the Reynolds number within the interval Re $\in ]0,{10}^4]$. The results presented in such diagram show that for values of Re $\le {10}^4$, the transition point is always located in the interval St $\in ]0.018,0.03[$, suggesting that St $=0.022$ is the critical value to generate the optimal vortex. Please note that the equations presented in this section related to the formation number were formulated for Re $={10}^4$; however, the Reynolds number of the results presented in the numerical simulations is Re $={10}^3$, which confirms that these equations are valid for smaller Reynolds numbers, showing that the optimal vortex is generated for values of St $=0.022$ and Re $\le {10}^4$.

The artificial generation of the optimal vortex has been studied in detail comparing numerical results with the previous theoretical approach. Figure 2 compares the topology patterns of three different test cases calculated at Re $=1000$, these are St $=0.02$ (regime (ii)), St $=0.022$ (critical value) and St $=0.03$ (regime (i)). The figure shows contours of instantaneous spanwise vorticity and the streamlines describing the vortex ring in the jet flow at a representative time instant of the injection phase, when the vortex ring leaves the cavity. As seen, at St $=0.02$ it is possible to distinguish two different vortex rings ejected from the cavity, while at St $=0.03$ a single vortex ring followed by a laminar jet with a well-defined structure is clearly identified. The transition between these two regimes is clearly shown at St $=0.022$, where a single vortex ring leaving the cavity is followed by a small vortex in a continuous structure, suggesting the presence of a laminar jet.

6. A ROM to Predict Spatio-Temporal Flow Structures in Synthetic Jets

This section presents a data-driven ROM based on HODMD that predicts spatio-temporal structures in synthetic jets. The main goal of this model is to approximate solutions converged in time from a reduced number of snapshots collected in the transitory region of a numerical simulation, when the solution is still continuously changing due to the presence of many non-permanent modes, decreasing in time. The model has been tested using four different configurations of the database generated from the numerical simulations calculated at St $=0.03$. Due to the presence of flow instabilities in the laminar jet stream combined with the flow instabilities related to the periodic vortex ring ejected from the cavity [30,31], this test case is the most complex of the three test cases simulated. Nevertheless, as it will be shown below, the versatility of this method makes that this tool can be successfully applied to construct ROMs in the data obtained at St $=0.02$ and $0.022$ obtaining similar results (similar speed-up factor and RRMS error, not shown on this paper for the sake of brevity). Two different models will be presented below. In the first case, the standard DMD-based ROM is applied to predict the flow features in synthetic jets, with the aim at showing the good performance of this technique to predict spatio-temporal data in this type of problems. In the second case, the novel approach DMDlr-ROM is applied to predict the previous spatio-temporal flow structures, using only the 30% of the total memory necessary to solve the initial problem with the DMD-based ROM, and reducing the computational time in the generation of the model by a factor of ~3.

6.1. Predictions Using the DMD-Based ROM

In all the cases presented, a group of K snapshots, equidistant in time with $\Delta t=0.053$ are collected starting at time $t_0$, thus they are defined in the time interval $[t_0,t_1]$, with $t_1=K·\Delta t$. The number of piston cycles contained in each snapshot group is given by $c_n=\frac{t_1-t_0}{T}$, with $T=\frac{1}{St}\simeq 33.33$ being the period of each oscillation. The parameters used for each HODMD analysis are $\epsilon ={\epsilon }_1=1.2\times {10}^{-2}$ and the value of d is adjusted proportional to the number of snapshots collected in each case as $d=344·c_n$ (the value of the tolerances and d is based on the results presented in [30,31]). In all cases, the DMD expansion (4) is constructed to create a ROM. The temporal term is then adjusted to time $t_2=333.3=10T$, representing the 10th cycle of oscillation. Then, the solution is predicted 15 more cycles, up to time $t_f=833.25=25T$. The speed-up factor of the numerical code is measured with the ratio $sp=t_f/t1$.

Table 1. Parameters set to generate 4 different DMD-based ROMs without rank reduction. Data collected using K snapshots, in the time interval $[t_0,t_1]$, representing $c_n$ piston cycles. The speed-up factor in the predictions is presented as $sp$. The general solution obtained by DMD-d retains $TM$ modes. The DMD-based ROM is constructed retaining M permanent modes, identified as $\left|\delta \right|<\gamma$, with ${\gamma }_1={10}^{-2}$ and ${\gamma }_2=2\times {10}^{-3}$ in brackets. Each model is constructed enforcing $\delta =0$ and leaving $\delta$ to its value approximated. RRMSE-N, RRMSE-F, and RRMSE represent the RRMSE made in the predictions of the near field, $x/D<5$, the far field $x/D\le 5$, and in the entire flow field, respectively.

Model	K	${\mathit{t}}_0$	${\mathit{t}}_1$	${\mathit{c}}_{\mathit{n}}$	$\mathit{sp}$	d	$\mathit{TM}$	M	$\mathit{\delta }$ = 0	RRMSE-N	RRMSE-F	RRMSE
1	625	$T+\Delta t$	$2T$	1	$12.5$	344	13	0 (0)	Yes	$9.50\times {10}^{-1}$	$9.45\times {10}^{-1}$	$9.50\times {10}^{-1}$
									No	1	1	1
2	1249	$T+\Delta t$	$3T$	2	$8.3$	688	21	7 (0)	Yes	$6.73\times {10}^{-1}$	$8.28\times {10}^{-1}$	$7.82\times {10}^{-1}$
									No	$6.71\times {10}^{-1}$	$6.83\times {10}^{-1}$	$6.79\times {10}^{-1}$
3	1873	$0+\Delta t$	$3T$	3	$8.3$	1032	31	31 (13)	Yes	$2.32\times {10}^{-1}$	$3.88\times {10}^{-1}$	$3.46\times {10}^{-1}$
								31	No	$2.27\times {10}^{-1}$	$3.75\times {10}^{-1}$	$3.35\times {10}^{-1}$
								13	No	$2.31\times {10}^{-1}$	$3.82\times {10}^{-1}$	$3.41\times {10}^{-1}$
4	3121	$0+\Delta t$	$5T$	5	5	1720	57	57 (29)	Yes	$2.21\times {10}^{-1}$	$3.79\times {10}^{-1}$	$3.37\times {10}^{-1}$
								57	No	$4.08\times {10}^{-1}$	$6.29\times {10}^{-1}$	$5.68\times {10}^{-1}$
								29	No	$2.17\times {10}^{-1}$	$3.64\times {10}^{-1}$	$3.24\times {10}^{-1}$

summarizes the analyses carried out and includes the number of DMD modes retained in each case, being $TM$ and M the total number of modes and the number of permanent modes, respectively, which are defined in this example as the modes with growth rate $\left|\delta \right|<\gamma$, with ${\gamma }_1={10}^{-2}$ and ${\gamma }_2=2\times {10}^{-3}$. The RRMSE made in the prediction of the 15 time intervals in the near and far fields and in the complete field are represented by RRMSE-N, RRMSE-F, and RRMSE, respectively.

As seen, the number of permanent modes retained in the DMD expansion (4) increases with the number of snapshots used in the analysis. This explains that the model predicts the real solution with smaller errors, although it is necessary to assume a penalty in the speed-up factor of the code. Moreover, increasing the number of snapshots in the original data analyzed, also increases the complexity in the DMD analysis meaning that for large three-dimensional databases representing complex flows (i.e., transitional or turbulent flows), it could be necessary to combine HODMD with the low-rank strategy to reduce the spatial dimension of the data, as it will be presented in the next section. The error made in the predictions of the data is larger than 60% in Models 1 and 2 and ~34% in Models 3 and 4, with the exception of Model 4 using predictions with ${\gamma }_1$ and $\delta \ne 0$, where the error is larger than 55%. The RRMS error in the near field of Models 3 and 4 is ~23%, meaning that the performance of the model is better predicting the near field than the far field. The speed-up factor in Model 3 is $sp=8.3$, while in Model 4 $sp=5$, suggesting that Model 3 is the best model for temporal predictions in synthetic jets. It is remarkable the large RRMS error made in these calculations, ~34%, nevertheless, the main goal of this DMD-based ROM is to predict the global dynamics representing the main temporal patterns describing the flow and the evolution of the periodic ejection of vortices. As it will be seen below, the model perfectly predicts the general dynamics driving the flow in synthetic jets, represented by the large size flow structures. Increasing the snapshot number, and collecting data in the permanent region of the numerical simulations, will diminish the RRMS error providing information related to the small-size flow scales; however, this will substantially increase the computational cost of the simulations. This work presents a balanced relationship between minimizing the error in the data predicted while minimizing the computational cost of the simulations. Although it is possible to diminish the RRMS error of the data predicted using lager number of snapshots paying a small penalty in the speed-up factor of the code (i.e., Model 4 compared to Model 3), the accuracy in the prediction of the large size flow structures will be kept almost invariant (the error is slightly reduced, some clarifications will be presented below), suggesting that it is better to assume the error with the aim at obtaining large benefits related to the reduction of the computational time.

Figure 3 and Figure 4 compare the original solution at time $t_2=10.05T$ and $t_2=10.25T$ with the predictions carried out in Models 2, 3, and 4. The results for Model 1 are omitted, since they are completely spurious, nevertheless the remaining cases succeed in the prediction of the large flow structures. The smallest flow structures, identified in the original data as small spatial elements (that could be compared with noise), are missing in all the DMD-based models. The jet flow structures in the near and far field are better represented by Model 4 than Models 2 and 3, being the representation of Model 3 better than Model 2. In the results provided by Model 4, it is possible to distinguish (i) the small-size element leaving the jet cavity, related to the vortex ring, and (ii) a group of medium-size structures in the far field, $x/D\ge 5$, possibly related with high-frequency temporal patterns. The presence of these structures is attenuated in Model 3, especially in the case with $\delta =0$, and in Model 2, where the wake of the jet start vanishing at spatial positions $x/D>10$ and the near field of the jet is not properly represented.

Figure 5 compares the frequencies vs. amplitudes and growth rates calculated in Models 1, 2, 3, and 4. The method retains a group of low-amplitude high-frequency spurious modes in Model 1, on the contrary, the remaining models retain a group of periodic modes, with leading frequency St $=0.03$ (forcing frequency), representing the general dynamics of the flow (see [14] for more details). The main differences found in these analyses lies in the (i) calculations of the growth rates, defining the accuracy of the frequency values (the more accurate, the smaller growth rate), showing that the modes are better calculated in Model 4, followed by Model 3 and Model 2, and in (ii) the number of high-frequency modes retained, which represent the flow physics of the far field, related to the small-size periodic flow structures identified in Figure 4. Model 4 retains a larger number of high-frequency modes than Model 3 and Model 2, explaining the better performance of the former model. Nevertheless, the RRMS error in the total interval predicted $[10T,25T]$, is very similar in Models 3 and 4, suggesting that both models perfectly predict the general dynamics in synthetic jets for large time periods.

The temporal evolution of Models 2, 3, and 4 is compared with the original solution in Figure 6 and Figure 7 in two representative points of the computational domain, extracted in the near and far fields. As predicted by the RRMS error (>60%), Model 2 fails in the prediction of the 15 cycles, while models 3 and 4 perfectly predicts both the near and far fields up for 26 cycles (from cycle 10 to cycle 36), with the exception of Model 4 using ${\gamma }_1$ and $\delta \ne 0$, whose solution diverges in time. It is possible to overcome such problem using ${\gamma }_2$ instead, as shown in Figure 8, where Model 4 provides similar solutions as Model 4 with $\delta =0$ and Model 3.

The previous figures suggest that the qualitative error made in the calculations of the flow structures in the near field is smaller than 23%, and in the far field is smaller than ~34%, as predicted in the RRMSE calculations presented in Table 1. Thus, the high value presented in such calculations could be related with the small-size flow structures defining the flow. This fact is confirmed in Figure 9, comparing the real solution with the velocity fields predicted in the near and far fields in cycle 22 (predictions up to cycle 25). The general dynamics, the large size flow scales, are predicted with RRMS error smaller than 3% (calculated with the flow field averaged in time); however, the method is not able to predict some very high-frequency flow oscillations, related with low-amplitude/small-size flow structures or with the medium-size structures describing the wake of the jet. This fact confirms the good performance of Model 3 as a DMD-based ROM to predict large size flow structures in synthetic jets with speed-up factor of $8.3$ and errors smaller than 3% in the streamwise velocity. Please note that the smaller size structures could be predicted with smaller error (i) using larger number of snapshots, with smaller temporal distant $\Delta t$ to construct the snapshot matrix (5) and (ii) using information of the numerical simulations converged in time, but these two steps will increase the computational cost of the ROM and will not provide additional information of the general flow dynamics, which is the main scope of this paper, although it remains as open topic for future research. Finally, it is remarkable that similar results are found for the normal velocity component, not shown for the sake of brevity, since such vector field is one order of magnitude smaller than the streamwise velocity field.

6.2. Efficient Computations Using the DMDlr-ROM

The previous section shows that Model 3 provides a good general approximation of the flow field in the temporal predictions carried out, thus this is the model selected to perform the low-rank approximation to construct the DMDlr-ROM. The parameters set for the DMD analysis are listed in Table 1 (snapshot number K, d, $t_0$, and $t_f$ for Model 3). Before starting with the DMD analysis, it is necessary to set the value of N singular values that will be retained to construct the X-snapshot and Y-snapshot tensors (15) and (16). Figure 10 shows the evolution of the singular values calculated in 1873 snapshots. A clear change of tendency is observed in a few snapshots for values of $N\simeq 10$, while this change of tendency is observed in the remaining snapshots for values closer to $N=20$. To ensure the proper performance of this ROM, the number of singular values is set as $N=20$. The size of the snapshot matrix is reduced by ~70%, from $1.6$ Gb to 472 Mb and 330 Mb for the X-snapshot and Y-snapshot tensors, respectively, and the computational time required to construct the model is also reduced by ~70%, from $461.70$ to $82.02$ seconds for each case.

HODMD is applied using the tolerances $\epsilon ={\epsilon }_1=5\times {10}^{-2}$ retaining 19 and 0 modes in the expansions (18) and (17), respectively, suggesting that the X-snapshot tensor (15) is related with smaller amplitude modes. Thus, the same analysis is carried out to such tensor, using the tolerances $\epsilon =5\times {10}^{-3}$ and ${\epsilon }_1={10}^{-2}$, retaining 15 modes. The frequencies vs. amplitudes and growth rates calculated for each analysis are compared in Figure 11.

From the 15 modes calculated in the X-snapshot tensor (low tolerances), 12 of them represent the periodic modes with leading frequency St $=0.03$, one of them represent the steady mode with zero frequency, and the two remaining modes are related with high frequencies, decaying in time with growth rate ~10⁻² (note that the modes include their complex conjugate counterpart). On the contrary, in analysis of the Y-snapshot tensor (high tolerances), the method retains 10 periodic modes, the zero-frequency mode, and a group of high-frequency modes. This fact suggests that the high-frequency spatio-temporal structures found in the wake of the jet, shown in Figure 4 presents higher amplitude in the data represented in the Y-snapshot tensor.

The predictions of the flow structures in the time interval $[10T,25T]$ are carried out using expansions (17) and (18), and the data are reconstructed using the three different models defined by Equations (19), (21) and (22).

Table 2. Parameters set to generate 3 different DMDlr-ROMs with $N=20$, using Model 3 with ${\gamma }_1$ defined in Table 1. RRMSE-N, RRMSE-F, and RRMSE represent the RRMSE made in the predictions of the near field with $x/D<5$, the far field with $x/D\le 5$, and in the entire flow field.

Model 3	Eq.	$\mathit{\delta }=0$	RRMSE-N	RRMSE-F	RRMSE
A	(21)	Yes	$2.48\times {10}^{-1}$	$2.54\times {10}^{-1}$	$2.52\times {10}^{-1}$
		No	$2.92\times {10}^{-1}$	$3.02\times {10}^{-1}$	$2.99\times {10}^{-1}$
B	(22)	Yes	$2.62\times {10}^{-1}$	$3.94\times {10}^{-1}$	$3.58\times {10}^{-1}$
		No	$3.25\times {10}^{-1}$	$4.47\times {10}^{-1}$	$4.12\times {10}^{-1}$
C	(19)	Yes	$3.64\times {10}^{-1}$	$4.47\times {10}^{-1}$	$4.43\times {10}^{-1}$
		No	$4.44\times {10}^{-1}$	$5.24\times {10}^{-1}$	$5.01\times {10}^{-1}$

presents the RRMS error made in the predictions of the near and far fields and the complete flow field. The error made in the results presented in Model $3-A$ using $\delta =0$ is smaller than ~25% in both the near and far field calculations, improving the results obtained in the previous section, without any rank reduction. This model predicts the spatio-temporal flow structures only considering the temporal variations found in X, related to the low-amplitude modes, mainly represented by the periodic modes with leading frequency St $=0.03$. On the contrary, the largest errors are found in Model $3-C$, which considers the variations in time of the X and Y spatial components. This fact suggests that the complex flow patterns describing the vortex ring and the wake in synthetic jets is related to small amplitude variations in time along the streamwise direction, while the temporal variations along the radial component remain almost constant.

Figure 12 and Figure 13 compare the streamwise velocity contours of the original flow field with the solution obtained with Models $3-A$, $3-B$, and $3-C$. The figures show qualitatively the good performance of the three models to predict the large size flow structures, and especially Model $3-A$, which also predicts the medium-size structures in the wake of the jet.

The temporal evolution of these three models is shown in Figure 14 and Figure 15. In contrast to Model 3 presented in the previous section, the performance of these low-rank approximations is better using $\delta =0$ in the permanent modes. For values of $\delta \ne 0$, all the models diverge as time passes by. Regarding the models with $\delta =0$, the best performance is described by Model $3-A$, whose predictions could be extended up to time 1200 (36 cycles). The performance of Model $3-B$ is also good, although this model starts slowly diverging after cycle 25. Finally, Model $3-C$ slowly diverges after cycle 18. The worsen performance of Models $3-B$ and $3-C$ compared to Model $3-A$ could be related to the high amplitude high-frequency modes retained by these models associated with the radial spatial component (Y-snapshots). On the contrary, Model $3-A$ only considers the temporal evolution of the small amplitude periodic modes related to the streamwise direction. The high-frequency modes describing the medium-size flow structures forming the wake of the jet seems to be related with a group of very low-amplitude modes found in the radial component (velocity fluctuations).

This fact is confirmed in Figure 16 that shows a zoomed-in view of Models $3-A$ and $3-B$ in a representative point extracted in the far field in cycle 22. As seen, Model $3-A$ predicts some of the very high-frequency flow oscillations, while Model $3-B$ provides the average of the flow oscillation related to the large size flow structures. Model $3-A$ uses the full model for the Y-snapshot tensor (radial component), which also includes the small amplitude high-frequency modes, related to the medium-size flow structures of the wake of the jet. On the contrary, in Model $3-B$ only a few modes, related to the large size flow structures, approximates the radial component.

To sum up, the DMDlr-ROM predicts the general dynamics in synthetic jets providing information about the flow physics and the main global temporal patterns. On the one hand, the mechanism describing the periodic formation of the characteristic vortex ring, which travels upstream, is driven by small amplitude periodic modes with leading frequency St $=0.03$ (forcing frequency). These modes describe the global dynamics of the flow, composed by large size flow structures. On the other hand, a group of very small amplitude and very high-frequency modes related to the radial component, which are connected with the velocity fluctuations, describe the group of medium-size flow structures found in the wake of the jet, suggesting the connection of these modes with the main mechanism leading to the bifurcation process (from large size to small-size flow scales) carried out in the wake of the jet.

7. Conclusions

This article presents three strategies to predict spatio-temporal structures in synthetic jets. This first method is a simple theoretical approach, connecting experimental and numerical results, which predicts the perfect set-up in synthetic jets with circular cavity to generate artificially the optimal vortex: the vortex that maximizes thrust and system efficiency. This approach is confirmed with results obtained in numerical simulations. The second method is a data-driven ROM that uses HODMD to predict the temporal evolution of the main flow structures describing the flow, with the aim at reducing the computational time in numerical simulations. This method predicts the large size flow structures with a speed-up factor of $8.3$ in the simulations. Finally, the third strategy is an extension of the second method that uses a low-rank approximation of the spatial information, reducing the size of the data analyzed by 70% and reducing the computational time to create the model by 3. This novel approach perfectly predicts the large size flow structures maintaining the speed-up factor to $8.3$, and provides a deepen insight into the flow physics. The model shows that the periodic generation and convection of the vortex ring is driven by a group of low-amplitude periodic modes in the streamwise direction, with leading frequency the forcing frequency of the jet. Moreover, the medium-size flow structures describing the wake of the jet are connected with small amplitude and very high-frequency modes in the radial spatial component.