Esercizio in Matlab su sistemi lineari e fattorizzazione

[Alessio_Ale]

Un altro esercizio che non riesco a risolvere.





Il codice:

n=18;

A1=6*ones(n,1);    %vettori per la creazione di A
A2=3*ones(n-1,1);
A3=6*ones(n-1,1);

    A=diag(A1)+diag(A2,1)+diag(A3,-1); %creo A
    
    b=linspace(5,8,n);

    tic %calcolo costo computazionale  LUP
    [L,U,P]= lu(A);
    y=L\(P*b');
    x=U\y;
    
    t1=toc

    tic %calcolo costo computazionale  QR
    [Q,R]=qr(A);
    x=R\(Q*b');
    
    t2=toc

Non riesco ad arrivare al risultato, ma prima ancora mi rendo conto di avere dei valori stranissimi derivanti sicuramente a qualche errore nel codice. Riuscite ad aiutarmi a correggerlo? Grazie!

[feddy]

Ciao Alessio,

le prossime volte inserisci il codice tra gli appositi tag, come ho fatto io nel tuo post, che ora risulta modificato. Aiuta la lettura.

[feddy]

Non hai calcolato gli autovalori. Se questi sono positivi, allora puoi applicare la fattorizzazione di Cholesky, in quanto la matrice è SPD (simmetrica e definita positiva). Il fatto che sia simmetrica è evidente da come è definita A.

Qualora non fosse SPD, allora va usato un altro metodo.

Il codice che hai scritto è molto breve, l'errore lo puoi trovare in pochissimo tempo. Inizia intanto a verificare che la tua incognita *risolva* il sistema lineare, cioè che $Ax-b$ sia piccolo in una qualche norma...

[feddy]

Siccome a quest'ora non ho nulla da fare...

Errori:
1) La matrice $A$ che hai creato non è quella richiesta. In particolare, la tua non è nemmeno simmetrica. Infatti,

A3

deve essere moltiplicata per $3$, invece che per $6$.

A dire il vero, non ti serve nemmeno definire

A3

come variabile, perché è uguale alla diagonale sopra. Dunque, ti basta un solo vettore per le sopra e sotto diagonali.

2) La tua risoluzione con LU è corretta, ma la QR no. Dalla teoria sai che $x=R^{-1}(Q^t b)$, ma nel tuo codice tu stai usando $Q$, non $Q^t$. Basta aggiungere un apice e risulta corretto.

3) Devi calcolare gli autovalori, notare che sono tutti positivi e, siccome la matrice è anche simmetrica, utilizzare l'algoritmo di Cholesky per fattorizzare $A$ in $A=R^tR$. Questi sono due sistemi triangolari, proprio quelli richiesti dall'esercizio.





Di seguito l'esercizio corretto: ho definito

b

direttamente come vettore colonna, per evitare il trasposto. Ho corretto la risoluzione con $QR$, che di fatto non era richiesta, e aggiunto la parte che mancava riguardo agli autovalori.

n=18;
d=6*ones(n,1);    %vettori per la creazione di A
sd=3*ones(n-1,1);
A=diag(d)+diag(sd,1)+diag(sd,-1); %creo A
b=linspace(5,8,n)';

tic %calcolo costo computazionale  LUP
[L,U,P]= lu(A);
y=L\(P*b);
x=U\y;
toc

tic %calcolo costo computazionale  QR
[Q,R]=qr(A);
x=R\(Q'*b);
toc

% Autovalori
eig(A)>0 %risulta un vettore di 1, cioè true, dunque la matrice è definita positiva

% Poichè sappiamo quindi che A è simmetrica e definita positiva, allora possiamo utilizzare
% tranquillamente la decomposizione di Cholesky:
tic
R = chol(A);
y = R'\b;
x = R\y;
toc

norm(x+y,1)

[Alessio_Ale]

Ciao! Ho visto che hai risposto a tutti i post che ho scritto e inizio ringraziandoti per la correzione e ancora più per avermi spiegato gli errori.
Prendo inoltre i consigli per scrivere i futuri post nel modo migliore!

[feddy]

Prego!