K. Kumerički, verzija 1.1, 2011-10-25

Linearna algebra

O nizovima objekata (brojeva, simbola, ...) biti će više govora u slijedećem poglavlju, a ovdje nam treba samo nekoliko osnovnih ideja. Često ćemo sretati dvije vrste nizova objekata: liste (engl. list) i tuplove (engl. tuple). Liste su nizovi elemenata odvojeni zarezom u uglatim zagradama, a tuplovi su to isto samo u okruglim zagradama. Razlike između liste i tupla će biti razjašnjene kasnije; u ovom poglavlju će obje strukture funkcionirati bez razlika.

{{{id=24| a = [1,1,2] # lista b = (2, 2, 4) # tupl /// }}}

Pristup pojedinim elementima se izvodi indeksom u uglatim zagradama (i za liste i za tuplove!), gdje je brojanje kao u C-u, tj. počinje od nule: prvi element liste a je a[0]:

{{{id=54| print a[2] b[0] /// }}}

Ukoliko želimo napraviti pridruživanje gdje se vrijednosti uzimaju iz liste na slijedeći način je to moguće izvesti za cijelu listu odjednom. (Tzv. raspakiravanje.)

{{{id=55| b1, b2, b3 = b; b2 /// }}}

Naravno, elementi listi mogu biti druge liste:

{{{id=56| c = [a, b, [a]] avec, bvec, [cvec] = c; cvec # ovo je korisno za donji zadatak sa svojstvenim vrijednostima matrice /// }}}

Vektori i matrice se konstruiraju pomoću funkcija vector() i matrix() kojima kao argument dajemo listu elemenata, odnosno listu listi elemenata (listu redak-vektora)

{{{id=18| vec1 = vector([1,1,2]) vec2 = vector(b) # moze i tupl /// }}}

Množenje skalarom je prirodno:

{{{id=66| 3*vec1 /// }}}

Skalarni produkt vektora ...

{{{id=21| vec1.inner_product(vec2) /// }}}

... se može zapisivati kao i obično množenje:

{{{id=19| vec1*vec2 /// }}}

Vektorski produkt (nema zapisa s $\times$!)

{{{id=22| vec1.cross_product(vec2) /// }}}

Norma ("duljina") vektora:

{{{id=102| vec2.norm() /// }}}

Svi vektori iste vrste i dimenzije su elementi apstraktnog vektorskog prostora koji je dostupan putem metode parent():

{{{id=106| vecspace=vec1.parent(); vecspace /// }}}

Ovdje je važno uočiti da je je prostor definiran na "prstenu" cijelih brojeva. Naravno, vektori mogu biti i nad drugim prstenovima/poljima:

{{{id=107| print vector([1/3, 1/4]).parent() vector([1., 2.]).parent() /// }}}

Baza vektorskog prostora:

{{{id=46| vecspace.basis() /// }}}

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

{{{id=20| mat = matrix([[1, 2, 1], [4, 3, 3], [9, 1, 7]]); mat /// }}} {{{id=27| mat*mat /// }}} {{{id=17| mat*vec1 /// }}} {{{id=29| ~mat # inverz matrice, moze i m.inverse() ili m^-1 /// }}} {{{id=1| ~mat*mat # provjera /// }}}

Matrice isto imaju svoje apstraktne prostore kojima pripadaju:

{{{id=45| print mat.parent() matinv.parent() /// }}} {{{id=47| mat2 = matrix([[1., 2.], [3., 4.]]) print mat2.parent(); mat2 /// }}}

Ovo defaultno polje realnih brojeva nije sasvim pogodno za numeričke račune pa je potrebno matrice s realnim koeficijentima kreirati uz eksplicitnu deklaraciju polja RDF  ("real double field"):

{{{id=48| mat2 = matrix(RDF, [[1., 2.], [3., 4.]]) print mat2.parent(); mat2 /// }}}

Pristup pojedinim elementima matrice se isto izvodi indeksiranjem:

{{{id=30| mat[0, 0] = 0; mat /// }}}

Zadatak 2-4.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}\;.$$

Odredite svojstvene vrijednosti  i svojstvene vektore  matrice

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

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

Zadatak 2-4.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.

Elementi vektora i matrica mogu biti i simboli:

{{{id=136| var('x y z') /// }}} {{{id=125| vec3 = vector([x, y, z]) vec3.norm() /// }}}

Međutim, ukoliko pokušamo već stvorenom vektoru nad poljem nekih brojeva zamijeniti neki element simbolom, to ne ide:

{{{id=123| vec1[2] = x /// }}}

Riječ je o tome da su npr. čisto realni vektori (nad poljem RDF) interno reprezentirani kao C-polja radi optimizacije. Da bismo mogli napraviti ovo što želimo trebamo prvo konvertirati vektor u simbolički vektor. Tu konverziju radi odgovarajući vektorski prostor nad simboličkim prstenom (SR, symbolic ring). Taj prostor dobijemo pomoću metode parent() vektora istog ranga, ali koji već jest simbolički. (Za detalje o takvim konverzijama vidi prvih par odjeljaka ovdje.)

{{{id=126| vec3.parent() /// }}} {{{id=127| vec1_symb = vec3.parent()(vec1) # konverzija vec1 u simbolicki vektor /// }}} {{{id=128| vec1_symb[1] = x; vec1_symb # sad ide /// }}}

Zadatak 2-4.3 (*): Isto kao zadatak 2-4.3, ali prije invertiranja matrice zamijenite njen središnji element simbolom x. Nakon invertiranja i množenja uvrstite za x slučajni broj (random()) i uvjerite se da dobivate jediničnu matricu.

Dijagonalizacija matrice A  je pronalaženje njenog rastava oblika

$$

A = P D P^{-1}

$$

gdje je $D$ dijagonalna matrica. To se izvodi metodom eigenmatrix_right() koja vraća matrice $P$ i $D$:

{{{id=135| A = matrix(QQ, [[3, 1], [1, 3]]) # metoda is_diagonalizable ne radi nad ZZ print A.is_diagonalizable() D, P = A.eigenmatrix_right(); (D, P) /// }}} {{{id=139| A == P*D*(~P) # provjera /// }}}

Treba uočiti da su elementi dijagonalne matrice $D$ i  stupci matrice $P$ upravo svojstvene vrijednosti odnosno svojstveni vektori od $A$.

{{{id=150| mat3.eigenvectors_right() /// }}}

Neke matrice nisu dijagonalizabilne, u slučaju čega će  matrice $P$ i $D$ koje vraća eigenmatrix_right() i dalje zadovoljavati

$$

AP = PD

$$

ali $P$ neće biti invertibilna:

{{{id=140| A = matrix(QQ, [[1, 1], [0, 1]]) print A.is_diagonalizable() D, P = A.eigenmatrix_right(); (D, P) /// }}} {{{id=144| A*P == P*D /// }}} {{{id=152| ~P /// }}}

U takvim slučajevima od koristi može biti i vrlo popularni rastav na singularne vrijednosti (singular value decomposition, SVD):

$$

A = U S V^{\dagger}

$$

gdje su $U$ i $V$ unitarne, a $S$ dijagonalna matrica.  (Jedino treba imati na umu da je odgovarajuća metoda SVD() trenutno implementirana samo za matrice nad realnim RDF poljem pa je po potrebi potrebno prvo provesti konverziju matrice.)

{{{id=143| U, S, V = matrix(RDF, A).SVD() /// }}} {{{id=147| U*S*V.conjugate_transpose() /// }}}

Zadatak 2-4.4Stupci matrice $U$ u SVD rastavu su svojstveni vektori matrice $A A^\dagger$, a odgovarajuće svojstvene vrijednosti su kvadrati elemenata dijagonale matrice $S$. Uvjerite se u to eksplicitno na gornjem primjeru.

Literatura

Robert A. Beezer, A First Course in Linear Algebra, (Sage-enhanced textbook)

{{{id=153| /// }}}