*2.3. Dynamic Mode Decomposition*

DMD is an equation-free, data-driven method for data analysis and behavior prediction of complex dynamical systems. It was first proposed by Schmid [34] to analyze the high-dimensional fluid dynamics data by decomposing it into coherent spatial structures that oscillate at distinct frequencies. Thanks to its ability both to analyze and to predict, it has gained successes not only in fluid mechanics, but also in other fields, including video processing and finance, where high-dimensional complex dynamic systems are involved [35].

For fluid mechanics applications, the input of DMD analysis is a sequence of snapshots of the flow field. The snapshot **x***<sup>i</sup>* is a column vector of dimension *n* containing all the interested variables at measure points in the flow field at time *t* = *ti*. The snapshots are taken at a fixed time interval Δ*t*. With total *m* snapshots, the snapshot matrix of the dataset can be written as:

$$\mathbf{X}\_1^{\mathrm{pr}} = \begin{bmatrix} | & | & & | \\ \mathbf{x}\_1 & \mathbf{x}\_2 & \dots & \mathbf{x}\_m \\ | & | & & | \end{bmatrix}\_{n \times m} \tag{9}$$

where the sub and sup indexes indicate the starting and the ending index of snapshots. The DMD approximates the dynamic system by a time-independent linear matrix **A** representing the temporal evolution from one snapshot to the next, as follows,

$$\mathbf{X}\_2^{\mathfrak{m}} = \mathbf{A}\_{\mathfrak{m}\times\mathfrak{n}} \mathbf{X}\_1^{\mathfrak{m}-1}.\tag{10}$$

In practice, the spatial dimension *n* is usually very large, which makes the direct solution of **A***n*×*<sup>n</sup>* difficult. Instead, the reduced-order Proper Orthogonal Decomposition (POD) projection matrix **A˜** *<sup>r</sup>*×*<sup>r</sup>* with *<sup>r</sup> <sup>n</sup>* is used. *<sup>r</sup>* is the number of the most energetic modes kept in the singular value decomposition of **<sup>X</sup>***m*−<sup>1</sup> <sup>1</sup> , defined as follows,

$$\mathbf{X}\_1^{m-1} \approx \mathbf{U}\_{n \times r} \boldsymbol{\Sigma}\_{r \times r} \mathbf{V}\_{r \times (m-1)}^T \tag{11}$$

where the left singular vectors in **U** are POD modes, the right singular vectors in **V** are time coefficient of these modes, and **Σ** is a diagonal matrix containing the first *r* largest singular values. **U** and **V** are orthonormal. **A˜** is defined as

$$\tilde{\mathbf{A}} = \mathbf{U}^T \mathbf{A} \mathbf{U} = \mathbf{U}^T \mathbf{X}\_2^m \mathbf{V} \mathbf{Z}^{-1}. \tag{12}$$

The next step is to calculate the eigenvalues and eigenvectors of **A˜** ,

$$
\mathbf{\tilde{A}W} = \mathbf{W}\Lambda\_\prime \tag{13}
$$

where **Λ** is a diagonal matrix containing *r* eigenvalues and **W** contains *r* eigenvectors of unit length. The snapshots **<sup>x</sup>***<sup>k</sup>* <sup>=</sup> **Ax***k*−<sup>1</sup> <sup>=</sup> **<sup>A</sup>***k*−1**x**<sup>1</sup> can be written as

$$\mathbf{x}\_k = \mathbf{U}\mathbf{W}\boldsymbol{\Lambda}^{k-1}\mathbf{W}^{-1}\mathbf{U}^T\mathbf{x}\_1 = \boldsymbol{\Phi}\boldsymbol{\Lambda}^{k-1}\mathbf{b},\tag{14}$$

where **Φ** = **UW** contains the DMD modes (*φi*) and **b** = **W**−1**U***T***x**<sup>1</sup> contains the amplitudes (*bi*) of these modes in the first snapshot. Each mode *φ<sup>i</sup>* has a corresponding eigenvalue *λ<sup>i</sup>* in **Λ**. The oscillating frequency (*fi*) of mode *φ<sup>i</sup>* is equal to

$$f\_i = |\frac{\heartsuit (\log(\lambda\_i))}{2\pi\Delta t}|,\tag{15}$$

and the growth rate *gi* is equal to

$$\mathbf{g}\_i = \frac{\Re(\log(\lambda\_i))}{\Delta t}.\tag{16}$$

For a stable oscillating system, all the growing rate should be equal to 0. A positive growing rate (*gi* > 0) indicates that the mode *φ<sup>i</sup>* has an amplitude *bi* increasing with time; a negative growing rate (*gi* < 0), in contrast, indicates a mode damping with time.

With the amplitudes *bi*, the most energetic DMD modes in the first snapshot can be easily identified and remain the same for other snapshots when the system is stable. In contrast, when growing and damping modes are concerned, the modes should be ordered by the time-averaged amplitudes *b*- *i* . To calculate the time-averaged amplitudes, the eigenvalue-weighted method proposed by Kou and Zhang [36] is selected for its simplicity and computational efficiency among others [37]. The eigenvalue-weighted amplitudes *b*- *<sup>i</sup>* are calculated as,

$$b\_i' = \frac{\sum\_{j=1}^{m} |b\_i \lambda\_i^{j-1}|}{m}.\tag{17}$$
