Metodi diretti: generalitàriccardo.broglia/Didattica/AN2004-05... · Metodo di sostituzione: costo...
28
Metodi diretti: generalità • Sono basati su una trasformazione del sistema algebrico iniziale in uno equivalente dalla struttura più semplice. • La soluzione è ottenuta in un numero finito di passi ed in assenza di errori di arrotondamento si otterrebbe la soluzione esatta. • Sono applicati a problemi “piccoli” e “densi”. • Efficienza stimata in termini di costo computazionale. B X A B AX ~ ~ = ⇔ =
Transcript of Metodi diretti: generalitàriccardo.broglia/Didattica/AN2004-05... · Metodo di sostituzione: costo...
Metodi diretti: generalità
• Sono basati su una trasformazione del sistema algebrico iniziale in uno equivalente dalla struttura più semplice.
• La soluzione è ottenuta in un numero finito di passi ed in assenza di errori di arrotondamento si otterrebbe lasoluzione esatta.
• Sono applicati a problemi “piccoli” e “densi”.• Efficienza stimata in termini di costo computazionale.
BXABAX ~~=⇔=
Metodi diretti: triangolare superiore• Sistema algebrico triangolare superiore regolare:
• Algoritmo di sostituzione all’indietro:
=
=++
=++++=+++++
⇔=
nnnn
knknkkk
nnkk
nnkk
bxu
bxuxu
bxuxuxubxuxuxuxu
LLLL
LLLLLLLLL
222222
111212111
BUX
( )
−=
−=
=
∑+=
1,,111
Knkxubu
x
ubx
n
kiikik
kkk
nn
nn
Metodi diretti: triangolare inferiore• Sistema algebrico triangolare inferiore regolare:
• Algoritmo di sostituzione in avanti:
=+++
=+++
=+=
⇔=
nnnnnn
kkkkkk
bxlxlxl
bxlxlxl
bxlxlbxl
LLLLLLLL
LLLL
2211
2211
2222121
1111
BLX
=
−=
=
∑−
=
nkxlbl
x
lbx
k
iikik
kkk ,,21 1
1
11
11
K
Metodo di sostituzione: costo• Ad ogni passo del metodo di sostituzione in avanti o
all’indietro ci sono:– n-k moltiplicazioni– 1 divisione
• Costo computazionale:
( ) flops22
122
2
1
nnnnknCn
kc ≈+−≈+−= ∑
=
• Perché è richiesto che U e L siano regolari?
Metodi diretti: eliminazione di Gauss• Trasforma il sistema algebrico originale in un sistema
equivalente triangolare superiore:
• Operazioni “lecite”, cioè che trasformano il sistema originale in uno equivalente:
– scambio tra equazioni– sostituzione di una equazione con una combinazioni
lineare tra la stessa e altre equazioni
• Il metodo opera sulle righe della matrice completa:
BUXBAX ~=⇔=
[ ]
=
nnnnn
n
n
b
bb
aaa
aaaaaa
ML
MLMMLL
2
1
21
22221
11211
B|A
Metodo di eliminazione di Gauss• Il metodo può essere così descritto: detta [A(1)|B(1)] la
matrice completa originale, al passo (k+1) la matrice [A(k+1)|B(k+1)] è ottenuta dalla matrice [A(k)|B(k)] annullando (per mezzo di operazioni “lecite”) tutti gli elementi sottodiagonali dalla colonna 1 alla colonna k.
• Al passo n-esimo il sistema è tridiagonale superiore:
• Il sistema A(n)X=B(n) può essere risolto tramite l’algoritmo delle sostituzioni successive.
[ ]
=
)(
)2(2
)1(1
)(
)2(2
)2(22
)1(1
)1(12
)1(11
00
0
nn
nnn
n
n
b
bb
a
aaaaa
M
L
MLMM
L
L
B|A
[ ]
=
)1(4
)1(3
)1(2
)1(1
)1(44
)1(43
)1(42
)1(41
)1(34
)1(33
)1(32
)1(31
)1(24
)1(23
)1(22
)1(21
)1(14
)1(13
)1(12
)1(11
11
bbbb
aaaaaaaaaaaaaaaa
)()( B|A
Eliminazione di Gauss: esempio con n=4; k=1
==
njni
,1,2
[ ]
=
)2(4
)2(3
)2(2
)1(1
)2(44
)2(43
)2(42
)2(34
)2(33
)2(32
)2(24
)2(23
)2(22
)1(14
)1(13
)1(12
)1(11
22
bbbb
aaaaaaaaaaaaa
)()( B|A
−=−=
)1(11
)1()2(
)1(11
)1()2(
bmbbamaa
iii
jiijij)1(
11
)1(1
1 aam i
i =con:per:
000
⇓
⇓
Eliminazione di Gauss: esempio con n=4; k=2
==
njni
,2,3
[ ]
=
)3(4
)3(3
)2(2
)1(1
)3(44
)3(43
)3(34
)3(33
)2(24
)2(23
)2(22
)1(14
)1(13
)1(12
)1(11
33
000
bbbb
aaaaaaaaaaa
)()( B|A
−=−=
)2(22
)2()3(
)2(22
)2()3(
bmbbamaa
iii
jiijij)2(
22
)2(2
2 aam i
i =con:per:
00
⇓
⇓
[ ]
=
)2(4
)2(3
)2(2
)1(1
)2(44
)2(43
)2(42
)2(34
)2(33
)2(32
)2(24
)2(23
)2(22
)1(14
)1(13
)1(12
)1(11
22
000
bbbb
aaaaaaaaaaaaa
)()( B|A
Eliminazione di Gauss: esempio con n=4; k=3
==
njni
,3,4
[ ]
=
)4(4
)3(3
)2(2
)1(1
)4(44
)3(34
)3(33
)2(24
)2(23
)2(22
)1(14
)1(13
)1(12
)1(11
44
0000
0
bbbb
aaaaaaaaaa
)()( B|A
−=−=
)3(33
)3()4(
)3(33
)3()4(
bmbbamaa
iii
jiijij)3(
33
)3(3
3 aam i
i =con:per:
⇓
⇓
[ ]
=
)3(4
)3(3
)2(2
)1(1
)3(44
)3(43
)3(34
)3(33
)2(24
)2(23
)2(22
)1(14
)1(13
)1(12
)1(11
33
0000
0
bbbb
aaaaaaaaaaa
)()( B|A
0
Eliminazione di Gauss: esercizio
[ ]
−−−
−−−
=
1224
041044/32/31
52/1322524
B|A
• Risolvere il sistema seguente tramite il metodo di eliminazione di Gauss:
• Soluzione:
[ ]
−−
−−−
=
10/21104
20/510002/1500
43202524
)4()4( B|A
Eliminazione di Gauss: applicabilità
−=−=
+
+
)()()1(
)()()1(
kkik
ki
ki
kkjik
kij
kij
bmbbamaa
=+=nkj
nki,
,1per:
• Le quantità: mik=aik(k)/akk
(k) sono detti moltiplicatori;• gli elementi akk
(k) sono detti elementi pivotali;• le condizioni akk
(k)≠0 garantiscono l’applicabilità del metodo, equivalgono alle condizioni det(Ak)≠0.
• se akk(k)=0 , poiché la matrice A è regolare (altrimenti il
sistema sarebbe irrisolvibile), esiste sempre la possibilità di uno scambio di righe per portare un elemento non nullo nella posizione diagonale (pivoting parziale).
1,1 −= nkPer: )(
)(
kkk
kik
ik aam =con:
Algoritmo
==
ii
ijij
bbaa
)1(
)1(
nji ,1, =per:
Eliminazione di Gauss: pivoting, scaling
• Pivoting parziale: akk(k)≠0 ma piccolo in modulo, per la stabilità del
metodo può comunque essere utile portare in diagonale un terminemaggiore in modulo, individuato quell’indice r≥k tale che:
si scambiano le righe r e k.
• Pivoting totale: individuata quella coppia di indici r,q ≥k tali che:
si scambiano le righe di indici r e k e le colonne di indici q e k.
• Scaling e pivoting scalato: divisione di ciascuna riga per per la norma del corrispondente vettore riga (norma uno o infinito):
)()( max ksknsk
krk aa
≤≤=
)(
,
)( max kspnpsk
krq aa
≤≤=
( ) sssk
spnpskrk
rq aaaaa A==≤≤
/max/ )(
,
)(
• Al k-esimo passo del metodo di Gauss sono necessarie:– n-k divisioni per il calcolo dei moltiplicatori– (n-k)(n-k+1) moltiplicazioni (per il calcolo di [aij
(k)])– (n-k) moltiplicazioni (per il calcolo di [bi
(k)])
• Costo computazionale:
Eliminazione di Gauss: costo
( ) ( ) ( )[ ] flops3
1)()(31
1GAUSS
nknknknknCn
kc ≈−++−−+−= ∑
−
=
• A cui vanno aggiunte le operazioni per la soluzione del sistema triangolare:
flops23
23 nnCc +=
• Al k-esimo passo del metodo di Gauss sono azzerati gli elementi sotto la diagonale principale, nel metodo di riduzione di Gauss-Jordan vengono annullati anche gli elementi sopra la diagonale principale.
• Il problema equivalente è diagonale• Costo computazionale:
Metodi diretti: riduzione di Gauss-Jordan
−=−=
+
+
)()()1(
)()()1(
kkik
ki
ki
kkjik
kij
kij
bmbbamaa
=≠=
nkjki
ni
,
,1per:1,1 −= nkPer: )(
)(
kkk
kik
ik aam =con:
flops2
3
nnCc +≈
Algoritmo
Esempio 4.9.1: risolvere il seguente sistema con il metodo di eleminazione Gauss con e senza pivoting parziale operando con una risoluzione a quattro cifre decimali:
Soluzione esatta: X= [10.0,1.0]T
Soluzione senza pivoting parziale: X= [6.667,1.001]T
Soluzione con pivoting parziale: X= [10.0,1.0]T
Eliminazione di Gauss: esempio
=−=+
4640.13160.23780.03690.13660.10003.0
21
21
xxxx
• Esercizi consigliati [GL] 2.12, 2.13
• Il metodo di eliminazione di Gauss può essere visto come una procedura per la fattorizzazione della matrice A in due matrici triangolari, inferiore con elementi diagonali unitari e superiore:
• Con L matrice dei moltiplicatori e U matrice finale dell’al-goritmo di Gauss.
Metodi diretti: fattorizzazione LU
==
)(
)3(3
)3(33
)2(2
)2(23
)2(22
)1(1
)1(13
)1(12
)1(11
321
3231
21
000
000
1
010010001
nnn
n
n
n
nnn a
aaaaaaaaa
mmm
mmm
L
MOMMM
L
L
L
L
MOMMM
L
L
L
LUA
• Matrice di trasformazione elementare Eij(c):
Metodi diretti: fattorizzazione LU
+
=
00
00
00
1000
010000100001
)(
LLL
MLLM
LLL
MLLLM
LLL
L
MOMMM
L
L
L
ccijE
• Moltiplicando a sinistra una matrice A per Eij(c) equivale a sommare alla riga i-esima (di A) la j-esima riga moltiplicata per c.
i
j
• Metodo di eliminazione di Gauss in forma matriciale:
Metodi diretti: fattorizzazione LU
)(
)(,1,1,2,2,,
)1( )()()(k
k
kkkkkkkkkknkn
k mmm
AL
AEEEA
=
=−−−= +++++ L
−
−=−=
+
+
=∏
100
1
010001
)(
,
,1
1
,,
LL
MLMLMM
MOLMM
MLLMM
MLMOMM
LML
LML
kn
kk
k
nikikik
m
mmEL
−=−=
+
+
)()()1(
)()()1(
kkik
ki
ki
kkjik
kij
kij
bmbbamaa
=+=nkj
nki,
,1per:1,1 −= nkPer: )(
)(
kkk
kik
ik aam =con:
• Metodo di eliminazione di Gauss in forma matriciale:
Metodi diretti: fattorizzazione LU
( ) LUAULULA
ALUALAL
ALLLALLALU
=⇔=
=
=⇔==
====
−−
=
−
−=
−=−=−=
−−−
−−−
−
∏∏
∏∏∏11
1
11
1
1
1
1
1
)1(1
1
)1(121
)2(21
)1(1
n
kk
nkk
nkk
nkk
nkk
nnn
nnn
n L
( ) ( )
=
= ∏−
=
−
+
−
1
010010001
100
1
010001
321
3231
211
1
1
,
,1
1
L
MOMMM
L
L
L
LL
MLMLMM
MOLMM
MLLMM
MLMOMM
LML
LML
nnn
n
kk
kn
kk
k
mmm
mmm
m
mLL
• Matrici di permutazione Pij di ordine n, matrice di identità di ordine con le righe i e j scambiate:
Metodi diretti: fattorizzazione LU
=
0001001000001001000000001
25P
• Moltiplicando a sinistra una matrice A per Pij equivale a scambiare la riga i-esima con la riga j-esima (di A).
2
5
1)det( −=ijP
Teorema: se A∈ℜnxn è tale che det(Ak)≠0, per k=1…n, ossia se tutti i suoi determinanti principali di testa sono non nulli,allora esistono L,U∈ℜnxn, matrici triangolare inferiore con elementi diagonali unitari e triangolare superiore tali che:
Teorema: se A∈ℜnxn è regolare, ma non soddisfa la condi-zione det(Ak)≠0, per k=1…n, esiste una matrice di permu-tazione P di ordine n, tale che la matrice: soddisfa le condizioni: det( )≠0, per k=1…n.
Nota: il costo computazionale della fattorizzazione LU è pari a quello del metodo di Gauss.
Fattorizzazione LU: applicabilità
LUA =
PAA =~
kA~
• Soluzione di un sistema algebrico lineare:
con un costo computazionale pari a n3/3+2n2/2.• Soluzione di più sistemi algebrici al variare del t.n.:
con un costo pari a n3/3+Nn2 operazioni floating point, con Gauss si avrebbe un costo totale di N(n3/3+n2/2)
Fattorizzazione LU: applicazioni
==
⇔=⇔== YUX
BLYBLUXBAX
LUA
Sostituzioni in avanti
Sostituzioni all’indietro
NiNiNi ii
iiiiii ,1,1,1
=
==
⇔
==
⇔
==
= YUXBLYBLUXBAX
LUA
• Calcolo del determinante:
dove con s si è indicato il numero di scambi di righe.
• Calcolo dell’inversa: detti Xi, i=1…n i vettori colonna di A-1:
Fattorizzazione LU: applicazioni
)det()det()det()det()det( ULLUAP ==
ni
ni
ii
ii
iinn
,,1
,,1][][ 2121
1
K
KKK
=
==
⇔
==
⇔=⇔=
=
−
YUXELY
EAXEEEXXXAIAA
LUA
( )( ) )()2(
22)1(
1111
111)det( n
nns
n
kkk
n
kkks aaaul L−=
−
= ∏∏==
A
Calcolo del rango di una matrice:• se l’algoritmo di eliminazione applicato alla matrice quadrata A∈ℜnxn
termina regolarmente dopo n-1 passi, allora la matrice A ha rango massimo pari ad n (la matrice è regolare).
• se l’algoritmo di eliminazione applicato alla matrice quadrata A∈ℜnxn
termina dopo q passi poiché ark(k)=0, r=k,…,n, allora la matrice A ha
rango q≤n (la matrice è singolare).
• l’algoritmo di eliminazione applicato alla matrice rettangolare A∈ℜmxn
termina necessariamente dopo q passi o perché q=m o perché tutti gli elementi dalla riga q+1-esima alla riga m sono nulli; la matrice A ha rango q.
Fattorizzazione LU: applicazioni
Esempio 4.10.2: fattorizzare la matrice
Fattorizzazione LU: esempio
−−
−−
=
1001121100211121
A
• Esercizi consigliati [GL] 2.16,2.17
• Banachiewicz-Doolittle: permette di fattorizzare la matrice A senza utilizzare il metodo di eliminazione di Gauss (stesso Cc≈ n3/3):
Formule compatte di fattorizzazione
≤<≤
−=
≤≤≤−=
∑
∑−
=
−
=
nijuulal
nijulau
jj
j
kkjikijij
j
kikjkjiji
1/
1
1
1
1
1
• Cholesky: se la matrice A è definita positiva può essere fattorizzata in due matrici una trasposta dall’altra: A=LLT con L triangolare inferiore (costo computazionale Cc≈ n3/6):
−=
≥>≥
−=
∑
∑−
=
−
=
1
1
2
1
11/
j
kjkjjjj
jj
j
kjkikijij
lal
jinlllal
• Esercizio consigliato [GL] 7.24
• Sistema sovradeterminato: sistema algebrico lineare di mequazioni in n incognite con m>n;.
Sistemi sovradeterminati
• Soluzione cercata: determinare la soluzione X che minimizza una norma del vettore dei residui R=AX-B∈ ℜm, ad esempio la norma 2:
BAX = con A∈ ℜmxn, B∈ ℜm e X∈ ℜn
In generale la soluzione non esiste; ciò accade per esempio nell’interpolazione ai minimi quadrati (retta “passante” per più di due punti).
2
1 1
2
2
2
2 ∑ ∑= =
−=−=
m
i
n
kikik bxaBAXR
• Teorema: Se la matrice dei coefficienti A ha rango n (numero delle incognite), esiste un solo vettore X che minimizza il residuo in norma 2, ed è dato dalla soluzione del sistema algebrico quadrato:
Sistemi sovradeterminati
BAAXA TT =
• Il vettore X che minimizza (||R||2)2 è detto soluzione del sistema nel senso dei minimi quadrati.
