LU rozklad

Tato práce vznikla jako zápočtový projekt na předmětu Algoritmy a datové struktury 1, MFF UK. Pojí se k ní ještě tato implementace v C#.

Efektivní maticové operace s využitím LU rozkladu

LU rozklad nám říká, že čtvercovou regulární matici lze rozložit na součin matic $A = L * U$, kde L je dolní a U je horní trojúhelníkova matice. Zatím nezkoumejme tuto větu do hloubky a pojďme na to postupně.

Využití:

Řešení soustavy rovnic

Mějme soustavu lineárních rovnic $A * x = b$, která se dá řešit gaussovou eliminací. Pokud se nám však podaří převést matici A na součin $L * U$, tedy $L * U * x = b$ (mějme na paměti, že matice L a U jsou trojúhelníkové), vektor $U*x$ snadno spočteme dopřednou a vektor x zpětnou substitucí.

Výpočet determinantu

Mějme za úkol určit determinant matice A, $det(A)$. Pokud se nám podaří převést A na součin $L*U$,  pak $det(A) = det(L*U) = det(L) * det(U)$, což je součin prvků na diagonále L a U, jelikož to jsou matice trojúhelníkové.

Inverze matice

Mějme za úkol určit inverzní matici k A, $A^{-1}$. Hledáme tedy řešení rovnice $A * A^{-1} = I$. K tomuto můžeme přistupovat jako k řešení n soustav rovnic, tj. $A * x_i = e_i$, kde $e_i$ je jednotkový vektor s 1 na i-tém místě a $x_i$ je sloupec námi hledané matice.

Jak je vidět, pro takovýto rozklad matice bychom našli spoustu využití.

Algoritmus:

Idea

Uvažujme normální gaussovu eliminaci. V první fázi je zde snaha převést matici na horní trojúhelníkovou, v našem případě by to mohla být matice U. Zde vyplývají na povrch dvě zajímavé otázky - zda by se nedaly elementární úpravy reprezentovat násobením maticí zleva a zda by součin těchto matic elementárních úprav nebyl dolní trojúhelníkovou maticí. Odpověď je ano. Převod na horní trojúhelníkovou matici je možné realizovat přičítáním alfa-násobku prvního vhodného řádku k řádkům pod ním, abychom vynulovali daný sloupec pod pivotem. Tato elementární úprava se dá reprezentovat násobením zleva maticí, kterou vyrobíme z jednotkové tak, že pod pivot umístíme potřebné koeficienty.

Příklad

Úskalí

Při tomto postupu se může stát, že se na diagonále objeví nula. Při běžné gaussově eliminaci tyto dva řádky (i. a j.) prohodíme. Tuto úpravu můžeme emulovat prohozením i. a j. řádku v matici L, ovšem potom L už nebude trojúhelníková. Vyřešíme to tzv. “permutační maticí” P, která bude na začátku jednotková a které budeme prohazovat řádky. Rozklad tedy bude $A = P L U$.

Vliv permutační matice na operace:

• řešení soustavy rovnic - $PLUx = b$ - ve vstupním vektoru b pouze prohodíme hodnoty podle matice $P$.
• determinant - $det(PLU)$ - někde si evidujeme počet prohození řádků v P, podle kterého výsledný determinant vynásobíme +1 / -1.
• inverze - ve výsledné inverzní matici prohodíme řádky podle matice P

Složitost:

V cyklu musíme probrat řádky matice a od každého z nich odečteme násobek prvního řádku, to vše pro každý ze sloupců. Jsou to tedy 3 vnořené cykly a složitost LU rozkladu je tedy $O(n^3)$.

Při řešení soustavy rovnic je třeba dopočítat vektor Ux dopřednou a x zpětnou substitucí (obojí má složitost $O(n^2)$). Nezískali jsme tedy žádné zrychlení oproti Gaussově eliminaci (složitost $O(n^3)$). Zrychlení ale získáme ve chvíli, kdy máme vyřešit stejnou soustavu A pro více pravých stran. Když už jednou spočítáme LU rozklad, vyřešení soustavy pro libovolnou pravou stranu má složitost $O(n^2)$.

Efektivnost výpočtu determinantu dle definice musíme sečíst n! členů, v každém z nich vynásobit n čísel, složitost je tedy $O(n*n!)$. Determinant lze vypočítat pomocí gaussovy eliminace v čase $O(n^3)$. Náš algoritmus funguje také v čase $O(n^3)$, jelikož LU rozklad je jiný způsob přístupu ke gaussově eliminaci.

Výpočet inverzní matice gaussovou eliminací má složitost $O(n^3)$. Nezískali jsme tedy žádné zrychlení.