Dynamic mode decomposition

In fluids applications, the size of a snapshot, ${\displaystyle M}$, is assumed to be much larger than the number of snapshots ${\displaystyle N}$, so there are many equally valid choices of ${\displaystyle A}$. The original DMD algorithm picks ${\displaystyle A}$ so that each of the snapshots in ${\displaystyle V_{2}^{N}}$ can be expressed as linear combinations of the snapshots in ${\displaystyle V_{1}^{N-1}}$. Because most of the snapshots appear in both data sets, this representation is error free for all snapshots except ${\displaystyle v_{N}}$, which is written as

${\displaystyle v_{N}=a_{1}v_{1}+a_{2}v_{2}+\dots +a_{N-1}v_{N-1}+r=V_{1}^{N-1}a+r,}$

where ${\displaystyle a={a_{1},a_{2},\dots ,a_{N-1}}}$ is a set of coefficients DMD must identify and ${\displaystyle r}$ is the residual. In total,

${\displaystyle V_{2}^{N}=AV_{1}^{N-1}+re_{N-1}^{T}=V_{1}^{N-1}S+re_{N-1}^{T},}$

where ${\displaystyle S}$ is the companion matrix

${\displaystyle S={\begin{pmatrix}0&0&\dots &0&a_{1}\\1&0&\dots &0&a_{2}\\0&1&\dots &0&a_{3}\\\vdots &\vdots &\ddots &\vdots &\vdots \\0&0&\dots &1&a_{N-1}\end{pmatrix}}.}$

The vector ${\displaystyle a}$ can be computed by solving a least squares problem, which minimizes the overall residual. In particular if we take the QR decomposition of ${\displaystyle V_{1}^{N-1}=QR}$, then ${\displaystyle a=R^{-1}Q^{T}v_{N}}$.

In this form, DMD is a type of Arnoldi method, and therefore the eigenvalues of ${\displaystyle S}$ are approximations of the eigenvalues of ${\displaystyle A}$. Furthermore, if ${\displaystyle y}$ is an eigenvector of ${\displaystyle S}$, then ${\displaystyle V_{1}^{N-1}y}$ is an approximate eigenvector of ${\displaystyle A}$. The reason an eigendecomposition is performed on ${\displaystyle S}$ rather than ${\displaystyle A}$ is because ${\displaystyle S}$ is much smaller than ${\displaystyle A}$, so the computational cost of DMD is determined by the number of snapshots rather than the size of a snapshot.

Instead of computing the companion matrix ${\displaystyle S}$, the SVD-based approach yields the matrix ${\displaystyle {\tilde {S}}}$ that is related to ${\displaystyle A}$ via a similarity transform. To do this, assume we have the SVD of ${\displaystyle V_{1}^{N-1}=U\Sigma W^{T}}$. Then

${\displaystyle V_{2}^{N}=AV_{1}^{N-1}+re_{N-1}^{T}=AU\Sigma W^{T}+re_{N-1}^{T}.}$

Equivalent to the assumption made by the Arnoldi-based approach, we choose ${\displaystyle A}$ such that the snapshots in ${\displaystyle V_{2}^{N}}$ can be written as the linear superposition of the columns in ${\displaystyle U}$, which is equivalent to requiring that they can be written as the superposition of POD modes. With this restriction, minimizing the residual requires that it is orthogonal to the POD basis (i.e., ${\displaystyle U^{T}r=0}$). Then multiplying both sides of the equation above by ${\displaystyle U^{T}}$ yields ${\displaystyle U^{T}V_{2}^{N}=U^{T}AU\Sigma W^{T}}$, which can be manipulated to obtain

${\displaystyle U^{T}AU=U^{T}V_{2}^{N}W\Sigma ^{-1}\equiv {\tilde {S}}.}$

Because ${\displaystyle A}$ and ${\displaystyle {\tilde {S}}}$ are related via similarity transform, the eigenvalues of ${\displaystyle S}$ are the eigenvalues of ${\displaystyle A}$, and if ${\displaystyle y}$ is an eigenvector of ${\displaystyle {\tilde {S}}}$, then ${\displaystyle Uy}$ is an eigenvector of ${\displaystyle A}$.

In summary, the SVD-based approach is as follows:

1. Split the time series of data in ${\displaystyle V_{1}^{N}}$ into the two matrices ${\displaystyle V_{1}^{N-1}}$ and ${\displaystyle V_{2}^{N}}$.
2. Compute the SVD of ${\displaystyle V_{1}^{N-1}=U\Sigma W^{T}}$.
3. Form the matrix ${\displaystyle {\tilde {S}}=U^{T}V_{2}^{N}W\Sigma ^{-1}}$, and compute its eigenvalues ${\displaystyle \lambda _{i}}$ and eigenvectors ${\displaystyle y_{i}}$.
4. The ${\displaystyle i}$-th DMD eigenvalue is ${\displaystyle \lambda _{i}}$ and ${\displaystyle i}$-th DMD mode is the ${\displaystyle Uy_{i}}$.

The advantage of the SVD-based approach over the Arnoldi-like approach is that noise in the data and numerical truncation issues can be compensated for by truncating the SVD of ${\displaystyle V_{1}^{N-1}}$. As noted in [3] accurately computing more than the first couple modes and eigenvalues can be difficult on experimental data sets without this truncation step.

Fig 1 Trailing edge Vortices (Entropy)

The wake of an obstacle in the flow may develop a Kármán vortex street. The Fig.1 shows the shedding of a vortex behind the trailing edge of a profile. The DMD-analysis was applied to 90 sequential Entropy fields (animated gif (1.9MB)) and yield an approximated eigenvalue-spectrum as depicted below. The analysis was applied to the numerical results, without referring to the governing equations. The profile is seen in white. The white arcs are the processor boundaries since the computation was performed on a parallel computer using different computational blocks.

Fig.2 DMD-spectrum

Roughly a third of the spectrum was highly damped (large, negative ${\displaystyle \lambda _{r}}$) and is not shown. The dominant shedding mode is shown in the following pictures. The image to the left is the real part, the image to the right, the imaginary part of the eigenvector.

Again, the entropy-eigenvector is shown in this picture. The acoustic contents of the same mode is seen in the bottom half of the next plot. The top half corresponds to the entropy mode as above.

The DMD analysis assumes a pattern of the form ${\displaystyle q(x_{1},x_{2},x_{3},\ldots )=e^{cx_{1}}{\hat {q}}(x_{2},x_{3},\ldots )}$ where ${\displaystyle x_{1}}$ is any of the independent variables of the problem, but has to be selected in advance. Take for example the pattern

${\displaystyle q(x,y,t)=e^{-i\omega t}{\hat {q}}(x,t)e^{-(y/b)^{2}}\Re \left\{e^{i(kx-\omega t)}\right\}+{\text{random noise}}}$

With the time as the preselected exponential factor.

A sample is given in the following figure with ${\displaystyle \omega =2\pi /0.1}$, ${\displaystyle b=0.02}$ and ${\displaystyle k=2\pi /b}$. The left picture shows the pattern without, the right with noise added. The amplitude of the random noise is the same as that of the pattern.

A DMD analysis is performed with 21 synthetically generated fields using a time interval ${\displaystyle \Delta t=1/90{\text{ s}}}$, limiting the analysis to ${\displaystyle f=45{\text{ Hz}}}$.

The spectrum is symmetric and shows three almost undamped modes (small negative real part), whereas the other modes are heavily damped. Their numerical values are ${\displaystyle \omega _{1}=-0.201,\omega _{2/3}=-0.223\pm i62.768}$ respectively. The real one corresponds to the mean of the field, whereas ${\displaystyle \omega _{2/3}}$ corresponds to the imposed pattern with ${\displaystyle f=10{\text{ Hz}}}$. Yielding a relative error of −1/1000. Increasing the noise to 10 times the signal value yields about the same error. The real and imaginary part of one of the latter two eigenmodes is depicted in the following figure.