2.    Implementazione a blocchi della fattorizzazione LU



Abbiamo visto, all'Esercitazione n.2, come sia importante utilizzare il più possibile operazioni di tipo BLAS3;  questo significa riformulare gli algoritmi "a blocchi", cioè effettuare singole operazioni su blocchi (= sotto-matrici) e non su scalari. Per la fattorizzazione LU si può operare nel seguente modo:

-    si partiziona la matrice A nel modo seguente :        A_1_1    A_1_2
                                                                                 A_2_1    A_2_2

-    si calcola la fattorizzazione L_1_1 * U_1_1 = A_1_1    con l'algoritmo scalare;

-    si risolve il sistema con membro destro multiplo L_1_1 * U_1_2 = A_1_2 ;

-    si risolve il sistema con membro destro multiplo  (U_1_1)' * (L_2_1)' = (A_2_1)'

-    si calcola  S = A_2_2  -  L_2_1 * U_1_2

e si ottiene :

        A_1_1    A_1_2                        L_1_1            0                       I            0                             U_1_1        U_1_2
                                            =                                                *                                    *               
        A_2_1    A_2_2                        L_2_1            I                       0            S                                  0                I

Notare che se si calcola la fattorizzazione  L_2_2 * U_2_2  =  S  , si ottiene

        A_1_1    A_1_2                        L_1_1            0                       I            0                             U_1_1        U_1_2
                                            =                                                *                                    *               
        A_2_1    A_2_2                        L_2_1        L_2_2                  0            I                                  0             U_2_2

e dunque per calcolare la fattorizzazione LU di A è sufficiente ripetere iterativamente il procedimento sopra esposto, suddividendo S in blocchi, ecc....


Esercizio: implementare l'algoritmo di fattorizzazione LU a blocchi in linguaggio Matlab/Octave e sperimentarlo con una matrice a scelta.