Crout-decompositie

De Crout-decompositie is een algoritme voor de LU-decompositie van een vierkante niet-singuliere matrix in een benedendriehoeksmatrix $L$ en een bovendriehoeksmatrix $U.$ In de matrix $U$ zijn de elementen op de hoofddiagonaal gelijk aan 1.

De methode is genoemd naar Prescott Durand Crout, een wiskundige aan het Massachusetts Institute of Technology die de methode in 1941 heeft beschreven.

Beschrijving

Via Gauss-eliminatie

Zij $A$ een $n\times n$-matrix met elementen $a_{i,j}.$ De Crout-decompositie is formeel gezien een Gauss-eliminatie. Het algoritme begint dan met als $L$ de $n\times n$-eenheidsmatrix en $U$ gelijk aan $A.$ In $n-1$ stappen wordt $U$ herleid tot een bovendriehoeksmatrix en $L$ tot een benedendriehoeksmatrix, zodanig dat het matrixproduct $LU$ steeds gelijk blijft aan $A.$

In de $k$-de stap worden het elementen op de hoofddiagonaal in de $k$-de kolom van $U$ op een gebracht en de elementen eronder op nul. Dit gebeurt door rij $k$ te delen door het pivot-element $u_{k,k}.$ Deze waarde kopiëren we in de corresponderende positie in matrix $L.$ Daarna vermenigvuldigen we de $k$-de rij met $u_{j,k}$ en trekken ze dan af van rij $j,$ wat een nieuwe rij $j$ geeft met nul in de $k$-de kolom ($j>k$). De factoren $u_{j,k}$ kopiëren we eveneens in de corresponderende kolom van $L.$

Indien $U$ een element $u_{k,k}$ op de hoofddiagonaal heeft dat nul is, gaat dit natuurlijk niet op. In dat geval moet de rij $k$ verwisseld worden met een onderliggende rij $j$ waarin $u_{j,k}$ niet nul is. De verwisseling gebeurt ook in $L.$ De verwisselingen kunnen bijgehouden worden in een permutatiematrix $P.$ In het begin is $P$ de eenheidsmatrix. Wanneer twee rijen in $U$ en $L$ verwisseld worden, worden dezelfde rijen in $P$ verwisseld. De decompositie is dan niet langer $A=LU$ maar $A=PLU.$

Rechtstreekse berekening

Indien men rekening houdt met de specifieke structuur van $L$ en $U,$ kan men de coëfficiënten ervan een voor een berekenen. Immers als men het matrixproduct $A=LU$ uitschrijft, verkrijgt men een stelsel van $n\times n$ lineaire vergelijkingen in de $n\times n$ onbekende coëfficiënten van $L$ en $U.$ Deze vergelijkingen kan men zo rangschikken dat elke volgende vergelijking een nieuwe coëfficiënt geeft; bijvoorbeeld:

$a_{1,1}=l_{1,1},a_{2,1}=l_{2,1}$ enz. (dit is de eerste kolom van $L$)

$a_{1,2}=l_{1,1}{u}_{1,2}$ of $u_{1,2}=\frac{a_{1,2}}{a_{1,1}}$

$a_{1,3}=l_{1,1}{u}_{1,3}$ of $u_{1,3}=\frac{a_{1,3}}{a_{1,1}}$ enz. Dit maakt de eerste rij van $U$ volledig.

$a_{2,2}=l_{2,1}{u}_{1,2}+l_{2,2}$ of $l_{2,2}=a_{2,2}-l_{2,1}{u}_{1,2}$

$a_{3,2}=l_{3,1}{u}_{1,2}+l_{3,2}$ of $l_{3,2}=a_{3,2}-l_{3,1}{u}_{1,2}$ enz. (dit vervolledigt de tweede kolom van $L$)

Op deze manier kan men afwisselend een kolom van $L$ en een rij van $U$ rechtstreeks berekenen. Dit is de "compacte" vorm van de Crout-decompositie.

Voorbeeld

Gegeven is de matrix:

$A=\left(\begin{array}{cccc}3& -7& -2& 2\\ -3& 5& 1& 0\\ 6& -4& 0& -5\\ -9& 5& -5& 12\end{array}\right)$

Neem als matrix $L$ de eenheidsmatrix en $U$ gelijk aan $A$.

$LU=\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right)\left(\begin{array}{cccc}3& -7& -2& 2\\ -3& 5& 1& 0\\ 6& -4& 0& -5\\ -9& 5& -5& 12\end{array}\right)$

Eerste stap

Deel de eerste rij van $U$ door 3, en vermenigvuldig de eerste rij achtereenvolgens met −3, 6, 9 en trek ze af van de tweede, derde en vierde rij zodat de elementen onder de diagonaal in de eerste kolom van $U$ nul worden. De getallen 3, −3, 6 en 9 kopiëren we naar de eerste kolom van $L,$ die dus een kopie wordt van de eerste kolom van $U.$ Dit geeft als tussenstand:

$LU=\left(\begin{array}{cccc}3& 0& 0& 0\\ -3& 1& 0& 0\\ 6& 0& 1& 0\\ -9& 0& 0& 1\end{array}\right)\left(\begin{array}{cccc}1& -7/3& -2/3& 2/3\\ 0& -2& -1& 2\\ 0& 10& 4& -9\\ 0& -16& -11& 18\end{array}\right)$

Men kan nagaan dat het product van deze twee matrices de oorspronkelijk matrix $A$ is.

Tweede stap

Om het element op de hoofddiagonaal in de tweede kolom van $U$ op een te krijgen, moet de tweede rij van $U$ gedeeld worden door −2. Om de rest van de tweede kolom van $U$ op nul te krijgen, moet de tweede rij van $U$ achtereenvolgens vermenigvuldig worden met 10 en −16, en afgetrokken van de derde, respectievelijk de vierde rij. De waarden −2, 10 en −16 gaan naar de tweede kolom van $L:$, die is dus een kopie van de tweede kolom van $U$ vanaf de tweede rij. Dit geeft:

$LU=\left(\begin{array}{cccc}3& 0& 0& 0\\ -3& -2& 0& 0\\ 6& 10& 1& 0\\ -9& -16& 0& 1\end{array}\right)\left(\begin{array}{cccc}1& -7/3& -2/3& 2/3\\ 0& 1& 1/2& -1\\ 0& 0& -1& 1\\ 0& 0& -3& 2\end{array}\right)$

Derde stap

In de derde stap wordt de derde rij van $U$ gedeeld door −1, waarna de derde rij vermenigvuldig wordt met −3 en afgetrokken van de vierde rij. De factoren −1 en −3 gaan naar de derde kolom van $L$:

$LU=\left(\begin{array}{cccc}3& 0& 0& 0\\ -3& -2& 0& 0\\ 6& 10& -1& 0\\ -9& -16& -3& 1\end{array}\right)\left(\begin{array}{cccc}1& -7/3& -2/3& 2/3\\ 0& 1& 1/2& -1\\ 0& 0& 1& -1\\ 0& 0& 0& -1\end{array}\right)$

Vierde stap

In de vierde en laatste stap moet alleen nog het element in de linker benedenhoek van $U$ op een gebracht worden, door het te delen door −1, en die waarde te kopiëren naar $L.$ Dit geeft als eindresultaat:

$LU=\left(\begin{array}{cccc}3& 0& 0& 0\\ -3& -2& 0& 0\\ 6& 10& -1& 0\\ -9& -16& -3& -1\end{array}\right)\left(\begin{array}{cccc}1& -7/3& -2/3& 2/3\\ 0& 1& 1/2& -1\\ 0& 0& 1& -1\\ 0& 0& 0& 1\end{array}\right)$

Daarmee is de decompositie van $A$ in $LU$ beëindigd. In dit geval was er geen verwisseling van rijen nodig.

Opmerkingen

Om een stelsel van lineaire vergelijkingen $Ax=B$ op te lossen gaat men als volgt te werk:

1. Vorm de matrices $L$ en $U$: $LUx=B$. We noemen $y$ het (voorlopig onbekende) product $Ux$.
2. Los $Ly=B$ op door voorwaartse substitutie.
3. Los $Ux=y$ op door achterwaartse substitutie.

Als bijvoorbeeld $B={\left(\begin{array}{cccc}-36& 20& 2& -34\end{array}\right)}^T,$ vinden we achtereenvolgens:

• $3y_1=-36$ of $y_1=-12$
• $-3y_1-2y_2=20$ waaruit $y_2=8$
• $6y_1+10y_2-y_3=2$ waaruit $y_3=6$
• $-9y_1-16y_2-3y_3-y_4=-34$ waaruit $y_4=-4$

en:

• $x_4=y_4=-4$
• $x_3-x_4=y_3=6$ waaruit $x_3=2$
• $x_2+\frac{x_3}{2}-x_4=y_2=8$ waaruit $x_2=3$
• $x_1-\frac{7x_2}{3}-\frac{2x_3}{3}+\frac{2x_4}{3}=y_1=-12$ waaruit $x_1=-1$

Indien deze methode in een computerprogramma wordt geïmplementeerd, is het niet nodig om twee afzonderlijke matrices voor $L$ en $U$ te gebruiken. Ze kunnen compact in een enkele matrix bewaard worden (daarvoor kan zelfs de matrix $A$ dienen als deze mag overschreden worden). Het is immers niet nodig om de elementen op de hoofddiagonaal van $U$ te bewaren, die per definitie gelijk zijn aan 1. De rij-permutaties kunnen bijgehouden worden in een permutatievector van $n$ elementen 1 tot en met $n;$ als twee rijen verwisseld worden worden de corresponderende elementen in de permutatievector verwisseld.

De Crout-decompositie gelijkt sterk op de Doolittle-decompositie. Dit is ook een LU-decompositie met dit verschil dat in een Doolittle-decompositie de elementen op de hoofddiagonaal van de benedendriehoeksmatrix $L$ gelijk zijn aan 1.

Websites

• Cholesky, Doolittle and Crout Factorization