Link Search Menu Expand Document

Lineare Algebra mit NumPy

Matrix-Operationen

NumPy bietet alle Operationen an, die man auf Matrizen anwenden kann. Diese werden im folgenden in ein paar Beispielen gezeigt, ohne groß erläutert zu werden. Die mathematische Bedeutung der einzelnen Operationen sei den entsprechenden Mathematikvorlesungen vorbehalten bzw. den einschlägigen Lehrbüchern.

  • Definition einiger Matrizen mit dem Wert 1.0 an jeder Stelle
  • Transponieren der Matrizen mit .T
a = np.ones((3, 2))
# [[ 1.,  1.],
#  [ 1., 1.],
#  [ 1.,  1.]]

a.T # Transponierte Matrix
# [[ 1.,  1.,  1.],
#  [ 1.,  1.,  1.]]

b = np.ones((2, 3))
# [[ 1.,  1.,  1.],
#  [ 1.,  1.,  1.]]
  • Matrix-Multiplikation mit dot
  • Dimensionen der Matrizen müssen passen
np.dot(a, b)
# [[ 2.,  2.,  2.],
#  [ 2.,  2.,  2.],
#  [ 2.,  2.,  2.]]

np.dot(b, a)
# [[ 3.,  3.],
#  [3., 3.]]

np.dot(A, B.T)
# ValueError: shapes (3,2) and (3,2) not aligned:  ...
# ... 2 (dim 1) != 3 (dim 0)

Für die Multiplikation der Matrizen müssen die Dimensionen der beiden Operanden passen. Kann auch kein Broadcasting erfolgen (siehe oben), kommt es zu einem ValueError.

  • Summieren aller Elemente oder entlang einer Achse mit sum
  • Kumulative Summe mit cumsum
a = np.array([[1, 2, 3], [4, 5, 6]])

# Spalten und Zeilensumme
a.sum(axis=0) # -> [5, 7, 9]
a.sum(axis=1) # -> [ 6, 15]
a.sum() # -> 21

# Kummulative Summe
a.cumsum(axis=0) # -> [[1, 2, 3], [5, 7, 9]]

a.cumsum(axis=1) # -> [[ 1,  3,  6], [ 4,  9, 15]]
  • Maximum und Minimum mit max und min
  • Auch entlang der Achsen möglich
a = np.array([[1, 2, 3], [4, 5, 6]])

# Minumum und Maximum
a.min() # -> 1
a.max() # -> 6

# Maximum und Minimum entlang einer Achse
a.min(axis=1) # -> [1, 4]
a.max(axis=1) # -> [3, 6]
  • eye(n): Erzeugt die Identitätsmatrix
  • trace(a): Spur berechnen
np.eye(3) # 3x3 Idenititätsmatrix
# [[ 1.,  0.,  0.],
#  [ 0.,  1.,  0.],
#  [ 0.,  0.,  1.]]

# Spur berechnen
a = np.array([[1, 2, 3], [4, 5, 6], [7, 6, 8]])
np.trace(a) # -> 14
  • column_stack((a, b, ..)): Erzeugt aus Vektoren eine Matrix (Spalten)
  • row_stack((a, b, ..)): Erzeugt aus den Vektoren eine Matrix (Zeilen)
# Zwei Vektoren
a = np.array([1, 2, 3])
b = np.array([4, 5, 6])

np.column_stack([a, b])
# [[1, 4],
#  [2, 5],
#  [3, 6]]

np.row_stack((a, b))
# [[1, 2, 3],
# [4, 5, 6]]

Lineare Algebra mit numpy

NumPy ist so umfangreich, dass die erweiterten Funktionen in Untermodulen untergebracht sind. Diese muss man mit einem weiteren import laden, z. B. import numpy.linalg für Funktoonen zur linearen Algebra. Verwendet man keinen Alias (siehe oben), dann sind die Funktionen mit dem Prefix numpy.linalg. ansprechbar, bzw. np.linalg., wenn NumPy selbst mit dem Alias np geladen wurde.

Um sich die vielen Punkte im Aufruf zu sparen, wird man diese Zusatzmodule häufig mit einem neuen Alias laden, z. B. import numpy.linalg as lin.

import numpy as np
import numpy.linalg as lin

a = eye(3)
# [[1., 0., 0.],
#  [0., 1., 0.],
#  [0., 0., 1.]]

n = lin.norm(a)
print(n) # -> 1.7320508075688772

Importieren der Funktionen für lineare Algebra mit import numpy.linalg

  • norm(A): Norm der Matrix (des Vektors) berechnen
  • inv(A): Inverse Matrix berechnen
  • solve(A, b): Gleichung Ax = b lösen
  • lstsq(A, b): Kleinste Fehlerquadrate berechnen
  • eig(A): Eigenwert-Vektor berechnen
  • eigvals(A): Eigenwert berechnen
  • pinv(A): Pseudo-inverse (Moore-Penrose) Matrix

Zufallszahlen

Für Simulationen oder zum Test von Funktionen benötigt man Zufallswerte. Auch hier bietet NumPy eine Reihe von nützlichen Funktionen, um Matrizen mit entsprechenden, zufälligen Werten zu erzeugen.

Paket numpy.random

  • rand(d0, d1, ..., dn): Matrix mit Zufallszahlen (Größe nach Dimensionen)
  • randn(d0, d1, ...,dn): Matrix mit normierten Zufallszahlen (Größe nach Dimensionen)
  • randint(lo, hi, size): Liste von size zufälligen Integers zwischen lo und hi
  • shuffle(a): Mischen der Matrix a. Matrix wird verändert
  • permutation(a): Mischen der Matrix a und zurückgeben einer neuen

Beispiel: Tragwerk

In einem komplexeren Beispiel wollen wir zeigen, wie man mit NumPy eine statische Berechnung vereinfachen kann. Hierzu gehen wir von einem Tragwerk aus, das in folgender Zeichnung visualisiert ist.

Am Punkt F ist das Tragwerk drehbar gelagert, am Punkt A ist es fest mit dem Untergrund verbunden.

Hierbei sind die Längen a und die Kräfte Q1 und Q2 dem Betrag nach gegeben.

\[ \begin{align*} a = 3m, Q_1 = 1kN, Q_2 = 2kN \end{align*} \]

In einem ersten Schritt wird das Tragwerk von dem Fundament freigeschnitten und die entsprechenden Kräfte (und Drehmomente) werden notiert.

\[ \begin{align*} a = 3m, Q_1 = 1kN, Q_2 = 2kN \end{align*} \]

Im nächsten Schritt werden die einzelnen Elemente voneinander freigeschnitten und die Kräfte werden eingetragen.

Da das Tragwerk im Stillstand ist, müssen die Kräfte sich gegenseitig aufheben. Dies notieren wir in Form einer Reihe von Gleichungen.

\[ \begin{align*} H_A + H_C &= 0 \\ V_A + V_C &= 0 \\ M_A + a \cdot V_C - a \cdot H_C &= 0 \\ -H_C + H_D &= 0 \\ -V_C - Q_1 + V_D - Q_2 &= 0 \\ a \cdot V_D - a \cdot Q_2 &= 0 \\ -H_D + H_F &= 0 \\ -V_D + V_F &= 0 \\ a \cdot H_D + a \cdot V_D &= 0 \\ a = 3m, Q_1 = 1kN, Q_2 = 2kN \end{align*} \]

Wir setzen die bekannten Werte für a, Q1 und Q2 ein.

\[ \begin{align*} H_A + H_C &= 0 \\ V_A + V_C &= 0 \\ M_A + 3 \cdot V_C - 3 \cdot H_C &= 0 \\ -H_C + H_D &= 0 \\ -V_C - 1 + V_D - 2 &= 0 \\ 3 \cdot V_D - 3 \cdot 2 &= 0 \\ -H_D + H_F &= 0 \\ -V_D + V_F &= 0 \\ 3 \cdot H_D + 3 \cdot V_D &= 0 \\ \end{align*} \]

Das Gleichungssystem lässt sich als Produkt von Matrix und Vektor schreiben.

\[ \begin{align*} \begin{bmatrix} 1 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & -3 & 3 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & -1 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & -1 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 3 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & -1 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & -1 & 0 & 1 \\ 0 & 0 & 0 & 0 & 0 & 3 & 3 & 0 & 0 \\ \end{bmatrix} \begin{bmatrix} H_A \\ V_A \\ M_A \\ H_C \\ V_C \\ H_D \\ V_D \\ H_F \\ V_F \\ \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ 0 \\ 0 \\ 1 + 2 \\ 3 \cdot 2 \\ 0 \\ 0 \\ 0 \\ \end{bmatrix} \end{align*} \]

Mit der Darstellung als Matrix haben wir alles, um das Problem von NumPy lösen zu lassen. D. h. wir formulieren die Matrix als NumPy-Array und lassen das System von NumPy über die solve-Funktion lösen.

import numpy as np
import numpy.linalg as lin

A = np.array(
    [[ 1, 0, 0, 1, 0, 0, 0, 0, 0 ],
     [ 0, 1, 0, 0, 1, 0, 0, 0, 0 ],
     [ 0, 0, 1, -3, 3, 0, 0, 0, 0 ],
     [ 0, 0, 0, -1, 0, 1, 0, 0, 0 ],
     [ 0, 0, 0, 0, -1, 0, 1, 0, 0 ],
     [ 0, 0, 0, 0, 0, 0, 3, 0, 0 ],
     [ 0, 0, 0, 0, 0, -1, 0, 1, 0 ],
     [ 0, 0, 0, 0, 0, 0, -1, 0, 1 ],
     [ 0, 0, 0, 0, 0, 3, 3, 0, 0 ]])

b = np.array([ 0, 0, 0, 0, 3, 6, 0, 0, 0])

x = lin.solve(A, b)

print(x) # => [ 2.  1. -3. -2. -1. -2.  2. -2.  2.]

Damit haben wir alle Kräfte aus dem Tragwerk berechnet.

  • HA = 2 kN
  • VA = 1 kN
  • MA = -3000 Nm
  • HC = -2 kN
  • VC = -1 kN
  • HD = -2 kN
  • VD = 2 kN
  • HF = -2 kN
  • VF = 2 kN

Copyright © 2022 Thomas Smits