DMDBase¶
Base module for the DMD: fit method must be implemented in inherited classes
- class DMDBase(svd_rank=0, tlsq_rank=0, exact=False, opt=False, rescale_mode=None, forward_backward=False, sorted_eigs=False, tikhonov_regularization=None)[source]
Bases:
object
Dynamic Mode Decomposition base class.
- Parameters:
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 no truncation.
exact (bool) – flag to compute either exact DMD or projected DMD. Default is False.
opt (bool or int) – If True, amplitudes are computed like in optimized DMD (see
_compute_amplitudes()
for reference). If False, amplitudes are computed following the standard algorithm. If opt is an integer, it is used as the (temporal) index of the snapshot used to compute DMD modes amplitudes (following the standard algorithm). The reconstruction will generally be better in time instants near the chosen snapshot; however increasing opt may lead to wrong results when the system presents small eigenvalues. For this reason a manual selection of the number of eigenvalues considered for the analyisis may be needed (check svd_rank). Also setting svd_rank to a value between 0 and 1 may give better results. Default is False.rescale_mode ({'auto'} or None or numpy.ndarray) – Scale Atilde as shown in 10.1016/j.jneumeth.2015.10.010 (section 2.4) before computing its eigendecomposition. None means no rescaling, ‘auto’ means automatic rescaling using singular values, otherwise the scaling factors.
forward_backward (bool) – If True, the low-rank operator is computed like in fbDMD (reference: https://arxiv.org/abs/1507.02264). Default is False.
sorted_eigs ({'real', 'abs'} or False) – Sort eigenvalues (and modes/dynamics accordingly) by magnitude if sorted_eigs=’abs’, by real part (and then by imaginary part to break ties) if sorted_eigs=’real’. Default: False.
tikhonov_regularization (int or float) – Tikhonov parameter for the regularization. If None, no regularization is applied, if float, it is used as the \lambda tikhonov parameter.
- Variables:
original_time (dict) –
dictionary that contains information about the time window where the system is sampled:
t0 is the time of the first input snapshot;
tend is the time of the last input snapshot;
dt is the delta time between the snapshots.
dmd_time (dict) –
dictionary that contains information about the time window where the system is reconstructed:
t0 is the time of the first approximated solution;
tend is the time of the last approximated solution;
dt is the delta time between the approximated solutions.
- _allocate_modes_bitmask_proxy()[source]
Utility method which allocates the activation bitmask proxy using the quantities that are currently available in this DMD instance. Fails quietly if the amplitudes are not set.
- _compare_data_shapes(y_snapshots)[source]
Method that ensures that the data inputs X and Y are the same shape if provided separately. Throws an error if the shapes do not agree.
- _compute_amplitudes()[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:
References for optimal amplitudes: Jovanovic et al. 2014, Sparsity-promoting dynamic mode decomposition, https://hal-polytechnique.archives-ouvertes.fr/hal-00995141/document
- _optimal_dmd_matrices()[source]
- _reset()[source]
Reset this instance. Should be called in
fit()
.
- _set_initial_time_dictionary(time_dict)[source]
Set the initial values for the class fields time_dict and original_time. This is usually called in fit() and never again.
- Parameters:
time_dict (dict) – Initial time dictionary for this DMD instance.
- _translate_eigs_exponent(tpow)[source]
Transforms the exponent of the eigenvalues in the dynamics formula according to the selected value of self._opt (check the documentation for opt in
__init__
).
- 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:
- property dmd_time
A dictionary which contains information about the time window used to reconstruct/predict using this DMD instance. By default this is equal to
original_time()
.Inside the dictionary:
Key
Value
t0
Time of the first output snapshot.
tend
Time of the last output snapshot.
dt
Timestep between two snapshots.
- Returns:
A dict which contains info about the input time frame.
- Return type:
- property dmd_timesteps
Get the timesteps of the reconstructed states.
- Returns:
the time intervals of the original snapshots.
- Return type:
- 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:
- property eigs
Get the eigenvalues of A tilde.
- Returns:
the eigenvalues from the eigendecomposition of atilde.
- Return type:
- fit(X)[source]
Abstract method to fit the snapshots matrices.
Not implemented, it has to be implemented in subclasses.
- property fitted
Check whether this DMD instance has been fitted.
- Returns:
True is the instance has been fitted, False otherwise.
- Return type:
- property frequency
Get the amplitude spectrum.
- Returns:
the array that contains the frequencies of the eigenvalues.
- Return type:
- property growth_rate
Get the growth rate values relative to the modes.
- Returns:
the Floquet values
- Return type:
- static load(fname)[source]
Load the object from fname using the pickle module.
- Returns:
The ReducedOrderModel loaded
Example:
>>> from pydmd import DMD >>> dmd = DMD.load('pydmd.dmd') >>> print(dmd.reconstructed_data)
- property modes
Get the matrix containing the DMD modes, stored by column.
- Returns:
the matrix containing the DMD modes.
- Return type:
- 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 ofamplitudes()
.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:
- property operator
Get the instance of DMDOperator.
- Returns:
the instance of DMDOperator
- Return type:
DMDOperator
- property original_time
A dictionary which contains information about the time window used to fit this DMD instance.
Inside the dictionary:
Key
Value
t0
Time of the first input snapshot (0 by default).
tend
Time of the last input snapshot (usually corresponds to the number of snapshots).
dt
Timestep between two snapshots (1 by default).
- Returns:
A dict which contains info about the input time frame.
- Return type:
- property original_timesteps
Get the timesteps of the original snapshot.
- Returns:
the time intervals of the original snapshots.
- Return type:
- property reconstructed_data
Get the reconstructed data.
- Returns:
the matrix that contains the reconstructed snapshots.
- Return type:
- save(fname)[source]
Save the object to fname using the pickle module.
- Parameters:
fname (str) – the name of file where the reduced order model will be saved.
Example:
>>> from pydmd import DMD >>> dmd = DMD(...) # Construct here the rom >>> dmd.fit(...) >>> dmd.save('pydmd.dmd')
- property snapshots
Get the input data (space flattened).
- Returns:
the matrix that contains the flattened snapshots.
- Return type:
- property snapshots_shape
Get the original input snapshot shape.
- Returns:
input snapshots shape.
- Return type:
tuple
- property snapshots_y
Get the input left-hand side data (space flattened) if given.
- Returns:
matrix that contains the flattened left-hand side snapshots.
- Return type: