Matrixexponential

Grundidee und Definition

Das Matrixexponential ist eine Matrixfunktion, die ähnlich wie die gewöhnliche Exponentialfunktion definiert ist. Für eine reelle oder komplexe n×n-Matrix X wird ihr Exponential mit e^X oder exp(X) bezeichnet. Es ist die n×n-Matrix

e^X = ∑_{k=0}^{∞} X^k/k! = E + X + X^2/2! + X^3/3! + ...

Dabei ist E die Einheitsmatrix. Diese Potenzreihe ist die Maclaurin-Reihe beziehungsweise die Taylor-Entwicklung an der Stelle x = 0. Sie konvergiert immer, genau wie die Reihe der gewöhnlichen Exponentialfunktion. Deshalb ist e^X für jede reelle oder komplexe quadratische Matrix X wohldefiniert.

Wenn X eine 1×1-Matrix ist, stimmt das Matrixexponential mit der normalen Exponentialfunktion überein. Eine Verallgemeinerung für unendliche Matrizen ist die Exponentialfunktion auf Banachalgebren. Besonders wichtig ist das Matrixexponential, weil es eine Verbindung zwischen Lie-Algebren und den zugehörigen Lie-Gruppen herstellt.

Wichtige Eigenschaften

Das Matrixexponential besitzt viele Eigenschaften der gewöhnlichen Exponentialfunktion, aber nicht alle gelten ohne Einschränkung. Für die n×n-Nullmatrix 0 gilt

e^0 = E.

Für jede komplexe n×n-Matrix X und komplexe Zahlen a und b gilt

e^{aX} · e^{bX} = e^{(a+b)X}.

Daraus folgt insbesondere

e^X · e^{-X} = E,

also ist e^X immer invertierbar und (e^X)^{-1} = e^{-X}.

Für Zahlen gilt immer e^{x+y} = e^x e^y. Für Matrizen gilt die entsprechende Regel nur dann allgemein, wenn die Matrizen kommutieren, also XY = YX. Dann gilt

e^{X+Y} = e^X · e^Y.

Für nichtkommutierende Matrizen gilt diese Gleichung im Allgemeinen nicht. Dann kann e^{X+Y} mit der Baker-Campbell-Hausdorff-Formel berechnet werden.

Weitere wichtige Regeln sind:

• exp(X^T) = (exp X)^T. Deshalb werden symmetrische Matrizen auf symmetrische Matrizen und schiefsymmetrische Matrizen auf orthogonale Matrizen abgebildet.
• exp(X*) = (exp X)*. Deshalb werden hermitesche Matrizen auf hermitesche Matrizen und schiefhermitesche Matrizen auf unitäre Matrizen abgebildet.
• Wenn Y invertierbar ist, dann gilt e^{YXY^{-1}} = Y e^X Y^{-1}.
• det(e^X) = e^{tr(X)}, wobei tr(X) die Spur von X ist.
• Für Diagonalmatrizen gilt e^{diag(x_1,...,x_n)} = diag(e^{x_1},...,e^{x_n}).

Exponentialabbildung und Lie-Gruppen

Da e^X immer invertierbar ist, liefert das komplexe Matrixexponential eine Abbildung

exp: M_n(C) → GL(n,C).

Hier ist M_n(C) der Vektorraum aller komplexen n×n-Matrizen und GL(n,C) die allgemeine lineare Gruppe, also die Gruppe aller komplexen invertierbaren n×n-Matrizen. Diese Abbildung ist auf ganz M_n(C) kein Gruppenhomomorphismus, aber auf jeder Untergruppe, deren Matrizen multiplikativ miteinander kommutieren.

Die Abbildung ist surjektiv im Sinn des Artikels: Jede reelle oder komplexe invertierbare Matrix kann als Exponentialmatrix einer komplexen Matrix geschrieben werden. Urbilder beziehungsweise lokale Schnitte lassen sich durch Matrixlogarithmen bestimmen.

Für Matrizen X und Y gilt außerdem die Abschätzung

||e^{X+Y} - e^X|| ≤ ||Y|| e^{||X||} e^{||Y||},

wobei ||·|| eine Matrixnorm ist. Daraus folgt, dass die Exponentialabbildung stetig und auf kompakten Teilmengen sogar lipschitzstetig ist. Für e^X selbst gibt es die genauere Schranke

||e^X|| ≤ e^{μ(X)},

wobei μ die logarithmische Matrixnorm bezeichnet.

Für eine feste Matrix X definiert

t ↦ e^{tX}, t ∈ R,

eine glatte Kurve in GL(n,C), die bei t = 0 die Einheitsmatrix ergibt. Wegen e^{tX}e^{sX} = e^{(t+s)X} ist dies eine Einparameter-Untergruppe. Ihre Ableitung ist

d/dt e^{tX} = X e^{tX}.

Für t = 0 ist die Ableitung genau X; man sagt daher, X erzeugt diese Einparameter-Untergruppe. Allgemeiner gilt für eine matrixwertige Funktion X(t):

d/dt e^{X(t)} = ∫_0^1 e^{(1-α)X(t)} (dX(t)/dt) e^{αX(t)} dα.

Beispiele für Lie-Algebren und Lie-Gruppen sind unter anderem gl(n)=M_n(R) zur allgemeinen linearen Gruppe GL(n), schiefsymmetrische Matrizen zur orthogonalen Gruppe O(n), schiefhermitesche Matrizen zur unitären Gruppe U(n), sowie su(2), das durch exp surjektiv auf SU(2) abgebildet wird. Für so(n,R) wird exp surjektiv auf SO(n,R) abgebildet. Dagegen ist exp von sl(2,C) nicht surjektiv auf SL(2,C); als Gegenbeispiel wird die Matrix ((-1,a),(0,-1)) mit a ≠ 0 genannt.

Lineare Differentialgleichungen

Eine wichtige Anwendung des Matrixexponentials ist die Lösung linearer gewöhnlicher Differentialgleichungssysteme. Für das Anfangswertproblem

d/dt y(t) = A y(t), y(t_0) = y_0,

mit einer quadratischen Matrix A lautet die Lösung

y(t) = e^{A(t-t_0)} y_0.

Das folgt aus der Ableitungsregel d/dt e^{tA} = A e^{tA}. Das Matrixexponential wird deshalb zum Beispiel bei Zustandsraummodellen dynamischer Übertragungssysteme verwendet.

Auch inhomogene lineare Differentialgleichungen der Form

d/dt y(t) = A y(t) + z(t), y(t_0)=y_0,

können mit dem Matrixexponential behandelt werden. Für Differentialgleichungen mit nichtkonstanter Matrix A(t), also

d/dt y(t) = A(t)y(t),

gibt es im Allgemeinen keine geschlossenen Lösungen. Die Magnus-Reihe liefert aber eine allgemeine Lösung in Matrixschreibweise über die Matrix-Exponentialfunktion, wobei der Exponent als unendliche Reihe dargestellt wird.

Für homogene Systeme y' = Cy ist e^{Ct} die Lösung. Für ein gekoppeltes System y'(t)=A y(t)+b kann man den Integrationsfaktor e^{-tA} verwenden. Im inhomogenen Fall benutzt man die Variation der Konstanten und sucht eine Lösung der Form y_p(t)=e^{tA}z(t). Daraus ergibt sich bei Anfangsbedingung y_p(t_0)=y_0:

y_p(t)=e^{(t-t_0)A}y_0 + ∫_{t_0}^{t} e^{(t-u)A} b(u) du.

Berechnungsmethoden

Das Matrixexponential kann auf verschiedene Arten berechnet werden. Direkt möglich ist die Taylor-Reihe

exp(At)=∑_{k=0}^{∞} A^k t^k/k! = E + At + A^2t^2/2 + A^3t^3/6 + ...

Die Reihe ist absolut konvergent, aber bei großen Matrixeinträgen müssen viele Terme berechnet werden. Eine Verbesserung ist die Methode „Scaling & Squaring“ mit der Regel e^{At} = (e^{At/m})^m. Wenn ||At|| nicht zu groß ist, kann die Reihe auch über eine Padé-Approximation berechnet werden. Scaling & Squaring hat im Wesentlichen Aufwand O(n^3), weil hauptsächlich Matrizen multipliziert werden.

Eine Matrix N heißt nilpotent, wenn N^q = 0 für eine natürliche Zahl q gilt. Dann endet die Reihe nach endlich vielen Gliedern:

exp(N)=E+N+1/2 N^2+1/6 N^3+...+1/(q-1)! N^{q-1}.

Bei Diagonalmatrizen wird einfach die normale Exponentialfunktion auf die Diagonaleinträge angewendet. Ist A diagonalisierbar, also A = VDV^{-1}, dann gilt

exp(At)=V e^{Dt} V^{-1}=V diag(e^{λ_1t}, e^{λ_2t}, ..., e^{λ_nt}) V^{-1}.

Dabei sind λ_1,...,λ_n die Eigenwerte. Die Diagonalisierung und auch der QR-Algorithmus haben Aufwand O(n^3), sind aber im Vergleich zu Taylor-basierten Methoden unabhängig von ||At||. Der große Aufwand entsteht vor allem bei der Bestimmung von Eigenwerten und Eigenvektoren und muss bei mehreren Zeitpunkten t nur einmal geleistet werden.

Weitere Methoden sind die Splitting-Methode und die jordansche Normalform. Bei der Splitting-Methode wird X = A + N zerlegt, wobei A diagonalisierbar, N nilpotent ist und AN = NA gilt. Dann ist e^X = e^A e^N. Bei der jordanschen Normalform gilt für X mit J = P^{-1}XP:

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

Da ein Jordan-Block die Form J_a(λ)=λE+N hat, lässt sich sein Exponential als e^{λE+N}=e^λ e^N berechnen. Numerisch ist die Jordan-Normalform jedoch instabil, weil Rundungsfehler in den Eigenwerten die Gruppierung identischer Eigenwerte erschweren. In der Numerik verwendet man daher eher Padé-Approximation mit Skalieren und Quadrieren, Diagonalisierung oder für große Matrizen Krylowräume mit Arnoldi-Verfahren.

Für 2×2-Matrizen gibt es explizite Formeln über die Diskriminante D des charakteristischen Polynoms. Je nachdem, ob D > 0, D < 0 oder D = 0 ist, werden sinh/cosh, sin/cos oder f(t)=t und f'(t)=1 verwendet. Der Putzer-Algorithmus stellt exp(At) als Polynom dar:

exp(At)=∑{k=1}^{n} p_k(t) M{k-1}.

Er hat Aufwand O(n^4) und eignet sich daher eher für kleine Matrizen.

Typische Beispiele und Anwendungen

Ein homogenes Beispiel ist das Differentialgleichungssystem

y_1' = 3y_1 - y_2, y_2' = y_1 + y_2.

Es wird als y'(t)=Ay(t) geschrieben mit

A = ((3,-1),(1,1)).

Das zugehörige Matrixexponential ist

e^{tA}=e^{2t} ((1+t,-t),(t,1-t)).

Damit erhält man die allgemeine Lösung

(y_1(t),y_2(t))^T = C_1 e^{2t}(1+t,t)^T + C_2 e^{2t}(-t,1-t)^T.

Ein inhomogenes Beispiel ist

y_1' = 3y_1 - y_2 - 2e^t, y_2' = y_1 + y_2 - e^t.

Mit derselben Matrix A gilt y'(t)=Ay(t)+b(t) und

b(t)=e^t(-2,-1)^T.

Über Variation der Konstanten erhält man eine spezielle Lösung

y_p(t)=e^{tA}∫_0^t e^{-uA}b(u)du = (-e^{2t}(1+t)+e^t, -t e^{2t})^T.

Eine weitere Anwendung ist die Tensorproduktformel. Für die Kroneckersumme zweier quadratischer Matrizen A und B gilt

A ⊕̄ B = A⊗I_m + I_n⊗B

und

exp(A ⊕̄ B)=exp(A)⊗exp(B).

Hier bezeichnet ⊗ das Tensorprodukt beziehungsweise Kronecker-Produkt. Diese Formel wird in der Darstellungstheorie direkter Summen von Lie-Algebren verwendet. In der Physik tritt sie in der Quantenmechanik nichtwechselwirkender Systeme auf, weil deren Zeitentwicklung über die Exponentialabbildung unabhängig für die einzelnen Bestandteile verläuft. Außerdem wird die umgekehrte Form log(A) ⊕̄ log(B)=log(A⊗B) bei der Von-Neumann-Entropie genutzt, um die Additivität der Entropie für unabhängige Systeme zu zeigen.