DMD en Python 1D

introduction

La décomposition en mode dynamique (DMD) a été annoncée en 2008 [^ 1] [^ 2]. La DMD est une méthode de décomposition qui combine l'analyse des composantes principales et la transformée de Fourier, et est une méthode qui peut séparer les données en entités dans les directions dimensionnelle et temporelle.

Les articles de commentaire auxquels j'ai fait référence sont ci-dessous.

• http://www.pyrunner.com/weblog/2016/07/25/dmd-python/
• <https://iqujack-lequina.hatenablog.com/entry/2018/05/20/ Mémorandum sur la décomposition en mode dynamique>
• https://qiita.com/Miyabi1456/items/702f62c9bcd9c063d664
• https://qiita.com/harmegiddo/items/84552c32f4b75c26878a

Les valeurs propres et les vecteurs propres de DMD sont des nombres complexes, et les composantes complexes des valeurs propres représentent la décomposition dans la direction du temps. Cependant, il existe un problème en ce qu'il ne peut pas être bien séparé s'il s'agit d'une onde stationnaire unidimensionnelle. L'explication est dans l'article suivant.

• <https://iqujack-lequina.hatenablog.com/entry/2018/12/02/ DMD onde stationnaire>

Je ne pense pas que DMD sera utilisé pour les données unidimensionnelles, mais comme le code Python n'a pas été publié, je le décrirai pour référence. Voici quelques exemples qui ne peuvent pas être bien décomposés et d'autres qui le peuvent. Le document de base est Tu 2014 [^ 3], et le code est DMD Book. Il est basé sur le code Matlab de.

environnement

• Windows 10 home

• Anaconda(Python 3.7.6)

• Numpy(1.18.1)

Objectif de démontage

Trouvez la période à partir des données numériques de l'onde sinusoïdale unidimensionnelle. La période est de 1j $ pour les nombres complexes. $ y = \sin(x) $

Exemple DMD qui ne fonctionne pas

Dans l'exemple qui n'a pas abouti, puisque les données ne sont qu'unidimensionnelles, la valeur propre devient un nombre réel et aucun nombre complexe n'apparaît. Par conséquent, la composante périodique ne peut pas être décomposée.

import numpy as np
import numpy.linalg as LA
import matplotlib.pyplot as plt
%matplotlib inline

t = np.arange(0, 10.01, 0.01)
x = np.sin(t)
dt = t[2] - t[1]

#Matrice de données
X1 = x[:-1]
X2 = x[1:]

#SVD
U,S,Vh = LA.svd(X1[np.newaxis,:],False)
V = Vh.T

#Atilde avec A réduit en dimension par le vecteur singulier gauche
Atilde = np.dot(np.dot(np.dot(U.T, X2[np.newaxis,:]), V), LA.inv(np.diag(S)))

#Trouver la valeur propre et le vecteur propre d'Atilde
Lam, W = LA.eig(Atilde)

print("valeur propre:",Lam)#Lam:1.00025541 est un nombre réel

#Trouver le vecteur propre de A à partir du vecteur propre d'Atilde
Phi = np.dot(np.dot(np.dot(X2[np.newaxis,:], V), LA.inv(np.diag(S))), W)

#Exp discrète à continue(**)de**Cherchant
Omega = np.log(Lam)/dt

#Restaurez la fonction d'origine d'Omega et de Phi.
b = np.dot(LA.pinv(Phi * np.exp(Omega * t)).T, x)
x_dmd = b * Phi * np.exp(Omega * t)

plt.plot(t, x_dmd[0,:])
plt.show()

Exemple DMD qui fonctionne

Il est introduit dans Article explicatif, mais les données sont décalées dans le sens du temps. Créez une autre dimension et faites-en des données bidimensionnelles. En faisant cela, un composant complexe apparaît dans la valeur propre, et $ \ sin (x) $ peut être trouvé dans les données numériques. C'est tout.

import numpy as np
import numpy.linalg as LA
import matplotlib.pyplot as plt
%matplotlib inline

t = np.arange(0, 10.01, 0.01)
x = np.sin(t)
dt = t[2] - t[1]

#Matrice de données
#La principale différence est ici. C'est bidimensionnel.
X_aug1 = np.array([x[:-2], x[1:-1]])
X_aug2 = np.array([x[1:-1], x[2:]])

#SVD
U,S,Vh = LA.svd(X_aug1,False)
V = Vh.conj().T

#Atilde avec A réduit en dimension par le vecteur singulier gauche
#U prend un conjugué complexe pour créer une matrice de translocation.
Atilde = np.dot(np.dot(np.dot(U.conj().T, X_aug2), V), LA.inv(np.diag(S)))

#Trouver la valeur propre et le vecteur propre d'Atilde
Lam, W = LA.eig(Atilde)

print("valeur propre:", Lam)#Lam:[0.9995+0.00999983j 0.9995-0.00999983j]Nombre complexe

#Trouver le vecteur propre de A à partir du vecteur propre d'Atilde
Phi = np.dot(np.dot(np.dot(X_aug2, V), LA.inv(np.diag(S))), W)

#Exp discrète à continue(**)de**Cherchant
Omega = np.diag(np.log(Lam)/dt)

print("Omega:", Omega)#Juste au cas où, ce sera 1j

#Restaurez la fonction d'origine d'Omega et de Phi.
#Ce que vous faites est le même que dans 1D, mais le style d'écriture est différent.
b = np.dot(LA.pinv(Phi), X_aug1[:,0])
x_dmd = np.zeros([2,len(t)], dtype='complex')
for i, _t in enumerate(t):
    x_dmd[:,i] = np.dot(np.dot(Phi, np.exp(Omega * _t)), b)

plt.plot(t, x_dmd[0,:].real)
plt.show()

Les références