Matrixexponential

In der Mathematik ist das Matrixexponential, auch als Matrix-Exponentialfunktion bezeichnet, eine Funktion auf der Menge der quadratischen Matrizen, welche analog zur gewöhnlichen (skalaren) Exponentialfunktion definiert ist. Das Matrixexponential stellt die Verbindung zwischen Lie-Algebra und der zugehörigen Lie-Gruppe her.

Definition

Sei eine reelle oder komplexe $(n \times n)$ -Matrix. Das Exponential von , welches mit $e^{X}$ oder $\exp(X)$ bezeichnet wird, ist die $(n \times n)$ -Matrix, welche durch die folgende Potenzreihe definiert ist (Taylor-Entwicklung):

$e^{X}=\sum _{{k=0}}^{\infty }{\frac {X^{k}}{k!}}$ .

Diese Reihe konvergiert immer. Daher ist das Exponential von wohldefiniert. Wenn eine $(1 \times 1)$ -Matrix ist, entspricht das Matrixexponential von der gewöhnlichen Exponentialfunktion. Eine Verallgemeinerung, welche auch für unendliche Matrizen sinnvoll ist, ist die Exponentialfunktion auf beliebigen Banachalgebren.

Eigenschaften

Das Matrixexponential teilt eine Reihe der Eigenschaften der gewöhnlichen Exponentialfunktion. Beispielsweise ist das Exponential der $(n \times n)$ -Nullmatrix $0$ gleich der $(n \times n)$ -Einheitsmatrix :

$e^{0}=E$ .

Für beliebige komplexe $(n \times n)$ -Matrizen und beliebige komplexe Zahlen und gilt

$e^{{aX}}\cdot e^{{bX}}=e^{{(a+b)X}}$ .

Daraus folgt

$e^{{X}}\cdot e^{{-X}}=e^{{(1-1)X}}=e^{0}=E$ ,

das heißt

$\left(e^{X}\right)^{{-1}}=e^{{-X}}$ .

Dabei bezeichnet $\left(e^{X}\right)^{{-1}}$ die zu $e^{X}$ inverse Matrix.

Die Exponentialfunktion erfüllt $e^{{x+y}}=e^{x}\,e^{y}$ für alle Zahlen und . Dasselbe gilt für kommutierende Matrizen und , das heißt aus

$X\cdot Y=Y\cdot X$

folgt

$e^{{X+Y}}=e^{X}\cdot e^{Y}$ .

Für nichtkommutierende Matrizen stimmt diese Gleichung im Allgemeinen nicht. In diesem Fall kann man $e^{{X+Y}}$ mit Hilfe der Baker-Campbell-Hausdorff-Formel berechnen.

Das Exponential der zu transponierten Matrix ist gleich der Transposition des Exponentials von :

$\exp \left(X^{{\mathrm {T}}}\right)=\left(\exp X\right)^{{\mathrm {T}}}$

Daraus folgt, dass die Matrixexponentialfunktion symmetrische Matrizen auf symmetrische Matrizen und schiefsymmetrische Matrizen auf orthogonale Matrizen abbildet. Analog gilt zwischen Adjunktion und Exponentiation die Beziehung

$\exp \left(X^{*}\right)=\left(\exp X\right)^{*}$ ,

so dass die Matrixexponentialfunktion hermitesche Matrizen auf hermitesche Matrizen und schiefhermitesche Matrizenauf unitäre Matrizen abbildet.

Weiterhin gilt:

Wenn invertierbar ist, dann ist $e^{{YXY^{{-1}}}}=Ye^{X}Y^{{-1}}$ .
$\det(e^{X})=e^{{{\mbox{tr}}(X)}}$ , Hier bezeichnet ${\mbox{tr}}(X)$ die Spur der quadratischen Matrix .
$e^{{\operatorname {diag}(x_{{1}},\ldots ,x_{{n}})}}=\operatorname {diag}\left(e^{{x_{{1}}}},\ldots ,e^{{x_{{n}}}}\right)$ .

Die Exponentialabbildung

Das Exponential einer Matrix ist immer eine invertierbare Matrix. Die Inverse von $e^{X}$ ist durch $e^{{-X}}$ gegeben. Das (komplexe) Matrixexponential liefert somit eine Abbildung

$\exp \colon M_{n}({\mathbb C})\to {\mbox{GL}}(n,{\mathbb C})$

aus dem Vektorraum aller (komplexen) $(n \times n)$ -Matrizen in die allgemeine lineare Gruppe, die Gruppe aller (komplexen) invertierbaren Matrizen. Diese Abbildung ist surjektiv, das heißt jede (reelle oder komplexe) invertierbare Matrix kann als die Exponentialmatrix einer komplexen Matrix geschrieben werden. Urbilder (bzw. lokale Schnitte) lassen sich durch Matrixlogarithmen berechnen.

Für je zwei Matrizen und gilt

$\|e^{{X+Y}}-e^{X}\|\leq \|Y\|e^{{\|X\|}}e^{{\|Y\|}}$ ,

wobei $\|\cdot \|$ eine beliebige Matrixnorm bezeichnet. Daraus folgt, dass die Exponentialabbildung stetig und auf kompakten Teilmengen von $M_{n}({\mathbb C})$ sogar lipschitzstetig ist. Für die Norm des Matrixexponentials selbst gibt es aber eine präzisere Schranke

$\|e^{X}\|\leq e^{{\mu (X)}}$

mit der logarithmischen Matrixnorm $\mu (A)$ und dem numerischen Wertebereich.

Die Zuordnung

$t\mapsto e^{{tX}},\qquad t\in {\mathbb R},$

definiert eine glatte Kurve in der allgemeinen linearen Gruppe, welche für t=0 die Einheitsmatrix liefert. Dies liefert eine Einparameter-Untergruppe der allgemeinen linearen Gruppe, da

$e^{{tX}}e^{{sX}}=e^{{(t+s)X}}$

gilt. Die Ableitung dieser Funktion am Punkt ist durch

${\frac {d}{dt}}e^{{tX}}=Xe^{{tX}}\qquad (1)$

gegeben. Die Ableitung für t=0 ist gerade die Matrix , das heißt erzeugt diese Einparameter-Untergruppe.

Allgemeiner gilt:

${\frac {d}{dt}}e^{{X(t)}}=\int _{0}^{1}e^{{(1-\alpha )X(t)}}{\frac {dX(t)}{dt}}e^{{\alpha X(t)}}d\alpha$

Beispiele von Lie-Algebren und zugehörigen Lie-Gruppen

Lie-Gruppe	Beispiel
Allgemeine lineare Gruppe: $\mathrm {GL} (n,K)$	${\mathfrak {gl}}(n)=\{X\in M_{n}(\mathbb {R} )\}$ $\mathrm {GL} (n)=\{X\in M_{n}(\mathbb {R} )\vert \det(X)\neq 0\}$ .
Orthogonale Gruppe: ${\mathrm {O}}(n,K)$	${\mathfrak {o}}(n)=\{X\in M_{n}(\mathbb {R} )\vert X^{T}=-X\}$ $\mathrm {O} (n)=\{X\in M_{n}(\mathbb {R} )\vert X^{T}=X^{-1}\}$
Unitäre Gruppe: ${\mathrm {U}}(n)$	${\mathfrak {u}}(n)=\{X\in M_{n}(\mathbb {C} )\vert X^{}=-X\}$ $\mathrm {U} (n)=\{X\in M_{n}(\mathbb {C} )\vert X^{}=X^{-1}\}$
Spezielle unitäre Gruppe: ${\mathrm {SU}}(n)$	${\mathfrak {su}}(2)=\{X\in M_{2}({\mathbb C})\vert X^{}=-X,\operatorname {tr}(X)=0\}$ wird von $\exp$ surjektiv auf ${\mathrm {SU}}(2)=\{A\in M_{2}({\mathbb C})\vert A^{}=A^{{-1}},\det(A)=1\}$ abgebildet.
Spezielle orthogonale Gruppe: ${\mathrm {SO}}(n,K)$	${\mathfrak {so}}(n,\mathbb {R} )={\mathfrak {o}}(n,\mathbb {R} )=\{X\in M_{n}(\mathbb {R} )\vert X^{T}=-X\}$ (schiefsymmetrische Matrizen) wird von $\exp$ surjektiv auf ${\mathrm {SO}}(n,{\mathbb R})=\{A\in M_{n}({\mathbb R})\vert A^{T}=A^{{-1}},\det(A)=1\}$ abgebildet.
Spezielle lineare Gruppe: ${\mathrm {SL}}(n,K)$	${\mathfrak {sl}}(2,{\mathbb C})$ wird von $\exp$ nicht surjektiv auf ${\mathrm {SL}}(2,{\mathbb C})$ abgebildet. Notorisches Gegenbeispiel ${\begin{pmatrix}-1&a\\0&-1\end{pmatrix}}\in {\mathrm {SL}}(2,{\mathbb C})$ mit $a\neq 0$ liegt nicht im Bild von ${\mathfrak {sl}}(2,{\mathbb C})$ .

Aus dem letzten Beispiel ist ersichtlich, dass die Exponentialabbildung für die Erzeugung von Lie-Gruppen (je nach Lie-Algebra) im Allgemeinen nicht surjektiv ist.

Lineare Differentialgleichungen

Einer der Vorzüge des Matrixexponentials ist, dass man es benutzen kann, um Systeme von linearen gewöhnlichen Differentialgleichungen zu lösen. Aus Gleichung (1) oben folgt zum Beispiel, dass die Lösung des Anfangswertproblems

${\frac {d}{dt}}y(t)=Ay(t),\quad y(t_{0})=y_{0}$ ,

wobei eine Matrix ist, durch

$y(t)=e^{A(t-t_{0})}y_{0}$

gegeben ist.

Das Matrixexponential kann auch zur Lösung der inhomogenen Gleichung

${\frac {d}{dt}}y(t)=Ay(t)+z(t),\quad y(t_{0})=y_{0}$ ,

verwendet werden. Beispiele findet man unten im Kapitel Anwendungen.

Für Differentialgleichungen der Form

${\frac {d}{dt}}y(t)=A(t)\,y(t),\quad y(t_{0})=y_{0}$

mit nicht-konstantem gibt es im Allgemeinen keine geschlossenen Lösungen. Die Magnus-Reihe liefert jedoch eine allgemeine Lösung in Matrixschreibweise über die Matrix-Exponentialfunktion auch im Fall nicht-konstanter Koeffizienten (als unendliche Reihe des Exponenten).

Berechnung des Matrixexponentials

Taylor-Reihe

Die Matrix-Exponentialfunktion kann prinzipiell über ihre Taylor-Entwicklung berechnet werden

$\exp(X)=\sum _{k=0}^{\infty }{\frac {X^{k}}{k!}}=E+X+{\frac {X^{2}}{2}}+{\frac {X^{3}}{6}}+{\frac {X^{4}}{24}}+\cdots$

wobei die Fakultät von bezeichnet. Bei ausreichender Genauigkeit (Reihe ist absolut konvergent) soll die Reihe bei einer endlichen Zahl an Berechnungsschritten abbrechen. Je größer die Einträge der Matrix sind, desto mehr Glieder der Reihe müssen aber berechnet werden (z.B. für die Lösung der linearen DGL für einen großen Zeitschritt). Um den Lösungsalgorithmus dahingehend zu verbessern, kann man die Einträge der Matrix mittels der Rechenregel $e^{x}=(e^{x/n})^{n}$ elegant skalieren. Die Berechnung der Reihe an sich kann ebenfalls über die Padé-Approximation verbessert werden.

Diagonalisierbare Matrizen

Ist die Matrix eine Diagonalmatrix

$A={\begin{pmatrix}a_{1}&0&\ldots &0\\0&a_{2}&\ldots &0\\\vdots &\vdots &\ddots &\vdots \\0&0&\ldots &a_{n}\end{pmatrix}}$ ,

dann kann man ihr Exponential ermitteln, indem man die übliche Exponentialfunktion auf jeden Eintrag der Hauptdiagonalen anwendet:

$e^{A}={\begin{pmatrix}e^{a_{1}}&0&\ldots &0\\0&e^{a_{2}}&\ldots &0\\\vdots &\vdots &\ddots &\vdots \\0&0&\ldots &e^{a_{n}}\end{pmatrix}}$ .

Damit kann man auch das Exponential diagonalisierbarer Matrizen berechnen. Wenn

$A=UDU^{{-1}}$

mit einer Diagonalmatrix und einer Basiswechselmatrix ist, dann gilt:

$e^{A}=Ue^{D}U^{{-1}}$ .

Nilpotenter Fall

Eine Matrix ist nilpotent, wenn $N^{q}=0$ für eine geeignete natürliche Zahl gilt. In diesem Fall kann das Matrixexponential $e^{N}$ aus der Reihenentwicklung berechnet werden, wobei die Reihe nach einer endlichen Anzahl von Termen abbricht:

$e^{N}=E+N+{\frac {1}{2}}N^{2}+{\frac {1}{6}}N^{3}+\cdots +{\frac {1}{(q-1)!}}N^{{q-1}}$ .

Allgemeiner Fall (Eigenwertproblem)

Zerfällt das Minimalpolynom (bzw. das charakteristische Polynom) der Matrix in Linearfaktoren (über $\mathbb {C}$ ist das stets der Fall), dann kann eindeutig in eine Summe

zerlegt werden, wobei

diagonalisierbar ist,
nilpotent ist und
mit kommutiert (d.h. ).

Damit kann man das Exponential von berechnen, indem man es auf die vorgenannten Fälle reduziert: $e^{X}=e^{{A+N}}=e^{A}e^{N}$ . Im letzten Schritt benötigt man die Kommutativität von und .

Eine andere (nah verwandte) Methode ist die Verwendung der Jordanschen Normalform von . Sei die Jordansche Normalform von mit der Basiswechselmatrix , das heißt, es gilt

$e^{{X}}=Pe^{{J}}P^{{-1}}.$

Wegen

$J=J_{{a_{1}}}(\lambda _{1})\oplus J_{{a_{2}}}(\lambda _{2})\oplus \cdots \oplus J_{{a_{n}}}(\lambda _{n})$

gilt

${\begin{aligned}e^{J}&=\exp {\big (}J_{a_{1}}(\lambda _{1})\oplus J_{a_{2}}(\lambda _{2})\oplus \cdots \oplus J_{a_{n}}(\lambda _{n}){\big )}\\&=\exp {\big (}J_{a_{1}}(\lambda _{1}){\big )}\oplus \exp {\big (}J_{a_{2}}(\lambda _{2}){\big )}\oplus \cdots \oplus \exp {\big (}J_{a_{k}}(\lambda _{k}){\big )}.\end{aligned}}$

Daher muss man nur das Exponential eines Jordan-Blocks kennen. Nun ist jeder Jordan-Block von der Form

$J_{a}(\lambda )=\lambda I+N,$

wobei eine spezielle nilpotente Matrix ist. Das Exponential des Jordan-Blocks ist also

$e^{{\lambda I+N}}=e^{{\lambda }}e^{N}.$

Beispiel 1

Man betrachte die Matrix

$B={\begin{pmatrix}21&17&6\\-5&-1&-6\\4&4&16\end{pmatrix}}$ ,

welche die Jordansche Normalform

$J=P^{{-1}}BP={\begin{pmatrix}4&0&0\\0&16&1\\0&0&16\end{pmatrix}}$

mit der Übergangsmatrix

$P={\begin{pmatrix}-1&1&{\tfrac {5}{8}}\\1&-1&-{\tfrac {1}{4}}\\0&2&0\end{pmatrix}}$

hat. Dann gilt

$J=J_{1}(4)\oplus J_{2}(16)$

und

$e^{B}=Pe^{{J}}P^{{-1}}=P(e^{{J_{2}(16)}}\oplus e^{{J_{1}(4)}})P^{{-1}}$ .

Somit ist

$\exp \left[16E+{\begin{pmatrix}0&1\\0&0\end{pmatrix}}\right]=e^{16}\left[{\begin{pmatrix}1&0\\0&1\end{pmatrix}}+{\begin{pmatrix}0&1\\0&0\end{pmatrix}}+{1 \over 2!}{\begin{pmatrix}0&0\\0&0\end{pmatrix}}+\cdots \right]={\begin{pmatrix}e^{16}&e^{16}\\0&e^{16}\end{pmatrix}}$ .

Das Exponential einer 1×1-Matrix ist trivial. Mit $e^{{J_{1}(4)}}=e^{4}$ folgt

${\begin{aligned}\exp(B)&=P\exp(J)P^{-1}=P{\begin{pmatrix}e^{4}&0&0\\0&e^{16}&e^{16}\\0&0&e^{16}\end{pmatrix}}P^{-1}\\[6pt]&={1 \over 4}{\begin{pmatrix}13e^{16}-e^{4}&13e^{16}-5e^{4}&2e^{16}-2e^{4}\\-9e^{16}+e^{4}&-9e^{16}+5e^{4}&-2e^{16}+2e^{4}\\16e^{16}&16e^{16}&4e^{16}\end{pmatrix}}.\end{aligned}}$

Die Jordansche Normalform und daraus das Exponential zu berechnen, ist auf diesem Weg sehr mühsam. Meist reicht es, die Wirkung der Exponential-Matrix auf einige Vektoren zu berechnen.

Beispiel 2

Es soll folgende Matrix-Exponentialfunktion berechnet werden:

$\exp(A\,t)=\exp \left[{\begin{pmatrix}0&1\\a_{2,1}&a_{2,2}\end{pmatrix}}\,t\right]$

Hierzu wird die $(2\times 2)$ -Matrix zunächst mittels der Eigenwerte und den Eigenvektoren diagonalisiert. Mit der Diagonalmatrix und der Eigenbasis folgt:

$\exp(A\,t)=V\exp(D\,t)\,V^{-1}$

Die Eigenwerte werden aus dem charakteristischen Polynom 2ten Grades bestimmt zu

$\lambda _{1,2}={\frac {a_{2,2}}{2}}\pm {\sqrt {\left({\frac {a_{2,2}}{2}}\right)^{2}+a_{2,1}}}$ .

Für die beiden Eigenvektoren bzw. die Eigenbasis gilt damit:

$V={\begin{pmatrix}1&1\\\lambda _{1}&\lambda _{2}\end{pmatrix}}$

Einsetzen für die Matrix-Exponantialfunktion liefert schließlich

$\exp(A\,t)={\frac {1}{\lambda _{2}-\lambda _{1}}}{\begin{pmatrix}\lambda _{2}e^{\lambda _{1}t}-\lambda _{1}e^{\lambda _{2}t}&e^{\lambda _{2}t}-e^{\lambda _{1}t}\\\lambda _{1}\lambda _{2}(e^{\lambda _{1}t}-e^{\lambda _{2}t})&\lambda _{2}e^{\lambda _{2}t}-\lambda _{1}e^{\lambda _{1}t}\end{pmatrix}}$

als geschlossene analytische Lösung.

Numerische Verfahren

Die Jordan-Normalform-Zerlegung ist numerisch instabil, da aufgrund der Gleitkommaarithmetik Rundungsfehler in die Eigenwerte eingeführt werden, die eine Gruppierung der Eigenwerte in Gruppen identischer Eigenwerte unmöglich macht. Daher werden in der Numerik andere Techniken zur Berechnung des Matrixexponentials verwendet. Einer der effektivsten verfügbaren Algorithmen ist die Padé-Approximation mit Skalieren und Quadrieren (s. Berechnung mittels Taylorreihe). Bei großen Matrizen kann der Rechenaufwand zusätzlich reduziert werden, indem Krylowräume verwendet werden, deren Basisvektoren mit dem Arnoldi-Verfahren orthogonalisiert worden sind.

Anwendungen

Homogene lineare Differentialgleichungen

Das Matrixexponential kann für die Lösung eines Systems von linearen Differentialgleichungen verwendet werden. Eine Differentialgleichung der Form

hat die Lösung $e^{{Ct}}$ . Wenn man den Vektor

${\mathbf {y}}(t)={\begin{pmatrix}y_{1}(t)\\\vdots \\y_{n}(t)\end{pmatrix}}$

betrachtet, dann kann man ein System von gekoppelten linearen Differentialgleichungen betrachten als

${\mathbf {y}}'(t)=A{\mathbf {y}}(t)+{\mathbf {b}}$ .

Wenn man den Integrationsfaktor $e^{{-tA}}$ ansetzt und auf beiden Seiten multipliziert, erhält man

$e^{{-tA}}{\mathbf {y}}'(t)-e^{{-tA}}A{\mathbf {y}}=e^{{-tA}}{\mathbf {b}}$ ,

$D(e^{{-tA}}{\mathbf {y}})=e^{{-tA}}{\mathbf {b}}$ .

Wenn man $e^{tA}$ berechnet, erhält man eine Lösung des Differentialgleichungssystems.

Beispiel (homogen)

Gegeben sei das folgende Differentialgleichungssystem

${\begin{aligned}y_{1}'&=3y_{1}-y_{2},\\y_{2}'&=y_{1}+y_{2}.\end{aligned}}$

Es lässt sich schreiben als y'(t)=Ay(t) mit der Koeffizientenmatrix

$A={\begin{pmatrix}3&-1\\1&1\end{pmatrix}}$ .

Damit ergibt sich das zugehörige Matrixexponential zu

$e^{tA}=e^{2t}{\begin{pmatrix}1+t&-t\\t&1-t\end{pmatrix}}$ .

Als allgemeine Lösung des Differentialgleichungssystems erhält man somit

${\begin{pmatrix}y_{1}(t)\\y_{2}(t)\end{pmatrix}}=C_{1}e^{{2t}}{\begin{pmatrix}1+t\\t\end{pmatrix}}+C_{2}e^{{2t}}{\begin{pmatrix}-t\\1-t\end{pmatrix}}$ .

Inhomogener Fall – Variation der Konstanten

Für den inhomogenen Fall kann man eine Methode ähnlich der Variation der Konstanten benutzen. Es wird eine Lösung der Form ${\mathbf {y}}_{p}(t)=e^{{tA}}{\mathbf {z}}(t)$ gesucht:

${\mathbf {y}}_{p}'=(e^{{tA}})'{\mathbf {z}}(t)+e^{{tA}}{\mathbf {z}}'(t)=Ae^{{tA}}{\mathbf {z}}(t)+e^{{tA}}{\mathbf {z}}'(t)=A{\mathbf {y}}_{p}(t)+e^{{tA}}{\mathbf {z}}'(t)$

Um die Lösung ${\mathbf {y}}_{p}$ zu ermitteln, setzt man

$e^{{tA}}{\mathbf {z}}'(t)={\mathbf {b}}(t)$

${\mathbf {z}}'(t)=(e^{{tA}})^{{-1}}{\mathbf {b}}(t)$

${\mathbf {z}}(t)=\int _{0}^{t}e^{{-uA}}{\mathbf {b}}(u)\,du+{\mathbf {c}}$

Damit ergibt sich

${\mathbf {y}}_{p}=e^{{tA}}\int _{0}^{t}e^{{-uA}}{\mathbf {b}}(u)\,du+e^{{tA}}{\mathbf {c}}=\int _{0}^{t}e^{{(t-u)A}}{\mathbf {b}}(u)\,du+e^{{tA}}{\mathbf {c}}$ ,

wobei durch die Anfangsbedingungen bestimmt wird.

Beispiel (inhomogen)

Gegeben sei das Differentialgleichungssystem

${\begin{aligned}y_{1}'&=3y_{1}-y_{2}-2e^{t},\\y_{2}'&=y_{1}+y_{2}-e^{t}.\end{aligned}}$

Mit der Matrix von oben schreibt sich das System

mit

$b(t)=e^{{t}}{\begin{pmatrix}-2\\-1\end{pmatrix}}$ .

Die allgemeine Lösung der homogenen Gleichung wurde bereits oben berechnet. Die Summe aus homogenen und speziellen Lösungen ergibt die Lösung für das inhomogene Problem. Man muss jetzt nur noch eine spezielle Lösung $y_{p}(t)$ finden (über die Variation der Konstanten). Von der Gleichung oben erhält man:

$y_{p}(t)=e^{{tA}}\int _{0}^{t}e^{{-uA}}b(u)\,du$ ,

also

$y_{p}(t)=e^{2t}{\begin{pmatrix}1+t&-t\\t&1-t\end{pmatrix}}\int _{0}^{t}e^{-2u}{\begin{pmatrix}1-u&u\\-u&1+u\end{pmatrix}}{\begin{pmatrix}-2e^{u}\\-e^{u}\end{pmatrix}}\,du={\begin{pmatrix}-e^{2t}(1+t)+e^{t}\\-te^{2t}\end{pmatrix}}$ .

Siehe auch

Basierend auf einem Artikel in:

Wikipedia.de