Main Content

Die Übersetzung dieser Seite ist veraltet. Klicken Sie hier, um die neueste Version auf Englisch zu sehen.

Lineare Algebra

Matrizen in der MATLAB Umgebung

Dieser Abschnitt enthält eine Einführung in die Erstellung von Matrizen und die Ausführung grundlegender Matrixberechnungen in MATLAB®.

In der MATLAB Umgebung wird der Begriff Matrix verwendet, um eine Variable anzugeben, die in einem zweidimensionalen Raster angeordnete reelle oder komplexe Zahlen enthält. Ein Array ist, allgemeiner ausgedrückt, ein Vektor, eine Matrix oder ein Zahlenraster mit einer höheren Dimension. Alle Arrays in MATLAB sind rechteckig, das heißt, die Komponentenvektoren entlang einer beliebigen Dimension weisen alle dieselbe Länge auf. Die mathematischen Operationen, die für die Matrizen definiert sind, gehören zur linearen Algebra.

Erstellen von Matrizen

In MATLAB gibt es zahlreiche Funktionen, mit denen die unterschiedlichsten Matrizen erstellt werden können. Beispielsweise können Sie eine symmetrische Matrix mit Einträgen auf Basis des Pascalschen Dreiecks erstellen:

A = pascal(3)
A =
       1     1     1
       1     2     3
       1     3     6

Oder Sie können eine asymmetrische Matrix gemäß dem magischen Quadrat erstellen, die identische Zeilen- und Spaltensummen aufweist:

B = magic(3)
B =
       8     1     6
       3     5     7
       4     9     2

Ein weiteres Beispiel ist eine rechteckige 3x2-Matrix aus Zufallszahlen. In diesem Fall beschreibt die erste Eingabe für randi den Bereich der möglichen Werte für die Ganzzahlen, während die beiden zweiten Eingaben die Anzahl der Zeilen und Spalten beschreiben.

C = randi(10,3,2)
C =

     9    10
    10     7
     2     1

Ein Spaltenvektor ist eine mx1-Matrix, ein Zeilenvektor ist eine 1xn-Matrix und ein Skalar ist eine 1x1-Matrix. Verwenden Sie zum manuellen Definieren einer Matrix eckige Klammern [ ], um den Anfang und das Ende des Arrays anzugeben. Verwenden Sie innerhalb der Klammern einen Strichpunkt (;), um das Ende einer Zeile anzugeben. Im Falle eines Skalars (1x1-Matrix) sind die Klammern nicht erforderlich. Diese Anweisungen erstellen einen Spaltenvektor, einen Zeilenvektor und einen Skalar:

u = [3; 1; 4]

v = [2 0 -1]

s = 7
u =
       3
       1
       4

v =
       2     0    -1

s =
       7

Weitere Informationen zum Erstellen und Bearbeiten von Matrizen finden Sie unter Creating, Concatenating, and Expanding Matrices.

Addieren und Subtrahieren von Matrizen

Die Addition und Subtraktion von Matrizen und Arrays wird Element für Element, oder elementweise, ausgeführt. Wenn Sie beispielsweise A zu B addieren und anschließend A vom Ergebnis subtrahieren, wird B wiederhergestellt:

X = A + B
X =
       9     2     7
       4     7    10
       5    12     8
Y = X - A
Y =
       8     1     6
       3     5     7
       4     9     2

Für die Addition und Subtraktion müssen die Dimensionen beider Matrizen kompatibel sein. Sind die Dimensionen inkompatibel, ergibt sich ein Fehler:

X = A + C
Error using  + 
Matrix dimensions must agree.

Weitere Informationen finden Sie unter Array vs. Matrix Operations.

Vektorprodukte und Transposition

Ein Zeilenvektor und ein Spaltenvektor derselben Länge können in beliebiger Reihenfolge multipliziert werden. Das Ergebnis ist entweder ein Skalar, auch inneres Produkt genannt, oder eine Matrix, auch äußeres Produkt genannt:

u = [3; 1; 4];
v = [2 0 -1];
x = v*u
x =

     2
X = u*v
X =

     6     0    -3
     2     0    -1
     8     0    -4

Bei reellen Matrizen vertauscht die Rechenoperation für das Transponieren aij und aji. Bei komplexen Matrizen ist eine weitere Überlegung, ob das komplexe Konjugat komplexer Einträge im Array ermittelt wird, um die Transponierte des komplexen Konjugats zu bilden. MATLAB verwendet den Apostrophoperator ('), um das komplexe Konjugat zu transponieren, und den Punktoperator (.'), um das Transponieren ohne Konjugation auszuführen. Bei Matrizen, die ausschließlich reelle Elemente enthalten, führen beide Operatoren zum gleichen Ergebnis.

Die Beispielmatrix A = pascal(3) ist symmetrisch, sodass A' gleich A ist. Allerdings ist B = magic(3) nicht symmetrisch, sodass mit B' die Elemente entlang der Hauptdiagonalen widergespiegelt werden:

B = magic(3)
B =

     8     1     6
     3     5     7
     4     9     2
X = B'
X =

     8     3     4
     1     5     9
     6     7     2

Bei Vektoren wandelt die Transposition einen Zeilenvektor in einen Spaltenvektor um (und umgekehrt):

x = v'

x =
       2
       0
      -1

Wenn x und y beides reelle Spaltenvektoren sind, wird das Produkt x*y nicht definiert. Allerdings erstellen die beiden Produkte

x'*y

und

y'*x

dasselbe skalare Ergebnis. Diese Größe wird so häufig verwendet, dass sie drei verschiedene Namen hat: inneres Produkt, Skalarprodukt oder Punktprodukt. Es gibt sogar eine dedizierte Funktion für Punktprodukte mit dem Namen dot.

Für einen komplexen Vektor oder eine Matrix, z, transponiert die Größe z' nicht nur den Vektor oder die Matrix, sondern konvertiert auch alle komplexen Elemente in ihr komplexes Konjugat. Das bedeutet, dass sich das Vorzeichen des Imaginärteils der einzelnen komplexen Elemente ändert. Angenommen, es liegt die folgende komplexe Matrix vor:

z = [1+2i 7-3i 3+4i; 6-2i 9i 4+7i]
z =

   1.0000 + 2.0000i   7.0000 - 3.0000i   3.0000 + 4.0000i
   6.0000 - 2.0000i   0.0000 + 9.0000i   4.0000 + 7.0000i

Die komplexe Transponierte des Konjugats von z ist:

z'
ans =

   1.0000 - 2.0000i   6.0000 + 2.0000i
   7.0000 + 3.0000i   0.0000 - 9.0000i
   3.0000 - 4.0000i   4.0000 - 7.0000i

Die nicht konjugierte komplexe Transponierte, bei der der komplexe Teil jedes Elements sein Vorzeichen beibehält, wird mit z.' angegeben:

z.'
ans =

   1.0000 + 2.0000i   6.0000 - 2.0000i
   7.0000 - 3.0000i   0.0000 + 9.0000i
   3.0000 + 4.0000i   4.0000 + 7.0000i

Bei komplexen Vektoren sind die beiden skalaren Produkte x'*y und y'*x komplexe Konjugate voneinander und das skalare Produkt x'*x eines komplexen Vektors mit sich selbst ist reell.

Multiplizieren von Matrizen

Die Multiplikation von Matrizen ist so definiert, dass die Zusammensetzung der zugrunde liegenden Transformationen reflektiert wird, und erlaubt die kompakte Darstellung von Systemen simultaner, linearer Gleichungen. Das Matrixprodukt C = AB wird definiert, wenn die Spaltendimension von A gleich der Zeilendimension von B ist oder wenn eine der beiden ein Skalar ist. Wenn A gleich m-mal-p und B gleich p-mal-n ist, ist ihr Produkt C gleich m-mal-n. Das Produkt kann mithilfe von MATLAB for-Schleifen, colon-Schreibweise und Vektor-Punkt-Produkten definiert werden:

A = pascal(3);
B = magic(3);
m = 3; 
n = 3;
for i = 1:m
     for j = 1:n
        C(i,j) = A(i,:)*B(:,j);
     end
end

MATLAB verwendet ein Sternchen, um die Matrixmultiplikation anzugeben, wie zum Beispiel in C = A*B. Die Matrixmultiplikation ist nicht kommutativ, das heißt, A*B ist typischerweise nicht gleich B*A:

X = A*B
X =
      15    15    15
      26    38    26
      41    70    39
Y = B*A
Y =
      15    28    47
      15    34    60
      15    28    43

Eine Matrix kann auf der rechten Seite mit einem Spaltenvektor und auf der linken Seite mit einem Zeilenvektor multipliziert werden:

u = [3; 1; 4];
x = A*u
x =

     8
    17
    30
v = [2 0 -1];
y = v*B
y =

    12    -7    10

Rechteckige Matrixmultiplikationen müssen die Bedingungen der Dimensionskompatibilität erfüllen. Da A 3x3 und C 3x2 ist, können Sie diese multiplizieren, um ein 3x2-Ergebnis zu erhalten (die allgemeine innere Dimension wird aufgehoben):

X = A*C
X =

    24    17
    47    42
    79    77

Allerdings funktioniert die Multiplikation nicht in umgekehrter Reihenfolge:

Y = C*A
Error using  * 
Incorrect dimensions for matrix multiplication. Check that the number of columns 
in the first matrix matches the number of rows in the second matrix. To perform 
elementwise multiplication, use '.*'.

Sie können alles mit einem Skalar multiplizieren:

s = 10;
w = s*y
w =

   120   -70   100

Wenn Sie ein Array mit einem Skalar multiplizieren, expandiert der Skalar implizit auf die Größe der anderen Eingabe. Dies wird häufig auch Skalarexpansion genannt.

Identitätsmatrix

In der allgemein akzeptierten mathematischen Schreibweise wird der Großbuchstabe I verwendet, um Identitätsmatrizen, Matrizen unterschiedlicher Größen mit Einsen in der Hauptdiagonalen und Nullen an anderer Stelle, anzugeben. Diese Matrizen haben die Eigenschaft, dass AI = A und IA = A, sofern die Dimensionen kompatibel sind.

Die ursprüngliche Version von MATLAB konnte I zu diesem Zweck nicht verwenden, da nicht zwischen Groß- und Kleinbuchstaben unterschieden wurde und i bereits als Indexzeichen und als komplexe Einheit diente. Daher wurde ein englisches Wortspiel eingeführt. Die Funktion

eye(m,n)

gibt eine rechteckige m-mal-n-Identitätsmatrix zurück und eye(n) gibt eine quadratische n-mal-n-Identitätsmatrix zurück.

Matrixinverse

Wenn eine Matrix A quadratisch oder nichtsingulär (Determinante ungleich null) ist, ergeben die Gleichungen AX = I und XA = I dieselbe Lösung X. Diese Lösung wird die Inverse von A genannt und wie folgt angegeben: A-1. Die Funktion inv und der Ausdruck A^-1 berechnen beide die Matrixinverse.

A = pascal(3)
A =
       1     1     1
       1     2     3
       1     3     6
X = inv(A)
X =

    3.0000   -3.0000    1.0000
   -3.0000    5.0000   -2.0000
    1.0000   -2.0000    1.0000
A*X
ans =

    1.0000         0         0
    0.0000    1.0000   -0.0000
   -0.0000    0.0000    1.0000

Die Determinante, die von det berechnet wird, ist ein Maß des Skalierungsfaktors der linearen Transformation, die von der Matrix beschrieben wird. Ist die Determinante exakt null, gilt die Matrix als singulär und es gibt keine Inverse.

d = det(A)
d =

     1

Einige Matrizen sind nahezu singulär und trotz der Tatsache, dass eine inverse Matrix vorliegt, ist die Berechnung anfällig für numerische Fehler. Die Funktion cond berechnet die Konditionszahl für die Inversion, die einen Anhaltspunkt für die Genauigkeit der Ergebnisse aus der Matrixinversion gibt. Die Konditionszahl liegt zwischen 1 für eine numerisch stabile Matrix und Inf für eine singuläre Matrix.

c = cond(A)
c =

   61.9839

Es ist nur selten erforderlich, die explizite Inverse einer Matrix zu bilden. Häufig wird inv falsch verwendet, wenn das System linearer Gleichungen gelöst wird Ax = b. Die beste Möglichkeit zur Lösung dieser Gleichung, vom Standpunkt der Ausführungszeit und der numerischen Genauigkeit aus betrachtet, ist die Verwendung des umgekehrten Schrägstrichs als Operator x = A\b für die Matrix. Weitere Informationen dazu finden Sie unter mldivide.

Kronecker-Tensor-Produkt

Das Kronecker-Produkt, kron(X,Y), zweier Matrizen ist die größere Matrix, die aus allen möglichen Produkten der Elemente von X mit denen von Y gebildet wird. Wenn X gleich m-mal-n und Y gleich p-mal-q ist, dann ist kron(X,Y) gleich mp-mal-nq. Die Elemente sind so angeordnet, dass jedes Element von X mit der gesamten Matrix Y multipliziert wird:

[X(1,1)*Y  X(1,2)*Y  . . .  X(1,n)*Y
                     . . .
 X(m,1)*Y  X(m,2)*Y  . . .  X(m,n)*Y]

Das Kronecker-Produkt wird häufig mit Matrizen aus Nullen und Einsen und mit Matrizen zum Erstellen wiederholter Kopien kleiner Matrizen verwendet. Wenn beispielsweise X die 2x2-Matrix

X = [1   2
     3   4]

und I = eye(2,2) die 2x2-Identitätsmatrix ist, gilt Folgendes:

kron(X,I)
ans =

     1     0     2     0
     0     1     0     2
     3     0     4     0
     0     3     0     4

und

kron(I,X)
ans =

     1     2     0     0
     3     4     0     0
     0     0     1     2
     0     0     3     4

Neben kron gibt es auch noch weitere Funktionen, die nützlich sind, um Arrays zu replizieren, beispielsweise repmat, repelem und blkdiag.

Vektor- und Matrixnormen

Die p-Norm eines Vektors x,

xp=(|xi|p)1p,

wird durch norm(x,p) berechnet. Diese Rechenoperation ist für einen beliebigen Wert von p > 1 definiert, doch die gängigsten Werte von p sind 1, 2 und ∞. Der Standardwert ist p = 2, was der euklidischen Länge oder dem Betrag des Vektors entspricht:

v = [2 0 -1];
[norm(v,1) norm(v) norm(v,inf)]
ans =

    3.0000    2.2361    2.0000

Die p-Norm einer Matrix A,

Ap=maxxAxpxp,

kann für p = 1, 2 und ∞ durch norm(A,p) berechnet werden. Und wieder ist der Standardwert p = 2:

A = pascal(3);
[norm(A,1) norm(A) norm(A,inf)]
ans =

   10.0000    7.8730   10.0000

Wenn Sie die Norm jeder Zeile oder Spalte einer Matrix berechnen möchten, können Sie vecnorm verwenden:

vecnorm(A)
ans =

    1.7321    3.7417    6.7823

Verwenden der Multithread-Berechnung mit Funktionen der linearen Algebra

MATLAB unterstützt Multithread-Berechnung für zahlreiche numerische Funktionen der linearen Algebra und elementweiser numerischer Funktionen. Diese Funktionen werden automatisch für mehrere Threads ausgeführt. Damit eine Funktion oder ein Ausdruck schneller auf mehreren CPUs ausgeführt wird, müssen verschiedene Bedingungen erfüllt sein:

  1. Die Funktion führt Operationen aus, die sich leicht in gleichzeitig ausführbare Abschnitte unterteilen lassen. Diese Abschnitte müssen in der Lage sein, mit möglichst wenig Kommunikation zwischen Prozessen ausgeführt zu werden. Sie sollten nur wenige aufeinanderfolgende Operationen erfordern.

  2. Die Datengröße muss so bemessen sein, dass alle Vorteile einer gleichzeitigen Ausführung die Zeit aufwiegen, die zum Unterteilen der Daten und zum Verwalten separater Ausführungsthreads erforderlich ist. So ist eine Beschleunigung der meisten Funktionen beispielsweise nur erforderlich, wenn das Array mehrere tausend Elemente oder mehr enthält.

  3. Die Rechenoperation ist nicht speichergebunden und die Verarbeitungszeit wird nicht von der Speicherzugriffszeit beeinträchtigt. Im Allgemeinen gilt, dass komplizierte Funktionen eher beschleunigt werden als einfache Funktionen.

Die Operatoren für Matrixmultiplikation (X*Y) und Matrixpotenz (X^p) werden bei großen Arrays mit doppelter Präzision (in der Größenordnung von 10.000 Elementen) wesentlich schneller. Die Matrixanalysefunktionen det, rcond, hess und expm sind bei großen Arrays mit doppelter Präzision wesentlich schneller.

Verwandte Themen

Systeme linearer Gleichungen

Rechnerische Überlegungen

Eines der wichtigsten Probleme bei der technischen Informatik ist die Lösung von Systemen simultaner linearer Gleichungen.

In der Matrixschreibweise sieht das allgemeine Problem wie folgt aus: Gibt es, wenn zwei Matrizen A und b vorliegen, eine eindeutige Matrix x, sodass Axb oder xAb?

Aufschluss gibt dabei ein 1x1-Beispiel. Hat beispielsweise die Gleichung

7x = 21

eine eindeutige Lösung?

Die Antwort lautet natürlich ja. Die Gleichung hat die eindeutige Lösung x = 3. Die Lösung lässt sich leicht durch Division ermitteln:

x = 21/7 = 3.

Die Lösung wird nicht einfach durch Berechnen der Inversen von 7, also 7–1 = 0,142857..., und anschließende Multiplikation von 7–1 mal 21 berechnet. Das würde mehr Aufwand bedeuten und wäre auch weniger genau wenn 7–1 für eine finite Anzahl von Ziffern stünde. Ähnliche Bedingungen gelten für Reihen linearer Gleichungen mit mehr als einer Unbekannten. MATLAB löst solche Gleichungen ohne Berechnung der Inversen der Matrix.

Auch wenn es sich nicht um eine mathematische Standardschreibweise handelt, verwendet MATLAB die vertraute Divisionsterminologie aus dem Skalarfall, um die Lösung eines allgemeinen Systems simultaner Gleichungen zu beschreiben. Die beiden Divisionssymbole, Schrägstrich, /, und umgekehrter Schrägstrich, \, entsprechen den beiden MATLAB Funktionen mrdivide und mldivide. Diese Operatoren werden für die beiden Situationen verwendet, in denen die unbekannte Matrix auf der linken oder rechten Seite der Koeffizientenmatrix angezeigt wird:

x = b/A

Gibt die Lösung für die Matrixgleichung xA = b an, die mithilfe von mrdivide ermittelt wurde.

x = A\b

Gibt die Lösung für die Matrixgleichung Ax = b an, die mithilfe von mldivide ermittelt wurde.

Stellen Sie sich die „Division“ beider Seiten der Gleichung Ax = b oder xA = b durch A vor. Die Koeffizientenmatrix A ist immer im „Nenner“.

Die Dimensionskompatibilitätsbedingungen für x = A\b erfordern, dass die beiden Matrizen A und b dieselbe Anzahl von Zeilen aufweisen. Anschließend weist die Lösung x dieselbe Anzahl von Spalten auf wie b und ihre Zeilendimension entspricht der Spaltendimension von A. Für x = b/A werden die Rollen der Zeilen und Spalten vertauscht.

In der Praxis treten lineare Gleichungen der Form Ax = b häufiger auf als solche der Form xA = b. Daher wird der umgekehrte Schrägstrich wesentlich häufiger verwendet als der Schrägstrich. Im übrigen Teil dieses Abschnitts geht es um den umgekehrten Schrägstrich als Operator. Die entsprechenden Eigenschaften des Schrägstrichs als Operator können von der Identität abgeleitet werden:

(b/A)' = (A'\b').

Die Koeffizientenmatrix A muss nicht quadriert sein. Wenn A die Größe m-mal-n aufweist, gibt es drei Fälle:

m = n

Quadratisches System. Sucht nach einer exakten Lösung.

m > n

Überbestimmtes System mit mehr Gleichungen als Unbekannten. Findet eine Kleinste-Quadrate-Lösung.

m < n

Unterbestimmtes System mit weniger Gleichungen als Unbekannten. Findet eine Basislösung mit höchstens m Komponenten ungleich null.

Der mldivide-Algorithmus.  Der Operator mldivide verwendet unterschiedliche Solver für verschiedene Arten von Koeffizientenmatrizen. Die verschiedenen Fälle werden automatisch durch Untersuchen der Koeffizientenmatrix diagnostiziert. Weitere Informationen finden Sie im Abschnitt „Algorithmen“ auf der Referenzseite zu mldivide.

Allgemeine Lösung

Die allgemeine Lösung für ein System linearer Gleichungen, Axb, beschreibt alle möglichen Lösungen. Sie können die allgemeine Lösung wie folgt finden:

  1. Lösen Sie das entsprechende homogene System Ax = 0. Verwenden Sie hierfür den Befehl null, indem Sie null(A) eingeben. Dadurch wird eine Basis für den Lösungsraum an Ax = 0 zurückgegeben. Alle Lösungen sind lineare Kombinationen von Basisvektoren.

  2. Suchen Sie eine bestimmte Lösung für das nicht homogene System Ax = b.

Anschließend können Sie eine beliebige Lösung für Axb als Summe der jeweiligen Lösung für Ax = b, aus Schritt 2, plus eine lineare Kombination der Basisvektoren aus Schritt 1 schreiben.

Im übrigen Teil dieses Abschnitts wird beschrieben, wie Sie MATLAB zum Suchen einer bestimmten Lösung für Ax = b wie in Schritt 2 verwenden.

Quadratische Systeme

Die gängigste Situation umfasst eine quadratische Koeffizientenmatrix A und einen einzelnen Spaltenvektor b auf der rechten Seite.

Nichtsinguläre Koeffizientenmatrix.  Wenn die Matrix A nichtsingulär ist, hat die Lösung x = A\b dieselbe Größe wie b. Beispiel:

A = pascal(3);
u = [3; 1; 4];
x = A\u

x =
      10
     -12
       5

Es kann bestätigt werden, dass A*x exakt mit u übereinstimmt.

Wenn A und b quadratisch sind und die gleiche Größe aufweisen, hat x= A\b ebenfalls diese Größe:

b = magic(3);
X = A\b

X =
      19    -3    -1
     -17     4    13
       6     0    -6

Es kann bestätigt werden, dass A*x exakt mit b übereinstimmt.

Beide Beispiele weisen exakte, ganzzahlige Lösungen auf. Dies liegt daran, dass die Koeffizientenmatrix als pascal(3) ausgewählt wurde, was einer Matrix entspricht, die vollen Rang hat (nicht singulär).

Singuläre Koeffizientenmatrix.  Eine quadratische Matrix A ist singulär, wenn sie keine linear unabhängigen Spalten aufweist. Wenn A singulär ist, ist die Lösung für Ax = b entweder nicht vorhanden oder nicht eindeutig. Der umgekehrte Schrägstrich als Operator, A\b, gibt eine Warnung aus, wenn A nahezu singulär ist oder wenn eine exakte Singularität erkannt wird.

Wenn A singulär ist und Ax = b eine Lösung hat, können Sie eine bestimmte, nicht eindeutige Lösung finden, indem Sie Folgendes eingeben:

P = pinv(A)*b

pinv(A) ist eine Pseudoinverse von A. Wenn Ax = b keine exakte Lösung aufweist, gibt pinv(A) eine Kleinste-Quadrate-Lösung zurück.

Beispiel:

A = [ 1     3     7
     -1     4     4
      1    10    18 ]

ist singulär, da eine Verifizierung möglich ist durch Eingabe von

rank(A)

ans =

     2

Da A keinen vollen Rang hat, weist sie einige singuläre Werte gleich null auf.

Exakte Lösungen. Für b =[5;2;12] weist die Gleichung Ax = b eine exakte Lösung auf, die wie folgt angegeben wird:

pinv(A)*b

ans =
    0.3850
   -0.1103
    0.7066

Verifizieren Sie, dass pinv(A)*b eine exakte Lösung ist, indem Sie Folgendes eingeben:

A*pinv(A)*b

ans =
    5.0000
    2.0000
   12.0000

Kleinste-Quadrate-Lösungen. Wenn jedoch b = [3;6;0], hat Ax = b keine exakte Lösung. In diesem Fall gibt pinv(A)*b eine Kleinste-Quadrate-Lösung zurück. Wenn Sie

A*pinv(A)*b

ans =
   -1.0000
    4.0000
    2.0000

eingeben, wird nicht der ursprüngliche Vektor b zurückgegeben.

Sie können bestimmen, ob Ax = b eine exakte Lösung aufweist, indem Sie nach der zeilenreduzierten Treppenform der erweiterten Matrix [A b] suchen. Geben Sie dazu für dieses Beispiel Folgendes ein:

rref([A b])
ans =
    1.0000         0    2.2857         0
         0    1.0000    1.5714         0
         0         0         0    1.0000

Da die unterste Zeile ausschließlich Nullen enthält, mit Ausnahme des letzten Eintrags, gibt es für die Gleichung keine Lösung. In diesem Fall gibt pinv(A) eine Kleinste-Quadrate-Lösung zurück.

Überbestimmte Systeme

Dieses Beispiel veranschaulicht, wie überbestimmte Systeme oft in den unterschiedlichsten Arten der Kurvenanpassung für experimentelle Daten vorkommen.

Eine Größe y wird an vielen verschiedenen Werten der Zeit t gemessen und führt zu den folgenden Beobachtungen. Sie können die Daten mithilfe der folgenden Anweisungen eingeben und in einer Tabelle anzeigen.

t = [0 .3 .8 1.1 1.6 2.3]';
y = [.82 .72 .63 .60 .55 .50]';
B = table(t,y)
B=6×2 table
     t      y  
    ___    ____

      0    0.82
    0.3    0.72
    0.8    0.63
    1.1     0.6
    1.6    0.55
    2.3     0.5

Versuchen Sie, die Daten mit einer abklingenden Exponentialfunktion zu modellieren

y(t)=c1+c2e-t.

Die vorherige Gleichung zeigt, dass der Vektor y durch eine lineare Kombination zweier anderer Vektoren approximiert werden sollte. Einer der Vektoren ist ein konstanter Vektor, der nur Einsen enthält, und der andere ist der Vektor mit den Komponenten exp(-t). Die unbekannten Koeffizienten, c1 und c2, können mit der Kleinste-Quadrate-Methode berechnet werden, bei der die Summe der Abweichungsquadrate der Daten aus dem Modell minimiert wird. Es gibt sechs Gleichungen in zwei Unbekannten, die durch eine 6x2-Matrix dargestellt werden.

E = [ones(size(t)) exp(-t)]
E = 6×2

    1.0000    1.0000
    1.0000    0.7408
    1.0000    0.4493
    1.0000    0.3329
    1.0000    0.2019
    1.0000    0.1003

Verwenden Sie den umgekehrten Schrägstrich als Operator, um die Kleinste-Quadrate-Lösung zu erhalten.

c = E\y
c = 2×1

    0.4760
    0.3413

Anders ausgedrückt, lautet die Kleinste-Quadrate-Methode für die Daten wie folgt:

y(t)=0.4760+0.3413e-t.

Die folgenden Anweisungen werten das Modell in Inkrementen mit regelmäßigen Abständen in t aus und plotten anschließend das Ergebnis zusammen mit den ursprünglichen Daten:

T = (0:0.1:2.5)';
Y = [ones(size(T)) exp(-T)]*c;
plot(T,Y,'-',t,y,'o')

Figure contains an axes object. The axes object contains 2 objects of type line.

E*c ist nicht exakt gleich y, doch die Differenz kann durchaus kleiner sein als Messfehler in den ursprünglichen Daten.

Eine rechteckige Matrix A hat einen unzureichenden Rang, wenn sie keine linear unabhängigen Spalten aufweist. Wenn A einen unzureichenden Rang hat, ist die Kleinste-Quadrate-Lösung für AX = B nicht eindeutig. A\B gibt eine Warnung aus, wenn A einen unzureichenden Rang hat, und erstellt eine Kleinste-Quadrate-Lösung. Sie können mithilfe von lsqminnorm die Lösung X finden, die die minimale Norm unter allen Lösungen aufweist.

Unterbestimmte Systeme

Dieses Beispiel veranschaulicht, wie die Lösung für unterbestimmte Systeme nicht eindeutig ist. Unterbestimmte lineare Systeme beinhalten weitaus mehr Unbekannte als Gleichungen. Die linke Divisionsoperation der Matrix in MATLAB findet eine grundlegende Kleinste-Quadrate-Lösung, die höchstens m Komponenten, die nicht null sind, für eine m-mal-n-Koeffizientenmatrix besitzt.

Im Folgenden ist ein kleines, willkürliches Beispiel abgebildet:

R = [6 8 7 3; 3 5 4 1]
rng(0);
b = randi(8,2,1)
R =

       6              8              7              3       
       3              5              4              1       


b =

       7       
       8      

Das lineare System Rp = b besitzt zwei Gleichungen in vier Unbekannten. Da die Koeffizientenmatrix kleine Ganzzahlen enthält, sollte der Befehl format die Lösung im rationalen Format anzeigen. Die jeweilige Lösung wird wie folgt ermittelt:

format rat
p = R\b
p =

       0       
      17/7     
       0       
     -29/7    

Eine der Komponenten, die nicht null sind, ist p(2), da R(:,2) die Spalte von R mit der größten Norm ist. Die andere Komponente, die nicht gleich null ist, ist p(4), da R(:,4) dominiert, nachdem R(:,2) eliminiert wurde.

Die vollständige allgemeine Lösung für das unterbestimmte System kann durch Hinzufügen von p zu einer beliebigen linearen Kombination der Nullraumvektoren charakterisiert werden, die mithilfe der Funktion null und einer Option zum Anfordern einer rationalen Basis gefunden werden kann.

Z = null(R,'r')
Z =

      -1/2           -7/6     
      -1/2            1/2     
       1              0       
       0              1       

Es kann bestätigt werden, dass R*Z null ist und dass das Residuum R*x - b für einen beliebigen Vektor x klein ist, wobei Folgendes gilt:

x = p + Z*q

Da die Spalten von Z die Nullraumvektoren sind, ist das Produkt Z*q eine lineare Kombination aus den folgenden Vektoren:

Zq=(x1x2)(uw)=ux1+wx2.

Wählen Sie zur Veranschaulichung eine beliebige Komponente q aus und konstruieren Sie x.

q = [-2; 1];
x = p + Z*q;

Berechnen Sie die Norm des Residuums.

format short
norm(R*x - b)
ans =

   2.6645e-15

Wenn unendlich viele Lösungen verfügbar sind, ist die Lösung mit der minimalen Norm von besonderem Interesse. Sie können lsqminnorm zum Berechnen der Kleinste-Quadrate-Lösung mit der minimalen Norm verwenden. Diese Lösung weist den kleinstmöglichen Wert für norm(p) auf.

p = lsqminnorm(R,b)
p =

    -207/137   
     365/137   
      79/137   
    -424/137  

Lösungen für verschiedene rechte Seiten

Einige Probleme betreffen die Lösung linearer Systeme, die dieselbe Koeffizientenmatrix A, doch unterschiedliche rechte Seiten b aufweisen. Wenn die verschiedenen Werte von b gleichzeitig verfügbar sind, können Sie b als Matrix mit mehreren Spalten konstruieren und alle Gleichungssysteme gleichzeitig mithilfe eines einzelnen Befehls mit umgekehrtem Schrägstrich lösen: X = A\[b1 b2 b3 …].

Manchmal sind jedoch die verschiedenen Werte von b nicht alle gleichzeitig verfügbar, sodass Sie verschiedene Gleichungssysteme nacheinander lösen müssen. Wenn Sie eines dieser Gleichungssysteme mit Schrägstrich (/) oder umgekehrtem Schrägstrich (\) lösen, faktorisiert der Operator die Koeffizientenmatrix A und verwendet diese Matrixzerlegung, um die Lösung zu berechnen. Allerdings berechnet der Operator jedes Mal, wenn Sie ein ähnliches Gleichungssystem mit einem anderen Wert von b lösen, dieselbe Zerlegung von A, bei der es sich um eine redundante Berechnung handelt.

Dieses Problem kann behoben werden, indem die Zerlegung von A vorab berechnet wird und anschließend die zu lösenden Faktoren für die verschiedenen Werte von b wiederverwendet werden. In der Praxis kann jedoch die Vorabberechnung der Zerlegung auf diese Weise schwierig sein, da Sie wissen müssen, welche Zerlegung zu berechnen ist (LU, LDL, Cholesky usw.) und wie Sie die Faktoren zum Lösen des Problems multiplizieren. Mit der LU-Zerlegung müssen Sie zwei lineare Systeme lösen, um das ursprüngliche System Ax = b zu lösen:

[L,U] = lu(A);
x = U \ (L \ b);

Stattdessen wird zum Lösen linearer Systeme mit verschiedenen aufeinanderfolgenden rechten Seiten die Verwendung von Zerlegungsobjekten (decomposition) empfohlen. Mit diesen Objekten können die Leistungsvorteile einer Vorabberechnung der Matrixzerlegung genutzt werden, doch es ist kein Wissen über die Verwendung der Matrixfaktoren erforderlich. Sie können die vorherige LU-Zerlegung wie folgt ersetzen:

dA = decomposition(A,'lu');
x = dA\b;

Wenn Sie nicht sicher sind, welche Zerlegung verwendet werden soll, wählt decomposition(A) den richtigen Typ basierend auf den Eigenschaften von A aus, ähnlich wie der umgekehrte Schrägstrich.

Hier ist ein einfacher Test der möglichen Leistungsvorteile dieses Ansatzes. Der Test löst dasselbe dünn besetzte lineare System 100 Mal mithilfe des umgekehrten Schrägstrichs (\) und decomposition.

n = 1e3;
A = sprand(n,n,0.2) + speye(n);
b = ones(n,1);

% Backslash solution
tic
for k = 1:100
    x = A\b;
end
toc
Elapsed time is 9.006156 seconds.
% decomposition solution
tic
dA = decomposition(A);
for k = 1:100
    x = dA\b;
end
toc
Elapsed time is 0.374347 seconds.

Bei diesem Problem ist die decomposition-Lösung wesentlich schneller als die alleinige Verwendung des umgekehrten Schrägstrichs – und doch bleibt die Syntax einfach.

Iterative Methoden

Wenn die Koeffizientenmatrix A groß und dünn besetzt ist, sind Faktorisierungsmethoden im Allgemeinen nicht effizient. Iterative Methoden generieren eine Reihe von Approximationslösungen. MATLAB stellt verschiedene iterative Methoden zur Handhabung großer, dünn besetzter Eingabematrizen zur Verfügung.

FunktionBeschreibung
pcg

Methode der vorausgesetzten konjugierten Gradienten. Diese Methode eignet sich für die hermitesche positiv definite Koeffizientenmatrix A.

bicg

Methode bikonjugierter Gradienten

bicgstab

Stabilisierte Methode bikonjugierter Gradienten

bicgstabl

BiCGStab(l)-Methode

cgs

Quadrierte Methode konjugierter Gradienten

gmres

Methode generalisierter minimaler Residuen

lsqr

LSQR-Methode

minres

Methode minimaler Residuen. Diese Methode eignet sich für die hermitesche Koeffizientenmatrix A.

qmr

Methode quasi-minimaler Residuen

symmlq

Symmetrische LQ-Methode

tfqmr

Transponiertenfreie QMR-Methode

Multithread-Berechnung

MATLAB unterstützt Multithread-Berechnung für zahlreiche numerische Funktionen der linearen Algebra und elementweiser numerischer Funktionen. Diese Funktionen werden automatisch für mehrere Threads ausgeführt. Damit eine Funktion oder ein Ausdruck schneller auf mehreren CPUs ausgeführt wird, müssen verschiedene Bedingungen erfüllt sein:

  1. Die Funktion führt Operationen aus, die sich leicht in gleichzeitig ausführbare Abschnitte unterteilen lassen. Diese Abschnitte müssen in der Lage sein, mit möglichst wenig Kommunikation zwischen Prozessen ausgeführt zu werden. Sie sollten nur wenige aufeinanderfolgende Operationen erfordern.

  2. Die Datengröße muss so bemessen sein, dass alle Vorteile einer gleichzeitigen Ausführung die Zeit aufwiegen, die zum Unterteilen der Daten und zum Verwalten separater Ausführungsthreads erforderlich ist. So ist eine Beschleunigung der meisten Funktionen beispielsweise nur erforderlich, wenn das Array mehrere tausend Elemente oder mehr enthält.

  3. Die Rechenoperation ist nicht speichergebunden und die Verarbeitungszeit wird nicht von der Speicherzugriffszeit beeinträchtigt. Im Allgemeinen gilt, dass komplizierte Funktionen eher beschleunigt werden als einfache Funktionen.

inv, lscov, linsolve und mldivide sind bei großen Arrays mit doppelter Präzision (in der Größenordnung von 10.000 Elementen oder mehr) erheblich schneller, wenn Multithreading aktiviert ist.

Siehe auch

| | | |

Verwandte Themen

Faktorisierungen

Einführung

Alle drei in diesem Abschnitt erläuterten Matrixfaktorisierungen verwenden Dreiecksmatrizen, bei denen alle Elemente entweder über oder unter der Diagonalen null sind. Systeme linearer Gleichungen, die Dreiecksmatrizen enthalten, werden mithilfe der Vorwärts- oder Rückwärtssubstitution schnell und einfach gelöst.

Cholesky-Faktorisierung

Die Cholesky-Faktorisierung drückt eine symmetrische Matrix als Produkt einer Dreiecksmatrix und ihrer Transponierten aus

A = RR,

wobei R eine obere Dreiecksmatrix ist.

Nicht alle symmetrischen Matrizen können auf diese Weise faktorisiert werden. Die Matrizen mit einer solchen Faktorisierung gelten als positiv definit. Dies impliziert, dass alle Diagonalelemente von A positiv sind und dass die Elemente, die nicht auf der Diagonalen liegen, „nicht zu groß“ sind. Die Pascal-Matrizen sind ein interessantes Beispiel. In diesem Kapitel war die Beispielmatrix A die 3-mal-3-Pascal-Matrix. Stellen Sie vorübergehend auf 6-mal-6 um:

A = pascal(6)

A =
       1     1     1     1     1     1
       1     2     3     4     5     6
       1     3     6    10    15    21
       1     4    10    20    35    56
       1     5    15    35    70   126
       1     6    21    56   126   252

Die Elemente von A sind Binomialkoeffizienten. Jedes Element ist die Summe seiner nördlichen und westlichen Nachbarn. Die Cholesky-Faktorisierung sieht wie folgt aus:

R = chol(A)

R =
     1     1     1     1     1     1
     0     1     2     3     4     5
     0     0     1     3     6    10
     0     0     0     1     4    10
     0     0     0     0     1     5
     0     0     0     0     0     1

Die Elemente sind auch hier Binomialkoeffizienten. Die Tatsache, dass R'*R gleich A ist, weist auf eine Identität hin, die Summen von Produkten der Binomialkoeffizienten umfasst.

Hinweis

Die Cholesky-Faktorisierung gilt auch für komplexe Matrizen. Alle komplexen Matrizen, die eine Cholesky-Faktorisierung aufweisen, erfüllen

A′ = A

und gelten als hermitesch positiv definit.

Die Cholesky-Faktorisierung erlaubt die Ersetzung des linearen Systems

Ax = b

durch

RRx = b.

Da der umgekehrte Schrägstrich als Operator Dreiecksysteme erkennt, kann dies in der MATLAB Umgebung schnell wie folgt gelöst werden:

x = R\(R'\b)

Wenn A gleich n-mal-n ist, ist die rechnerische Komplexität von chol(A) gleich O(n3), doch die Komplexität der nachfolgenden Lösungen mit umgekehrtem Schrägstrich ist nur O(n2).

LU-Faktorisierung

Die LU-Faktorisierung, oder das Gaußsche Eliminationsverfahren, drückt eine beliebige quadratische Matrix A als Produkt der Permutation einer unteren Dreiecksmatrix und einer oberen Dreiecksmatrix aus,

A = LU,

wobei L eine Permutation der unteren Dreiecksmatrix mit Einsen entlang der Diagonalen und U eine obere Dreiecksmatrix ist.

Die Permutationen sind aus theoretischen und rechnerischen Gründen erforderlich. Die Matrix

[0110]

kann nicht als Produkt von Dreiecksmatrizen ausgedrückt werden, wenn nicht die beiden Zeilen vertauscht werden. Auch wenn die Matrix

[ε110]

als Produkt der Dreiecksmatrizen ausgedrückt werden kann, wenn ε klein ist, sind die Elemente in den Faktoren groß und vergrößern Fehler, sodass die Permutationen zwar nicht unbedingt erforderlich, doch wünschenswert sind. Durch teilweise Pivotisierung kann sichergestellt werden, dass die Größe der Elemente von L durch 1 begrenzt ist und dass die Elemente von U nicht wesentlich größer sind als die von A.

Beispiel:

[L,U] = lu(B)

L =
    1.0000         0         0
    0.3750    0.5441    1.0000
    0.5000    1.0000         0

U =
    8.0000    1.0000    6.0000
         0    8.5000   -1.0000
         0         0    5.2941

Durch die LU-Faktorisierung von A kann das lineare System

A*x = b

wie folgt schnell gelöst werden:

x = U\(L\b)

Determinanten und Inverse lassen sich mit der LU-Faktorisierung wie folgt berechnen:

det(A) = det(L)*det(U)

und

inv(A) = inv(U)*inv(L)

Sie können die Determinanten auch mithilfe von det(A) = prod(diag(U)) berechnen, wobei jedoch die Vorzeichen der Determinanten umgekehrt sein können.

QR-Faktorisierungen

Eine orthogonale Matrix, oder eine Matrix mit orthonormalen Spalten, ist eine reelle Matrix, deren Spalten alle eine Einheitenlänge aufweisen und senkrecht zueinander sind. Wenn Q orthogonal ist, dann gilt Folgendes:

QTQ = I,

wobei I die Identitätsmatrix ist.

Die einfachsten orthogonalen Matrizen sind zweidimensionale Koordinatendrehungen:

[cos(θ)sin(θ)sin(θ)cos(θ)].

Für komplexe Matrizen lautet der entsprechende Begriff unitär. Orthogonale und unitäre Matrizen sind für die numerische Berechnung wünschenswert, da sie Länge und Winkel beibehalten und Fehler nicht vergrößern.

Die orthogonale Faktorisierung, oder QR-Faktorisierung, drückt eine beliebige rechteckige Matrix als Produkt einer orthogonalen oder unitären Matrix und einer oberen Dreiecksmatrix aus. Eine Spaltenpermutation kann ebenfalls beteiligt sein:

A = QR

oder

AP = QR,

wobei Q die orthogonale oder unitäre Matrix, R die obere Dreiecksmatrix und P eine Permutation ist.

Es gibt vier Varianten der QR-Faktorisierung – vollständig oder reduziert sowie mit und ohne Spaltenpermutation.

Überbestimmte lineare Systeme umfassen eine rechteckige Matrix mit mehr Zeilen als Spalten, also m-mal-n mit m > n. Bei der vollständigen QR-Faktorisierung ergibt sich eine quadratische, orthogonale m-mal-m-Q-Matrix und eine rechteckige, obere m-mal-n-R-Dreiecksmatrix:

C=gallery('uniformdata',[5 4], 0);
[Q,R] = qr(C)

Q =

    0.6191    0.1406   -0.1899   -0.5058    0.5522
    0.1506    0.4084    0.5034    0.5974    0.4475
    0.3954   -0.5564    0.6869   -0.1478   -0.2008
    0.3167    0.6676    0.1351   -0.1729   -0.6370
    0.5808   -0.2410   -0.4695    0.5792   -0.2207


R =

    1.5346    1.0663    1.2010    1.4036
         0    0.7245    0.3474   -0.0126
         0         0    0.9320    0.6596
         0         0         0    0.6648
         0         0         0         0

In vielen Fällen sind die letzten m – n-Spalten von Q nicht erforderlich, da sie mit den Nullen im unteren Teil von R multipliziert werden. Bei der reduzierten QR-Faktorisierung entsteht eine rechteckige m-mal-n-Q-Matrix mit orthonormalen Spalten und eine quadratische, obere n-mal-n-R-Dreiecksmatrix. Für das 5x4-Beispiel stellt dies keine besondere Einsparung dar, doch für größere, stark rechteckige Matrizen können die Einsparungen hinsichtlich Zeit und Speicher tatsächlich von Bedeutung sein:

[Q,R] = qr(C,0)	
Q =

    0.6191    0.1406   -0.1899   -0.5058
    0.1506    0.4084    0.5034    0.5974
    0.3954   -0.5564    0.6869   -0.1478
    0.3167    0.6676    0.1351   -0.1729
    0.5808   -0.2410   -0.4695    0.5792


R =

    1.5346    1.0663    1.2010    1.4036
         0    0.7245    0.3474   -0.0126
         0         0    0.9320    0.6596
         0         0         0    0.6648

Im Gegensatz zur LU-Faktorisierung sind für die QR-Faktorisierung weder Pivotisierung noch Permutationen erforderlich. Doch eine optionale Spaltenpermutation, die durch die Präsenz eines dritten Ausgabearguments ausgelöst wird, ist hilfreich, um Singularität oder einen nicht ausreichenden Rang zu erkennen. Bei jedem Schritt der Faktorisierung wird die Spalte der verbleibenden, nicht faktorisierten Matrix mit der größten Norm als Basis für den jeweiligen Schritt verwendet. So ist sichergestellt, dass die Diagonalelemente von R in absteigender Reihenfolge auftreten und eine lineare Abhängigkeit entlang der Spalten nahezu sicher bei der Überprüfung dieser Elemente erkannt wird. Im folgenden kleinen Beispiel weist die zweite Spalte von C eine größere Norm auf als die erste, sodass die beiden Spalten vertauscht werden:

[Q,R,P] = qr(C)

Q =
   -0.3522    0.8398   -0.4131
   -0.7044   -0.5285   -0.4739
   -0.6163    0.1241    0.7777

R =
  -11.3578   -8.2762
         0    7.2460
         0         0

P =
     0     1
     1     0

Wenn die reduzierte und die Spaltenpermutationen kombiniert werden, ist das dritte Ausgabeargument ein Permutationsvektor und keine Permutationsmatrix:

[Q,R,p] = qr(C,0)

Q =
   -0.3522    0.8398
   -0.7044   -0.5285
   -0.6163    0.1241

R =
  -11.3578   -8.2762
         0    7.2460


p =
     2     1

Die QR-Faktorisierung transformiert ein überbestimmtes lineares System in ein äquivalentes Dreieckssystem. Der Ausdruck

norm(A*x - b)

ist gleich

norm(Q*R*x - b)

Durch Multiplikation mit orthogonalen Matrizen bleibt die euklidische Norm beibehalten, sodass dieser Ausdruck gleich

norm(R*x - y)

ist, wobei y = Q'*b. Da die letzten m-n-Zeilen von R gleich null sind, besteht dieser Ausdruck aus zwei Teilen:

norm(R(1:n,1:n)*x - y(1:n))

und

norm(y(n+1:m))

Wenn A den vollen Rang aufweist, können Sie nach x auflösen, sodass der erste dieser Ausdrücke null ist. Dann gibt der zweite Ausdruck die Norm des Residuums an. Wenn A nicht den vollen Rang aufweist, kann mit der Dreiecksstruktur von R eine Basislösung für das Kleinste-Quadrate-Problem gefunden werden.

Verwenden der Multithread-Berechnung für die Faktorisierung

Die MATLAB Software unterstützt Multithread-Berechnung für zahlreiche numerische Funktionen der linearen Algebra und elementweiser numerischer Funktionen. Diese Funktionen werden automatisch für mehrere Threads ausgeführt. Damit eine Funktion oder ein Ausdruck schneller auf mehreren CPUs ausgeführt wird, müssen verschiedene Bedingungen erfüllt sein:

  1. Die Funktion führt Operationen aus, die sich leicht in gleichzeitig ausführbare Abschnitte unterteilen lassen. Diese Abschnitte müssen in der Lage sein, mit möglichst wenig Kommunikation zwischen Prozessen ausgeführt zu werden. Sie sollten nur wenige aufeinanderfolgende Operationen erfordern.

  2. Die Datengröße muss so bemessen sein, dass alle Vorteile einer gleichzeitigen Ausführung die Zeit aufwiegen, die zum Unterteilen der Daten und zum Verwalten separater Ausführungsthreads erforderlich ist. So ist eine Beschleunigung der meisten Funktionen beispielsweise nur erforderlich, wenn das Array mehrere tausend Elemente oder mehr enthält.

  3. Die Rechenoperation ist nicht speichergebunden und die Verarbeitungszeit wird nicht von der Speicherzugriffszeit beeinträchtigt. Im Allgemeinen gilt, dass komplizierte Funktionen eher beschleunigt werden als einfache Funktionen.

lu und qr sind bei großen Arrays mit doppelter Präzision (in der Größenordnung von 10.000 Elementen) erheblich schneller.

Siehe auch

| |

Verwandte Themen

Potenzen und Exponentialfunktionen

In diesem Abschnitt wird die Berechnung von Matrixpotenzen und Exponentialfunktionen mithilfe unterschiedlicher Methoden veranschaulicht.

Positive ganzzahlige Potenzen

Wenn A eine quadratische Matrix und p eine positive Ganzzahl ist, multipliziert A^p effektiv A p-1 Mal mit sich selbst. Beispiel:

A = [1 1 1
     1 2 3
     1 3 6];
A^2
ans = 3×3

     3     6    10
     6    14    25
    10    25    46

Inverse und Potenzen mit rationalem Exponent

Wenn A quadratisch und nicht singulär ist, multipliziert A^(-p) effektiv inv(A) p-1 Mal mit sich selbst.

A^(-3)
ans = 3×3

  145.0000 -207.0000   81.0000
 -207.0000  298.0000 -117.0000
   81.0000 -117.0000   46.0000

MATLAB® berechnet inv(A) und A^(-1) mit demselben Algorithmus, sodass die Ergebnisse exakt identisch sind. Sowohl inv(A) als auch A^(-1) geben Warnungen aus, wenn die Matrix nahezu singulär ist.

isequal(inv(A),A^(-1))
ans = logical
   1

Potenzen mit rationalem Exponent wie A^(2/3) sind ebenfalls zulässig. Ergebnisse bei Verwendung von Potenzen mit rationalem Exponent hängen von der Verteilung der Eigenwerte der Matrix ab.

A^(2/3)
ans = 3×3

    0.8901    0.5882    0.3684
    0.5882    1.2035    1.3799
    0.3684    1.3799    3.1167

Elementweise Potenzen

Der Operator .^ berechnet elementweise Potenzen. Um beispielsweise jedes Element in einer Matrix zu quadrieren, können Sie A.^2 verwenden.

A.^2
ans = 3×3

     1     1     1
     1     4     9
     1     9    36

Quadratwurzeln

Die Funktion sqrt ist eine komfortable Möglichkeit, die Quadratwurzel jedes Elements in einer Matrix zu berechnen. Alternativ kann auch A.^(1/2) verwendet werden.

sqrt(A)
ans = 3×3

    1.0000    1.0000    1.0000
    1.0000    1.4142    1.7321
    1.0000    1.7321    2.4495

Für andere Wurzeln können Sie nthroot verwenden. Berechnen Sie beispielsweise A.^(1/3).

nthroot(A,3)
ans = 3×3

    1.0000    1.0000    1.0000
    1.0000    1.2599    1.4422
    1.0000    1.4422    1.8171

Diese elementweisen Wurzeln unterscheiden sich von der Quadratwurzel einer Matrix, die eine zweite Matrix B so berechnet, dass A=BB. Die Funktion sqrtm(A) berechnet A^(1/2) mit einem genaueren Algorithmus. Das m in sqrtm unterscheidet diese Funktion von der Funktion sqrt(A), die wie A.^(1/2) elementweise ausgeführt wird.

B = sqrtm(A)
B = 3×3

    0.8775    0.4387    0.1937
    0.4387    1.0099    0.8874
    0.1937    0.8874    2.2749

B^2
ans = 3×3

    1.0000    1.0000    1.0000
    1.0000    2.0000    3.0000
    1.0000    3.0000    6.0000

Skalarbasen

Neben dem Erheben einer Matrix in eine Potenz können Sie auch einen Skalar in die Matrixpotenz erheben.

2^A
ans = 3×3

   10.4630   21.6602   38.5862
   21.6602   53.2807   94.6010
   38.5862   94.6010  173.7734

Wenn Sie einen Skalar in die Potenz einer Matrix erheben, verwendet MATLAB die Eigenwerte und Eigenvektoren der Matrix, um die Matrixpotenz zu berechnen. Wenn [V,D] = eig(A), dann 2A=V 2D V-1.

[V,D] = eig(A);
V*2^D*V^(-1)
ans = 3×3

   10.4630   21.6602   38.5862
   21.6602   53.2807   94.6010
   38.5862   94.6010  173.7734

Exponentialfunktionen einer Matrix

Die Exponentialfunktion einer Matrix ist ein Sonderfall beim Erheben eines Skalars in die Matrixpotenz. Die Basis für die Exponentialfunktion einer Matrix ist die Eulersche Zahl e = exp(1).

e = exp(1);
e^A
ans = 3×3
103 ×

    0.1008    0.2407    0.4368
    0.2407    0.5867    1.0654
    0.4368    1.0654    1.9418

Die Funktion expm ist eine komfortable Möglichkeit, Exponentialfunktionen einer Matrix zu berechnen.

expm(A)
ans = 3×3
103 ×

    0.1008    0.2407    0.4368
    0.2407    0.5867    1.0654
    0.4368    1.0654    1.9418

Die Exponentialfunktion einer Matrix kann auf unterschiedliche Weise berechnet werden. Weitere Informationen dazu finden Sie unter Matrix Exponentials.

Umgang mit kleinen Zahlen

Die MATLAB Funktionen log1p und expm1 berechnen log(1+x) und ex-1 für sehr kleine Werte von x genau. Wenn Sie beispielsweise versuchen, eine Zahl kleiner als die Maschinengenauigkeit zu 1 zu addieren, wird das Ergebnis auf 1 gerundet.

log(1+eps/2)
ans = 0

Allerdings kann log1p eine genauere Antwort zurückgeben.

log1p(eps/2)
ans = 1.1102e-16

Ähnliches gilt für ex-1, denn wenn x sehr klein ist, wird der Wert auf null gerundet.

exp(eps/2)-1
ans = 0

Und wiederum kann expm1 eine genauere Antwort zurückgeben.

expm1(eps/2)
ans = 1.1102e-16

Siehe auch

| | | | | | |

Verwandte Themen

Eigenwerte

Zerlegung von Eigenwerten

Eigenwert und Eigenvektor einer Quadratmatrix A sind ein Skalar λ bzw. ein Vektor υ ungleich null, die folgende Gleichung erfüllen:

= λυ.

Wenn die Eigenwerte entlang der Diagonalen einer Diagonalmatrix Λ und die entsprechenden Eigenvektoren die Spalten einer Matrix V bilden, liegt Folgendes vor:

AV = .

Wenn V nicht singulär ist, wird dies zur Eigenwertzerlegung

A = VΛV–1.

Ein gutes Beispiel ist die Koeffizientenmatrix der Differenzialgleichung dx/dt = Ax:

A =
     0    -6    -1
     6     2   -16
    -5    20   -10

Die Lösung für diese Gleichung wird mithilfe der Exponentialfunktion der Matrix ausgedrückt x(t) = etAx(0). Die Anweisung

lambda = eig(A)

ergibt einen Spaltenvektor, der die Eigenwerte von A enthält. Für diese Matrix sind die Eigenwerte komplex:

lambda =
     -3.0710         
     -2.4645+17.6008i
     -2.4645-17.6008i

Der Realteil der jeweiligen Eigenwerte ist negativ, sodass eλt gegen null geht, wenn t erhöht wird. Der Imaginärteil von zwei der Eigenwerte, ±ω, trägt die oszillatorische Komponente, sin(ωt), zur Lösung der Differenzialgleichung bei.

Mit zwei Ausgabeargumenten berechnet eig die Eigenvektoren und speichert die Eigenwerte in einer Diagonalmatrix:

[V,D] = eig(A)
V =
  -0.8326         0.2003 - 0.1394i   0.2003 + 0.1394i
  -0.3553        -0.2110 - 0.6447i  -0.2110 + 0.6447i
  -0.4248        -0.6930            -0.6930          

D =
  -3.0710                 0                 0         
        0           -2.4645+17.6008i        0         
        0                 0           -2.4645-17.6008i

Der erste Eigenvektor ist reell und die beiden anderen Vektoren sind komplexe Konjugate voneinander. Alle drei Vektoren sind normalisiert, sodass sie eine euklidische Länge, norm(v,2), von eins aufweisen.

Die Matrix V*D*inv(V), die prägnanter auch als V*D/V geschrieben werden kann, ist innerhalb des Rundungsfehlers von A. Und inv(V)*A*V, oder V\A*V, ist innerhalb des Rundungsfehlers von D.

Mehrere Eigenwerte

Einige Matrizen haben keine Eigenvektorzerlegung. Diese Matrizen sind nicht diagonalisierbar. Beispiel:

A = [ 1    -2    1 
      0     1    4 
      0     0    3 ]

Für diese Matrix ergibt

[V,D] = eig(A)

Folgendes:

V =

    1.0000    1.0000   -0.5571
         0    0.0000    0.7428
         0         0    0.3714


D =

     1     0     0
     0     1     0
     0     0     3

Es liegt ein doppelter Eigenwert bei λ = 1 vor. Die erste und zweite Spalte von V sind identisch. Für diese Matrix liegt kein vollständiger Satz linearer, unabhängiger Eigenvektoren vor.

Schur-Zerlegung

Viele erweiterte Matrixberechnungen erfordern keine Eigenwertzerlegungen. Sie basieren stattdessen auf der Schur-Zerlegung

A = USU ′ ,

wobei U eine orthogonale Matrix und S eine obere Dreiecks-Blockmatrix mit 1x1 und 2x2 Blöcken entlang der Diagonalen ist. Die Eigenwerte werden von den Diagonalelementen und Blöcken von S erkannt, während die Spalten von U eine orthogonale Basis bereitstellen, die wesentlich bessere numerische Eigenschaften hat als eine Eigenvektormenge.

Vergleichen Sie beispielsweise die Eigenwert- und die Schur-Zerlegungen dieser fehlerhaften Matrix:

A = [ 6    12    19 
     -9   -20   -33 
      4     9    15 ];

[V,D] = eig(A)
V =

  -0.4741 + 0.0000i  -0.4082 - 0.0000i  -0.4082 + 0.0000i
   0.8127 + 0.0000i   0.8165 + 0.0000i   0.8165 + 0.0000i
  -0.3386 + 0.0000i  -0.4082 + 0.0000i  -0.4082 - 0.0000i


D =

  -1.0000 + 0.0000i   0.0000 + 0.0000i   0.0000 + 0.0000i
   0.0000 + 0.0000i   1.0000 + 0.0000i   0.0000 + 0.0000i
   0.0000 + 0.0000i   0.0000 + 0.0000i   1.0000 - 0.0000i
[U,S] = schur(A)
U =

   -0.4741    0.6648    0.5774
    0.8127    0.0782    0.5774
   -0.3386   -0.7430    0.5774


S =

   -1.0000   20.7846  -44.6948
         0    1.0000   -0.6096
         0    0.0000    1.0000

Die Matrix A ist fehlerhaft, weil sie keinen vollständigen Satz linear unabhängiger Eigenvektoren hat (die zweite und dritte Spalte von V sind identisch). Da nicht alle Spalten von V linear unabhängig sind, hat sie eine hohe Bedingungszahl von etwa ~1e8. Allerdings kann schur drei verschiedene Basisvektoren in U berechnen. Da U orthogonal ist, gilt cond(U) = 1.

Die Matrix S weist den reellen Eigenwert als ersten Eintrag entlang der Diagonalen auf, während der wiederholte Eigenwert durch den 2x2-Block rechts unten dargestellt wird. Die Eigenwerte des 2x2-Blocks sind ebenfalls Eigenwerte von A:

eig(S(2:3,2:3))
ans =

   1.0000 + 0.0000i
   1.0000 - 0.0000i

Siehe auch

|

Verwandte Themen

Singuläre Werte

Ein singulärer Wert und entsprechende singuläre Vektoren einer Rechteckmatrix A sind ein Skalar σ bzw. ein Vektorpaar (u und v), die Folgendes erfüllen:

Av=σuAHu=σv,

wobei AH die hermitesche Transponierte von A ist. Die singulären Vektoren u und v werden typischerweise so skaliert, dass sie eine Norm von 1 aufweisen. Ebenfalls gilt, wenn u und v singuläre Vektoren von A sind, dann sind -u und -v auch singuläre Vektoren von A.

Die singulären Werte σ sind stets reell und nicht negativ, auch dann wenn A komplex ist. Mit den singulären Werten in einer Diagonalmatrix Σ und den entsprechenden singulären Vektoren, die die Spalten der beiden orthogonalen Matrizen U und V bilden, ergeben sich die folgenden Gleichungen:

AV=UΣAHU=VΣ.

Da U und V unitäre Matrizen sind, ergibt sich durch Multiplikation der ersten Gleichung mit VH auf der rechten Seite die Gleichung der Singulärwertzerlegung

A=UΣVH.

Die vollständige Singulärwertzerlegung einer m-mal-n-Matrix umfasst Folgendes:

  • m-mal-m-Matrix U

  • m-mal-n-Matrix Σ

  • n-mal-n-Matrix V

Anders ausgedrückt, U und V sind beide quadratisch und Σ hat dieselbe Größe wie A. Wenn A wesentlich mehr Zeilen als Spalten hat (m > n), ist die resultierende m-mal-m-Matrix U groß. Allerdings werden die meisten Spalten in U mit Nullen in Σ multipliziert. In dieser Situation spart die reduzierte Zerlegung Zeit und Speicherkapazität, indem sie eine m-mal-n-Matrix U, eine n-mal-n-Matrix Σ und dieselbe Matrix V erstellt:

In the economy-sized decomposition, columns in U can be ignored if they multiply zeros in the diagonal matrix of singular values.

Die Eigenwertzerlegung ist das richtige Tool für die Analyse einer Matrix, wenn sie eine Zuordnung von einem Vektorraum zu sich selbst darstellt wie für eine gewöhnliche Differenzialgleichung. Allerdings ist die Singulärwertzerlegung das richtige Tool für die Analyse einer Zuordnung von einem Vektorraum in einen anderen Vektorraum, möglicherweise mit einer anderen Dimension. Die meisten Systeme simultaner linearer Gleichungen fallen in diese zweite Kategorie.

Wenn die Matrix A quadratisch, symmetrisch und positiv definit ist, sind ihr Eigenwert und ihr Singulärwert identisch. Wenn jedoch A von der Symmetrie und positiven Definitheit abweicht, vergrößert sich die Differenz zwischen den beiden Zerlegungen. Insbesondere die Singulärwertzerlegung einer reellen Matrix ist immer reell, doch die Eigenwertzerlegung einer reellen, nicht symmetrischen Matrix kann komplex sein.

Für die Beispielmatrix

A =
     9     4
     6     8
     2     7

lautet die vollständige Singulärwertzerlegung wie folgt:

[U,S,V] = svd(A)

U =

    0.6105   -0.7174    0.3355
    0.6646    0.2336   -0.7098
    0.4308    0.6563    0.6194


S =

   14.9359         0
         0    5.1883
         0         0


V =

    0.6925   -0.7214
    0.7214    0.6925

Sie können verifizieren, dass U*S*V' gleich A innerhalb des Rundungsfehlers ist. Für dieses kleine Problem ist die reduzierte Zerlegung nur minimal kleiner.

[U,S,V] = svd(A,0)

U =

    0.6105   -0.7174
    0.6646    0.2336
    0.4308    0.6563


S =

   14.9359         0
         0    5.1883


V =

    0.6925   -0.7214
    0.7214    0.6925

Und auch hier ist U*S*V' gleich A innerhalb des Rundungsfehlers.

SWZ-Stapelberechnung

Wenn Sie eine große Sammlung von Matrizen mit derselben Größe zerlegen müssen, ist es ineffizient, alle Zerlegungen in einer Schleife mit svd durchzuführen. Stattdessen können Sie alle Matrizen in einem mehrdimensionalen Array verketten und pagesvd verwenden, um mit einem einzigen Funktionsaufruf Singulärwertzerlegungen auf allen Arrayseiten durchzuführen.

FunktionVerwendung
pagesvdVerwenden Sie pagesvd, um Singulärwertzerlegungen auf den Seiten eines mehrdimensionalen Arrays durchzuführen. Dies ist eine effiziente Methode zur Durchführung von SWZ an einer großen Sammlung von Matrizen, die alle dieselbe Größe haben.

Angenommen, es liegt eine Sammlung von drei 2x2-Matrizen vor. Verketten Sie die Matrizen mit der Funktion cat zu einem 2x2x3-Array.

A = [0 -1; 1 0];
B = [-1 0; 0 -1];
C = [0 1; -1 0];
X = cat(3,A,B,C);

Verwenden Sie nun pagesvd, um die drei Zerlegungen gleichzeitig durchzuführen.

[U,S,V] = pagesvd(X);

Für jede Seite von X gibt es entsprechende Seiten in den Ausgaben U, S und V. Zum Beispiel befindet sich die Matrix A auf der ersten Seite von X und ihre Zerlegung ist durch U(:,:,1)*S(:,:,1)*V(:,:,1)' gegeben.

SWZ-Niedrigrangapproximation

Bei großen, dünn besetzten Matrizen ist die Verwendung von svd zur Berechnung aller Singulärwerte und Singulärvektoren nicht immer praktikabel. Wenn Sie zum Beispiel nur einige der größten Singulärwerte kennen müssen, stellt die Berechnung aller Singulärwerte einer dünn besetzten 5000x5000-Matrix zusätzlichen Aufwand dar.

Wenn nur ein Teil der Singulärwerte und Singulärvektoren erforderlich ist, sind die Funktionen svds und svdsketch gegenüber svd zu bevorzugen.

FunktionVerwendung
svdsVerwenden Sie svds, um eine Approximation der SWZ mit dem Rang k zu berechnen. Sie können angeben, ob die Untermenge der Singulärwerte die größte oder die kleinste Zahl oder eine Zahl sein soll, die einer bestimmten Zahl am nächsten ist. svds berechnet in der Regel die bestmögliche Approximation mit dem Rang k.
svdsketchVerwenden Sie svdsketch, um eine teilweise SWZ der Eingabematrix zu berechnen, die innerhalb einer bestimmten Toleranz liegt. Während svds erfordert, dass Sie den Rang angeben, bestimmt svdsketch adaptiv den Rang der Matrixskizze basierend auf der angegebenen Toleranz. Die Approximation mit dem Rang k, die svdsketch schließlich verwendet, liegt innerhalb der Toleranz, doch im Gegensatz zu svds ist nicht garantiert, dass es sich um die bestmögliche Approximation handelt.

Angenommen, es liegt eine dünn besetzte 1000x1000-Zufallsmatrix mit einer Dichte von etwa 30 % vor.

n = 1000;
A = sprand(n,n,0.3);

Die sechs größten singulären Werte sind

S = svds(A)

S =

  130.2184
   16.4358
   16.4119
   16.3688
   16.3242
   16.2838

Und die sechs kleinsten singulären Werte sind

S = svds(A,6,'smallest')

S =

    0.0740
    0.0574
    0.0388
    0.0282
    0.0131
    0.0066

Für kleinere Matrizen, die vollständig im Speicher Platz finden, kann full(A), unter Verwendung von svd(full(A)), dennoch schneller sein als svds oder svdsketch. Wenn jedoch wirklich große und dünn besetzte Matrizen vorliegen, lässt sich die Verwendung von svds oder svdsketch nicht vermeiden.

Siehe auch

| | | |

Verwandte Themen