4.4. Linearna algebra

U Sageu postoje vektori i matrice kao specijalni objekti, ali mi ćemo koristiti NumPy polja koja nude sličnu funkcionalnost zahvaljujući rutinama NumPy modula za linearnu algebru. (Za “domaću” Sage linearnu algebru vidi ovaj odjeljak.)

sage: import numpy as np
sage: import numpy.linalg as la

sage: v1 = np.array([1,1,2])
sage: v2 = np.array((2,2,4))

Množenje skalarom je prirodno:

sage: 3*v1
array([3, 3, 6])

Skalarni i vektorski produkt vektora:

sage: np.dot(v1, v2)
12

sage: np.cross(v1, v2)
array([0, 0, 0])

Norma (“duljina”) vektora:

sage: la.norm(v2)
4.8989794855663558

Množenje matrica te množenje matrice i vektora ide na prirodan način:

sage: A = np.array([[1, 2, 1], [4, 3, 3], [9, 1, 7]]); A
array([[1, 2, 1],
       [4, 3, 3],
       [9, 1, 7]])

sage: np.dot(A, A)
array([[18,  9, 14],
       [43, 20, 34],
       [76, 28, 61]])

sage: np.dot(A, v1)
array([ 5, 13, 24])

Determinanta i inverz matrice:

sage: la.det(A)   # = 7
-6.9999999999999982

sage: la.inv(A)
array([[-2.57142857,  1.85714286, -0.42857143],
       [ 0.14285714,  0.28571429, -0.14285714],
       [ 3.28571429, -2.42857143,  0.71428571]])

Numerika po defaultu ide s double precision točnošću, pa kad provjeravamo da li je \(A^{-1} A = 1\) ne dobijemo egzaktno jediničnu matricu:

sage: should_be_unit = np.dot(la.inv(A), A); should_be_unit
array([[  1.00000000e+00,   1.38777878e-15,   8.32667268e-16],
       [  2.77555756e-17,   1.00000000e+00,  -2.77555756e-17],
       [ -1.22124533e-15,  -7.77156117e-16,   1.00000000e+00]])

Moguće je zatražiti od NumPyja da ispisuje ovakve male brojeve kao nule:

sage: np.set_printoptions(suppress=True)
sage: should_be_unit
array([[ 1.,  0.,  0.],
       [ 0.,  1., -0.],
       [-0., -0.,  1.]])

Zadatak 1

Za datu matricu A definiramo svojstvene vektore (eigenvectors) \({\bf v}\) i njima pripadajuće svojstvene vrijednosti \(\lambda\) (eigenvalues) kao rješenja matrične jednadžbe

\[A {\bf v} = \lambda {\bf v}\;.\]

Pomoću funkcije la.eig() odredite svojstvene vrijednosti i svojstvene vektore matrice

\[\begin{split}\left(\begin{array}{cc} 2.3 & 4.5 \\ 6.7 &-1.2 \end{array}\right) \;,\end{split}\]

i provjerite da dobivena rješenja zaista zadovoljavaju gornju jednadžbu.

Zadatak 2

Kreirajte 3x3 matricu sa slučajnim realnim brojevima između 0 i 10. Invertirajte je i pomnožite s originalnom matricom te se uvjerite da dobijete jediničnu matricu.

Dijagonalizacija matrice A je pronalaženje njenog rastava oblika

\[A = P D P^{-1}\]

gdje je \(D\) dijagonalna matrica. Takav rastav se također izvodi fukcijom la.eig(), koja daje svojstvene vrijednosti i svojstvene vektore matrice. Naime, stupci matrice \(P\) su upravo svojstveni vektori od \(A\), a dijagonalna matrica \(D\) na dijagonali ima upravo odgovarajuće svojstvene vrijednosti:

sage: A = np.array([[3, 1], [1, 3]]); A
array([[3, 1],
       [1, 3]])

sage: evs, P = la.eig(A)

sage: D = np.diag(evs); D
array([[ 4.,  0.],
       [ 0.,  2.]])
sage: np.dot(np.dot(la.inv(P), A), P)  # provjerevamo P^-1 A P = D
array([[ 4.,  0.],
       [ 0.,  2.]])

Neke matrice nije moguće dijagonalizirati, u slučaju čega će matrice \(P\) i \(D\) koje vraća la.eig() i dalje zadovoljavati

\[AP = PD\]

ali \(P\) neće biti invertibilna:

sage: A = np.array([[1, 1], [0, 1]]); A
array([[1, 1],
       [0, 1]])
sage: evs, P = la.eig(A)

sage: np.dot(A, P) - np.dot(P, np.diag(evs))
array([[ 0.,  0.],
       [ 0.,  0.]])

Neinvertibilnost od \(P\) se ogleda u ogromnim brojevima koje dobijemo kad “invertiramo” tu matricu:

sage: print la.inv(P)
[[  1.00000000e+00   4.50359963e+15]
 [  0.00000000e+00   4.50359963e+15]]