OptDMD

Derived module from pydmd.dmdbase() for the optimal closed-form solution to dmd.

Note

P. Heas & C. Herzet. Low-rank dynamic mode decomposition: optimal solution in polynomial time. arXiv:1610.02962. 2016.

class OptDMD(factorization='evd', svd_rank=0, tlsq_rank=0, opt=False)[source]

Bases: DMDBase

Dynamic Mode Decomposition

This class implements the closed-form solution to the DMD minimization problem. It relies on the optimal solution given by [HeasHerzet16].

[HeasHerzet16]

P. Heas & C. Herzet. Low-rank dynamic mode decomposition: optimal solution in polynomial time. arXiv:1610.02962. 2016.

Parameters:
  • factorization (str) – compute either the eigenvalue decomposition of the unknown high-dimensional DMD operator (factorization=”evd”) or its singular value decomposition (factorization=”svd”). Default is “evd”.

  • svd_rank (int or float) – the rank for the truncation; If 0, the method computes the optimal rank and uses it for truncation; if positive interger, the method uses the argument for the truncation; if float between 0 and 1, the rank is the number of the biggest singular values that are needed to reach the ‘energy’ specified by svd_rank; if -1, the method does not compute truncation.

  • tlsq_rank (int) – rank truncation computing Total Least Square. Default is 0, that means TLSQ is not applied.

  • opt (bool or int) – argument to control the computation of DMD modes amplitudes. See DMDBase. Default is False.

_compute_amplitudes(modes, snapshots, eigs, opt)[source]

Compute the amplitude coefficients. If self._opt is False the amplitudes are computed by minimizing the error between the modes and the first snapshot; if self._opt is True the amplitudes are computed by minimizing the error between the modes and all the snapshots, at the expense of bigger computational cost.

This method uses the class variables self.snapshots (for the snapshots), self.modes and self.eigs.

Returns:

the amplitudes array

Return type:

numpy.ndarray

References for optimal amplitudes: Jovanovic et al. 2014, Sparsity-promoting dynamic mode decomposition, https://hal-polytechnique.archives-ouvertes.fr/hal-00995141/document

property amplitudes

Get the coefficients that minimize the error between the original system and the reconstructed one. For futher information, see dmdbase._compute_amplitudes.

Returns:

the array that contains the amplitudes coefficient.

Return type:

numpy.ndarray

property dynamics

Get the time evolution of each mode.

\mathbf{x}(t) \approx \sum_{k=1}^{r} \boldsymbol{\phi}_{k} \exp \left( \omega_{k} t \right) b_{k} = \sum_{k=1}^{r} \boldsymbol{\phi}_{k} \left( \lambda_{k} \right)^{\left( t / \Delta t \right)} b_{k}

Returns:

the matrix that contains all the time evolution, stored by row.

Return type:

numpy.ndarray

property eigs

Get the eigenvalues of A tilde.

Returns:

the eigenvalues from the eigendecomposition of atilde.

Return type:

numpy.ndarray

property factorization
fit(X, Y=None)[source]

Compute the Dynamic Modes Decomposition to the input data.

Parameters:
  • X (numpy.ndarray or iterable) – the input snapshots.

  • Y (numpy.ndarray or iterable) – the input snapshots at sequential timestep, if passed. Default is None.

property fitted

Check whether this DMD instance has been fitted.

Returns:

True is the instance has been fitted, False otherwise.

Return type:

bool

property modes

Get the matrix containing the DMD modes, stored by column.

Returns:

the matrix containing the DMD modes.

Return type:

numpy.ndarray

property modes_activation_bitmask

Get the bitmask which controls which DMD modes are enabled at the moment in this DMD instance.

The DMD instance must be fitted before this property becomes valid. After fit() is called, the defalt value of modes_activation_bitmask is an array of True values of the same shape of amplitudes().

The array returned is read-only (this allow us to react appropriately to changes in the bitmask). In order to modify the bitmask you need to set the field to a brand-new value (see example below).

Example:

>>> # this is an error
>>> dmd.modes_activation_bitmask[[1,2]] = False
ValueError: assignment destination is read-only
>>> tmp = np.array(dmd.modes_activation_bitmask)
>>> tmp[[1,2]] = False
>>> dmd.modes_activation_bitmask = tmp
Returns:

The DMD modes activation bitmask.

Return type:

numpy.ndarray

predict(X)[source]

Predict the output Y given the input X using the fitted DMD model.

Parameters:

X (numpy.ndarray) – the input vector.

Returns:

one time-step ahead predicted output.

Return type:

numpy.ndarray