Odwracanie kolejności elementów w macierzy – algorytm
- Szczegóły
Odwracanie macierzy jest uogólnieniem dzielenia liczb. Macierz odwrotna do macierzy kwadratowej A (ozn. A-1) to taka macierz, która pomnożona przez A (z lewej lub prawej strony) daje macierz jednostkową. Innymi słowy, gdy pomnożymy macierz A przez macierz do niej odwrotną (lub wykonamy to mnożenie w odwrotnej kolejności, tzn. A-1 * A), otrzymamy macierz jednostkową.
Macierz odwrotna istnieje tylko dla tzw. macierzy nieosobliwych (czyli takich, których wyznacznik jest różny od zera). Macierz odwrotna nie istnieje dla macierzy osobliwych.
Chcesz sprawdzić czy dobrze obliczyłeś/aś macierz odwrotną?
Macierz odwrotna do macierzy transponowanej jest równa transponowanej macierzy odwrotnej (nie ważna jest kolejność, tzn. (AT)-1 = (A-1)T).
Dlaczego? Wzorem zapiszemy to tak: det(A-1). Żeby to obliczyć, musielibyśmy najpierw odwrócić macierz A i dopiero potem obliczyć wyznacznik macierzy A-1.
Przeczytaj także: Sterowniki i usterki ASUS K52J
\[\det(A^{-1})=(detA)^{-1}\]
która mówi, że zamiast męczyć się z odwracaniem samej macierzy, możemy najpierw obliczyć wyznacznik... Wystarczy przypomnieć sobie, że wyznacznik macierzy to nic innego jak zwykła liczba rzeczywista, więc (detA)-1 to po prostu jakaś liczba podniesiona do potęgi "-1", czyli odwrotność, tzn. 1/detA.
Eliminacja Gaussa
Na macierzach (podobnie jak na układach równań) można wykonywać trzy proste operacje: dodać do jednego wiersza macierzy inny wiersz pomnożony przez liczbę, zamienić dwa wiersze miejscami, mnożyć wiersz przez liczbę różną od zera. Te trzy opisane powyżej operacje nazywamy operacjami elementarnymi na wierszach. Operacje elementarne na wierszach macierzy odpowiadają takim samym operacjom na równaniach układu równań.
Rozważ układ równań:
Danemu układowi równań odpowiada następująca macierz:
Przeczytaj także: Zastosowanie wężyków do filtra osmozy
Wykonamy teraz kilka operacji elementarnych na równaniach układu i jednocześnie na wierszach odpowiadającej mu macierzy. Na początku pomnóżmy pierwszy wiersz przez liczbę \(2\):
Napis \(w_1 \rightarrow w_1\cdot 2\) jest symbolicznym zapisem tego, że "wiersz pierwszy przekształcamy w wiersz pierwszy pomnożony przez \(2\)". Teraz do drugiego wiersza dodajemy pierwszy wiersz pomnożony przez \(2\):
W wyniku powyższej operacji pozbyliśmy się jednej niewiadomej z drugiego równania. Teraz pozbędziemy się niewiadomej \(x_1\) z trzeciego równania. W tym celu do trzeciego równania dodajemy pierwsze równanie pomnożone przez \(3\):
Teraz będziemy pozbywać się kolejnych niewiadomych maksymalnie upraszczając układ równań (oraz macierz). Wykonywane operacje elementarne będą zapisywane przy strzałkach.
Wykonanie wszystkich powyższych operacji elementarnych doprowadziło nas do rozwiązania układu równań:
Przeczytaj także: Odwrócona osmoza: Twój przewodnik
Ciągłe przepisywanie całego układu równań jest dość czasochłonne. Dlatego warto wykonywać działania jedynie na macierzy odpowiadającej danemu układowi, aby na końcu tylko odczytać z niej rozwiązanie. Przy takim podejściu całe zadanie sprowadza się do przekształcenia macierzy na postać schodkową zredukowaną. Z takiej postaci od razu widać rozwiązania układu równań.
Podczas eliminacji Gaussa, możemy dostać algorytm (tzw. ijk, ikj, bikj() --- nazwy pochodzą od sposobu organizacji pętli, zobacz poniżej), zaimplementowanych w programie w języku C, wykorzystującym technikę pozwalającą przechowywać macierze w pamięci komputera kolumnowo. Dla porównania zmierzyliśmy czas wykonania tej samej operacji przy użyciu wyspecjalizowanych bibliotek z pakietów BLAS (algorytm DGEMM) i ATLAS (algorytm ATLAS DGEMM).
Okazuje się, że o ile liczba operacji arytmetycznych jest identyczna, powodem różnic musi być odmienne wykorzystanie pamięci cache, wynikające z odmiennej organizacji dostępu do pamięci w naszych algorytmach. Przedyskutujmy to dokładniej.
Algorytm DGEMM z pakietu BLAS --- to właśnie profesjonalny algorytm blokowy, ale niezoptymalizowany na naszą architekturę.
Organizacja pamięci macierzy
W Fortranie, elementy macierzy są przechowywane w pamięci kolumnami, tzn. Dla odmiany, C przechowuje w pamięci elementy macierzy wierszami, tzn. Co więcej, standard języka C nie gwarantuje, że kolejne wiersze macierzy będą przechowywane w przylegających do siebie obszarach pamięci. To zaś powoduje kolejne komplikacje. To jest jeden z licznych powodów, dla których, posługując się dwuwymiarowymi macierzami w C, będziemy stosowali pewien prosty trick. Będziemy przechowywać elementy macierzy w jednym, długim wektorze. Uzyskujemy w ten sposób zwarty obszar pamięci. C (np. (i,j)-ty element macierzy w wektorze przechowującym jej elementy.
Biblioteki BLAS i LAPACK
BLAS (Basic Linear Algebra Subprograms) oraz LAPACK (Linear Algebra PACKage) to wyspecjalizowane biblioteki. Istnieje inna wersja optymalizowanych BLASów, tzw. Goto BLAS.
Porównanie czasu działania kodu w C, implementującego algorytm rozkładu LU z wykładu, z czasem działania procedury DGESV z LAPACKa, niezoptymalizowanej i zoptymalizowanej (ATLAS) na daną architekturę.
BLAS to zestaw procedur do wykonywania podstawowych operacji algebry liniowej: skalarami, wektorami i macierzami. Istnieją różne poziomy BLAS (w dwóch precyzjach).
- BLAS Level 1 -- działania typu wektor--wektor, np. albo obliczanie normy wektora.
- BLAS Level 2 -- działania typu macierz--wektor, np. Zapis algorytmów z użyciem BLAS Level 2 umożliwia potencjalnie przyspieszenie programu, m.in. ze względu na to, że zoptymalizowane procedury BLAS Level 2 mogą np.
- BLAS Level 3 -- operacje typu macierz--macierz, np.
Nazwy procedur w BLAS są dość regularne i łatwo je odgadywać. RR oznacza rodzaj zadania, np. FF wreszcie określa samo zadanie, np. SV oznacza SolVe (w domyśle: układ równań), MV --- Matrix-Vector (w domyśle: mnożenie), EV --- EigenValues, czyli wartości własne, itp. Są też warianty trzyliterowe, np. poniżej. Każda z tych procedur ma swój wariant "ekspercki", np.
tags: #odwracanie #kolejności #elementów #w #macierzy #algorytm
