LUdecomposition

Il codice per la fattorizzazione LU della matrice $A$ segue la strategia ricorsiva descritta nel paragrafo , comunque, per ottimizzare la funzione, è stata sostituita la ricorsione con un ciclo iterativo. Poiché la matrice di output $L$ ha degli zero sopra la diagonale e ha dei valori uno nella diagonale, la funzione non si occuperà di riempire quelle posizioni. Analogamente, poiché la matrice $U$ ha tutti gli elementi sotto la diagonale uguali a zero neanche queste posizioni saranno riempite. Il codice calcola, pertanto, soltanto gli elementi significativi delle matrici $L$ ed $U$.

void LUdecomposition(float **A, int dim) {
  int k, i, j;
1   for (k = 0; k < dim; k++) {
2     for (i = k + 1; i < dim; i++) {
3       A[i][k] = A[i][k] / A[k][k];
4     }
5     for (i = k + 1; i < dim; i++)
6       for (j = k + 1; j < dim; j++)
7         A[i][j] = A[i][j] - A[i][k] * A[k][j];
8   }
9 }

La funzione, per ottimizzare la memoria utilizzata, memorizza le matrici $L$ ed $U$ nella matrice $A$ stessa, terminata la funzione l' elemento $a_{i,j}$ apparterrà a $L$ se $i > j$, oppure ad $U$ se $i \leq j$.

Il ciclo for esterno che comincia dalla riga 1, si ripete per ogni passo ricorsivo. Dentro il ciclo for alle linee 2-4, che non vengono eseguite quando $k = dim$, si pone in $a_{i,k}$ il valore $v_i$ (linea 3). L' elemento $w^T_i$ sarà invece il valore in $a_{k,i}$. Infine la sottomatrice $A'$ verrà memorizzata nella matrice $A$ stessa nelle linee 5-8. LUdecomposition() impiega tempo $\Theta(\frac{1}{3} n^3)$.