56
5
Metodo di Newton
5.1
Generalitµ
a
Il metodo di Newton µe stato introdotto originariamente come un metodo di soluzione di
un sistema di equazioni non lineari
F (x) = 0;
in cui F : Rn ! Rn si assume continuamente di®erenziabile. Indichiamo con J(x) la
matrice Jacobiana di F , ossia:
J(x) =
µ
@Fi (x)
@xj
¶
;
i; j = 1; : : : ; n;
o scrivere:
in cui Fi (x) sono le componenti del vettore F (x). Assegnato xk 2 Rn , si puµ
F (xk + s) = F (xk ) + J(xk )s + °(xk ; s)
dove °(xk ; s)=ksk ! 0 per s ! 0.
Il metodo di Newton consiste nel determinare xk+1 = xk + sk scegliendo sk in modo da
annullare l'approssimazione lineare di F (xk + s), ossia imponendo:
F (xk ) + J(xk )sk = 0:
Se J(xk ) µe invertibile, si ottiene:
sk = ¡[J(xk )]¡1 F (xk );
e quindi il metodo di Newton per la soluzione del sistema F (x) = 0 diviene:
xk+1 = xk ¡ [J(xk )]¡1 F (xk ):
(5.1)
In un problema di minimizzazione non vincolata, il metodo di Newton si puµo interpretare
come un algoritmo per la risoluzione del sistema di n equazioni in n incognite rf (x) = 0,
ottenute imponendo che il gradiente di f si annulli. Se f µe convessa ciµ
o equivale a costruire
una successione di punti minimizzando a ogni passo un'approssimazione quadratica della
funzione.
Se xk µe un punto assegnato, si puµo scrivere:
1
f (xk + s) = f (xk ) + rf (xk )0 s + s0 r2 f (xk )s + ¯(xk ; s);
2
57
in cui ¯(xk ; s)=ksk2 ! 0 per s ! 0. Per valori su±cientemente piccoli di ksk si puµ
o allora
pensare di approssimare f (xk + s) con la funzione quadratica
1
qk (s) = f (xk ) + rf (xk )0 s + s0 r2 f (xk )s
2
e determinare il punto successivo
xk+1 = xk + sk
scegliendo sk in modo da minimizzare (ove possibile) la funzione qk (s) rispetto a s. Poich¶e
rqk (s) = rf (xk ) + r2 f (xk )s;
se si suppone che r2 f (xk ) sia de¯nita positiva, il punto di minimo di qk (s) sarµ
a dato da:
sk = ¡[r2 f (xk )]¡1 rf (xk ):
Il metodo di Newton µe allora de¯nito dall'iterazione
xk+1 = xk ¡ [r2 f (xk )]¡1 rf (xk ):
(5.2)
µ facile rendersi conto che, se si identi¯ca F (x) con il gradiente rf (x) di una funzione
E
f : Rn ! R, la matrice Jacobiana di F coincide con la matrice Hessiana di f e quindi la
(5.1) coincide con la (5.2).
Nel seguito riportiamo dapprima alcuni risultati fondamentali sulla convergenza locale
del metodo di Newton, riferendoci, per semplicitµ
a di notazioni, al caso di un sistema di
equazioni non lineari F (x) = 0. Successivamente, con riferimento ai problemi di minimizzazione non vincolata, accenneremo ad alcuni dei metodi proposti per assicurare la
convergenza globale del metodo di Newton.
5.2
Convergenza locale
Le proprietµ
a di convergenza locale del metodo di Newton per la soluzione di sistemi di
equazioni non lineari sono state oggetto di studi approfonditi.
In particolare, uno dei risultati piµ
u importanti sull'argomento, noto anche come teorema di
1
Newton-Kantorovich, stabilisce condizioni su±cienti che assicurano l'esistenza di soluzioni
dell'equazione F (x) = 0 su spazi di funzioni e fornisce una stima della regione di convergenza.
In quel che segue, tuttavia, ci limiteremo a caratterizzare le proprietµ
a di convergenza
o consente di sempli¯locale in Rn assumendo come ipotesi l'esistenza di soluzioni. Ciµ
care notevolmente lo studio della convergenza; d'altra parte l'esistenza di punti stazionari
nei problemi di minimizazione non vincolata viene usualmente dedotta sulla base delle
ipotesi che assicurano l'esistenza di un punto di minimo e delle condizioni di ottimalitµ
a.
Stabiliamo, in particolare, il risultato seguente.
1
cfr. Kantorovich, L. and Akilov, G. Functional Analysis in Normed Spaces, Pergamon Press, Oxford,
1964.
58
Proposizione 5.1 (Convergenza locale del metodo di Newton)
Sia F : Rn ! Rn continuamente di®erenziabile su un insieme aperto D µ Rn .
Supponiamo inoltre che valgano le condizioni seguenti:
(i) esiste un x? 2 D tale che F (x? ) = 0;
(ii) la matrice Jacobiana J(x? ) µe non singolare;
(iii) esiste una costante L > 0 tale che, per ogni x; y 2 D, si abbia
kJ(x) ¡ J(y)k · Lkx ¡ yk:
Allora esiste una sfera aperta B(x? ; ") ½ D, tale che, se x0 2 B(x? ; "), la successione fxk g
a di
generata dal metodo di Newton (5.1) rimane in B(x? ; ") e converge a x? con rapiditµ
convergenza quadratica.
Dim. Poich¶e J(x? ) µe non singolare e J(x) µe continua su D, µe possibile trovare un "1 > 0
e un ¹ > 0 tali che B(x? ; "1 ) µ D e che, per ogni x 2 B(x? ; "1 ) risulti:
°
°
°
¡1 °
°J(x) ° · ¹:
Sia ora
" < min ["1 ; 2=¹L]
e supponiamo che sia xk 2 B(x? ; ").
Essendo per ipotesi F (x? ) = 0, possiamo riscrivere la (5.1) nella forma:
xk+1 ¡ x? = ¡J(xk )¡1 [ ¡ J(xk )(xk ¡ x? ) + F (xk ) ¡ F (x? )];
da cui segue:
°
°
kxk+1 ¡ x? k · °J(xk )¡1 ° k¡J(xk )(xk ¡ x? ) + F (xk ) ¡ F (x? )k
· ¹ k¡J(xk )(xk ¡ x? ) + F (xk ) ¡ F (x? )k :
Poich¶e F µe di®erenziabile, si ha:
F (xk ) ¡ F (x? ) =
Z
1
0
J(x? + ¸(xk ¡ x? ))(xk ¡ x? )d¸;
e quindi, per la (5.3) si puµ
o scrivere:
°Z
°
kxk+1 ¡ x k · ¹ °
°
1
?
0
°
°
[J(x + ¸(xk ¡ x )) ¡ J(xk )](xk ¡ x )d¸°
°:
?
?
?
Dalla diseguaglianza precedente si ottiene:
?
kxk+1 ¡ x k · ¹
Z
0
1
kJ(x? + ¸(xk ¡ x? )) ¡ J(xk )k d¸ kxk ¡ x? k
59
(5.3)
e quindi, tenendo conto dell'ipotesi (iii) si ha:
kxk+1 ¡ x? k · ¹L
Z
0
1
(1 ¡ ¸)d¸kxk ¡ x? k2 =
¹L
kxk ¡ x? k2 :
2
(5.4)
Poich¶e si µe assunto xk 2 B(x? ; ") e si µe scelto " < 2=¹L, dalla (5.4) segue
kxk+1 ¡ x? k ·
¹L
¹L"
kxk ¡ x? k2 ·
kxk ¡ x? k < kxk ¡ x? k:
2
2
(5.5)
Ciµo implica xk+1 2 B(x? ; ") e di conseguenza, per induzione, si ha xk 2 B(x? ; ") per ogni
k. Applicando ripetutamente la (5.5) si ha anche:
?
kxk ¡ x k ·
µ
¹L"
2
¶k
kx0 ¡ x? k;
da cui segue, essendo ¹L"=2 < 1, che xk ! x? . La (5.4) implica allora che la rapiditµ
a di
convergenza µe quadratica.
Il risultato espresso dalla proposizione precedente si puµ
o facilmente riportare a un risultato
sulla convergenza del metodo di Newton nella minimizzazione di una funzione f : Rn ! R;
basta infatti tener presente che le ipotesi su F e J si traducono in ipotesi su rf e r2 f .
In particolare, si puµ
o enunciare la proposizione seguente, che µe una diretta conseguenza
dalla Proposizione 5.1
Proposizione 5.2 (Convergenza locale del metodo di Newton)
Sia f : Rn ! R una funzione due volte continuamente di®erenziabile su un insieme aperto
D µ Rn . Supponiamo inoltre che valgano le condizioni seguenti:
(i) esiste un x? 2 D tale che rf (x? ) = 0;
(ii) la matrice Hessiana r2 f (x? ) µe non singolare;
(iii) esiste una costante L > 0 tale che, per ogni x; y 2 D, si abbia
°
°
° 2
°
°r f (x) ¡ r2 f (y)° · Lkx ¡ yk:
Allora esiste una sfera aperta B(x? ; ") ½ D, tale che, se x0 2 B(x? ; "), la successione fxk g
a di
generata dal metodo di Newton (5.2) rimane in B(x? ; ") e converge a x? con rapiditµ
convergenza quadratica.
Si puµ
o osservare che il risultato precedente caratterizza la convergenza locale del metodo
di Newton nell'intorno di un qualsiasi punto stazionario in cui la matrice Hessiana sia non
singolare; si puµ
o trattare quindi, in particolare, sia di un minimo che di un massimo locale.
60
5.3
Modi¯che globalmente convergenti
Nell'applicazione del metodo di Newton alla minimizzazione non vincolata occorre tener
conto dei problemi seguenti:
(i) la direzione di Newton puµo non essere de¯nita in xk (r2 f (xk ) µe singolare);
(ii) la successione prodotta dal metodo di Newton puµ
o non essere convergente;
(iii) si puµ
o avere convergenza verso massimi locali.
Per superare tali di±coltµ
a si rende necessario modi¯care, con opportuni criteri, il metodo
di Newton. Le modi¯che devono tuttavia essere tali da non far perdere le caratteristiche
di rapiditµ
a di convergenza stabilite per il metodo di Newton nella sua forma \pura". Per
precisare meglio il problema introduciamo la seguente de¯nizione.
De¯nizione 5.1 (Modi¯ca globalmente convergente del metodo di Newton)
Sia f : Rn ! R due volte continuamente di®erenziabile e supponiamo che l'insieme di
livello
L0 = fx 2 Rn : f (x) · f (x0 )g
sia compatto.
Diremo che l'algoritmo de¯nito dall'iterazione:
xk+1 = xk + sk
µe una modi¯ca globalmente convergente del metodo di Newton se valgono le seguenti proprietµ
a:
(i) se fxk g µe in¯nita ogni punto di accumulazione di fxk g µe un punto stazionario di f
appartenente a L0 ;
(ii) nessun punto di accumulazione di fxk g µe un punto di massimo locale di f ;
(iii) se fxk g converge a un punto di minimo locale x? di f e r2 f (x) soddisfa le ipotesi
della Proposizione 5.2, esiste un k? tale che, per ogni k ¸ k ? la direzione sk coincide
con la direzione di Newton, ossia: sk = ¡[r2 f (xk )]¡1 rf (xk ).
Le modi¯che globalmente convergenti, note anche come metodi tipo-Newton (Newton type
methods) possono essere ricondotte, essenzialmente a tre classi di metodi.
(a) metodi di ricerca unidimensionale applicati a una direzione di ricerca tipo-Newton,
ottenuta modi¯cando la matrice Hessiana in modo da soddisfare opportuni requisiti
di discesa;
(b) metodi di tipo ibrido, basati sulla combinazione del metodo di Newton con un metodo
globalmente convergente, come il metodo del gradiente
(c) metodi tipo trust region (regione di con¯denza), in cui il punto xk+1 viene generato risolvendo un sottoproblema vincolato consistente nella minimizzazione di una
approssimazione quadratica di f in un intorno sferico di xk .
61
Prima di illustrare alcuni esempi signi¯cativi, premettiamo, senza dimostrazione, il seguente
risultato in cui si forniscono condizioni sotto cui il passo ® = 1 µe accettabile in un algoritmo
tipo-Armijo lungo una direzione dk .
Proposizione 5.3 Sia f : Rn ! R una funzione due volte continuamente di®erenziabile
e siano fxk g e fdk g due successioni tali che valgano le condizioni seguenti:
(i) fxk g converge a un punto x? in cui rf (x? ) = 0 e r2 f (x? ) µe de¯nita positiva;
^ risulti:
(ii) esistono un indice k^ e un numero ½ > 0 tali che, per ogni k ¸ k,
(iii) vale il limite:
rf (xk )0 dk · ¡½kdk k2 ;
(5.6)
kxk + dk ¡ x? k
= 0:
k!1
kxk ¡ x? k
(5.7)
lim
Allora, per ogni ° < 1=2 assegnato, esiste un indice k? , tale che, per ogni k ¸ k? si ha:
f (xk + dk ) · f (xk ) + °rf (xk )0 dk :
5.3.1
Modi¯che della matrice Hessiana e ricerche unidimensionali
Una modi¯ca globalmente convergente del metodo di Newton puµ
o essere ottenuta modi¯cando la matrice Hessiana, in modo tale che la direzione di ricerca sia una direzione
di discesa opportuna, lungo cui puµ
o essere e®ettuata una ricerca unidimensionale tipoArmijo.
Il criterio piµ
u semplice puµ
o essere quello di sommare alla matrice Hessiana un'opportuna
matrice de¯nita positiva Dk in modo tale che la matrice r2 f (xk ) + Dk risulti \ su±cientemente" de¯nita positiva. In tal caso si puµ
o assumere come direzione di ricerca una
direzione del tipo
h
dk = ¡ r2 f (xk ) + Dk
i¡1
rf (xk );
ed e®ettuare poi lungo dk una ricerca unidimensionale.
Un metodo conveniente per determinare Dk µe quello di utilizzare un procedimento di
fattorizzazione di Cholesky modi¯cata, mediante il quale vengono de¯nite una matrice
diagonale Dk e una matrice triangolare inferiore Lk con elementi diagonali positivi lii , tali
che:
r2 f (xk ) + Dk = Lk L0k :
µ da notare inoltre che,
Ciµo assicura che la matrice r2 f (xk ) + Dk sia de¯nita positiva. E
0
disponendo della fattorizzazione Lk Lk µe immediato risolvere il sistema lineare che fornisce
dk . Basta infatti determinare la soluzione s del sistema triangolare
Lk s = ¡rf (xk )
e successivamente ricavare dk risolvendo il sistema (anch'esso triangolare)
L0k dk = s:
62
Sono stati sviluppati procedimenti di fattorizazione in cui la determinazione di Lk viene
e®ettuata alterando il meno possibile la matrice Hessiana (ossia scegliendo gli elementi di
Dk abbastanza piccoli) e in modo da soddisfare anche le condizioni della De¯nizione 5.1.
Una descrizione sempli¯cata di un procedimento di questo tipo µe riportata nel seguito. Si
dimostra che se vale la (i), allora la direzione di Newton
dk = ¡[r2 f (xk )]¡1 rf (xk )
soddisfa le (ii) (iii) in un intorno di x¤ ; di conseguenza la proposizione precedente assicura
che il Metodo di Armijo accetta il passo unitario alla prima iterazione.
Come esempio di modi¯ca globalmente convergente del metodo di Newton consideriamo
l'algoritmo seguente.
Metodo di Newton modi¯cato (MNM)
Dati: q ¸ 3, p ¸ 2, c1 > 0, c2 > 0.
Inizializzazione: Si sceglie x0 2 Rn e si pone k = 0.
For k=0,1,. . .
1. Calcola rf (xk ); se rf (xk ) = 0 stop; altrimenti calcola r2 f (xk ).
2. Se il sistema: r2 f (xk )s = ¡rf (xk ) non ha soluzioni poni
dk = ¡rf (xk ) e vai al passo 5; altrimenti calcola una soluzione sN .
3. Se jrf (xk )0 sN j < c1 krf (xk )kq oppure ksN kp > c2 krf (xk )k poni
dk = ¡rf (xk ) e vai al passo 5.
4. Se rf (xk )0 sN < 0 assumi dk = sN ; altrimenti assumi dk = ¡sN .
5. E®ettua una ricerca unidimensionale con il metodo di Armijo, con ° < 1=2
(assumendo come stima iniziale ¢k = 1 se la direzione di Newton sN non µe
stata modi¯cata).
6. Poni xk+1 = xk + ®k dk , k = k + 1.
End for
µ facile veri¯care che le direzioni calcolate per mezzo dell'algoritmo precedente o coinciE
dono con l'antigradiente oppure soddisfano le condizioni
rf (xk )0 dk · ¡c1 krf (xk )kq ;
kdk kp · c2 krf (xk )k;
il che assicura che l'algoritmo di Armijo con ¢k = 1 sia convergente e che siano soddisfatte
le condizioni di convergenza globale introdotte in precedenza.
Osservazione La scelta dk = ¡sN e®ettuata al passo 4 se rf (xk )0 s > 0, µe motivata
dal fatto che la direzione dk = ¡sN , oltre a essere una direzione di discesa, µe anche una
direzione a curvatura negativa. Infatti, µe facile veri¯care che, nelle ipotesi poste:
d0k r2 f (xk )dk = rf (xk )0 [r2 f (xk )]¡1 rf (xk ) < 0:
Una direzione a curvatura negativa puµ
o essere particolarmente vantaggiosa in quanto µe
presumibile che, riducendosi lungo di essa la derivata direzionale, si possa ottenere una
63
notevole riduzione della funzione obiettivo. In alternativa, si potrebbe anche assumere
dk = ¡rf (xk ) quando rf (xk )0 sN > 0.
Si puµ
o dimostrare il risultato seguente.
Proposizione 5.4 Sia f : Rn ! R due volte continuamente di®erenziabile e si supponga che l'insieme di livello L0 sia compatto. Allora l'algoritmo (MNM) µe una modi¯ca
globalmente convergente del metodo di Newton nel senso della de¯nizione 5.1.
L'algoritmo precedente µe basato sull'uso della direzione dell'antigradiente nei punti in cui
la direzione di Newton non µe de¯nita o non soddisfa opportune condizioni.
Un diverso punto di vista nella de¯nizione di modi¯che globalmente convergenti del metodo
di Newton µe quello noto come metodo della \regione di con¯denza" (trust region) che µe
stato oggetto di un'intensa attivitµa di ricerca. Esso si basa sulla determinazione della
direzione e dell'ampiezza dello spostamento da e®ettuare a partire da xk in modo da
minimizzare l'approssimazione quadratica di f in una regione sferica di centro xk . Si
parla perciµo anche di metodo a passo ristretto (restricted step).
Per illustrare le motivazioni di questa impostazione, ricordiamo che il metodo di Newton per la minimizzazione non vincolata puµ
o essere interpretato come una sequenza di
minimizzazioni dell'approssimazione quadratica qk (s) di f de¯nita da:
1
qk (s) = f (xk ) + rf (xk )0 s + s0 r2 f (xk )s:
2
µ tuttavia evidente che, se r2 f (xk ) non µe almeno semide¯nita positiva, la funzione qk (s)
E
non ammette minimo. In tal caso, invece di fornire un diverso criterio per il calcolo della
direzione, si puµ
o pensare di e®ettuare la minimizzazione di qk (s) in un intorno limitato di
xk , ossia per ksk · ak , assumendo:
xk+1 = xk + sk ;
dove sk 2 Rn µe la soluzione del problema
min qk (s)
ksk · ak :
Negli algoritmi ¯nora proposti il raggio ak che de¯nisce la regione sferica attorno a xk
viene determinato in modo da assicurare che sia f (xk+1 ) < f (xk ) e che la riduzione di
f sia prossima a quella che si dovrebbe ottenere se f (xk + s) fosse e®ettivamente una
funzione quadratica. Ciµ
o quindi corrisponde a de¯nire la regione attorno a xk come la
regione in cui si puµ
o ritenere ragionevolmente valida l'approssimazione quadratica. Di qui
il nome di metodo della \regione di con¯denza".
5.4
Metodi di stabilizzazione non monotoni
Tutte le modi¯che globalmente convergenti del metodo di Newton ¯nora considerate si
basano, in modo essenziale, sulla generazione di una successione di punti a cui corrispondono valori monotonamente decrescenti della funzione obiettivo.
64
Ci si puµ
o chiedere tuttavia se tale proprietµ
a sia una caratteristica irrinunciabile e se non
comporti, in molti casi, delle conseguenze sfavorevoli sull'e±cienza del processo di minimizzazione.
In e®etti, l'esperienza di calcolo mette in evidenza che il metodo di Newton, nella sua
versione originaria, µe spesso notevolmente piµ
u e±ciente delle varie modi¯che globalmente
convergenti e l'analisi di tali casi rivela che la maggiore e±cienza µe legata al fatto che le
tecniche utilizzate per assicurare la convergenza globale comportano spesso una riduzione
del passo lungo la direzione di Newton, proprio per imporre una riduzione monotona
dell'obiettivo.
Un caso tipico µe quello in cui la funzione presenta delle \valli ripide" ad andamento curvilineo, per cui la riduzione monotona della funzione implica spostamenti di piccola entitµ
a
lungo le direzioni di ricerca e quindi un elevato e spesso inaccettabile costo computazionale.
Si consideri, ad esempio, la funzione (nota come funzione di Rosenbrock), de¯nita da:
f (x) = c (x2 ¡ x21 )2 + (1 ¡ x1 )2 ;
u ripida, quanto piµ
u elevato µe
che presenta una valle lungo la parabola x2 = x21 , tanto piµ
il coe±ciente c. Se si e®ettua la minimizzazione con un algoritmo standard globalmente
convergente (routine E04LBF della libreria NAG, che utilizza una tecnica di ricerca unidimensionale) si trova, ad esempio, che a partire dal punto iniziale x1 = ¡1:2, x2 = 1 e
per c = 106 sono richieste 358 iterazioni e 536 valutazioni della funzione e del gradiente
per ottenere un'approssimazione dell'ottimo in cui f = :6 ¢ 10¡6 . Il metodo di Newton con
passo unitario richiede invece solo 5 iterazioni e fornisce una stima in cui f = :2 ¢ 10¡27 .
Osservazioni analoghe si possono e®ettuare nella minimizzazione di molte altre funzioni
\di±cili" che presentano zone a forte curvatura delle super¯ci di livello. Funzioni di questo
tipo sono originate, tipicamente, dall'uso di tecniche di penalizzazione per la soluzione di
problemi vincolati.
Per superare tali di±coltµ
a si puµ
o pensare di indebolire i requisiti imposti sulla riduzione
della funzione, pur continuando ad assicurare le proprietµ
a di convergenza globale desiderate.
Tecniche di stabilizzazione di tipo non monotono si possono de¯nire, ad esempio, accettando che la funzione possa aumentare durante un numero di passi pre¯ssato.
In particolare, in luogo del metodo di Armijo, si puµ
o considerare l'algoritmo seguente,
che ne de¯nisce una versione non monotona, in cui M µe un intero pre¯ssato e si suppone
m(0) = 0 e m(k + 1) = minfm(k) + 1; M g.
65
Metodo di Armijo non monotono2
Dati: ± 2 (0; 1), ° 2 (0; 1=2), m(k) · minfk; M g.
Passo 1. Si pone ® = 1 e si de¯nisce il valore di riferimento:
Wk =
max ff (xk¡j )g:
0·j·m(k)
Passo 2. Se risulta
f (xk + ®dk ) · Wk + °®rf (xk )0 dk
si assume ®k = ®. Stop.
Passo 3. Si pone ® = ±® e si ritorna al Passo 2.
Notiamo che, al Passo 2, in luogo di imporre la condizione di su±ciente riduzione usuale
ci si riferisce al valore Wk de¯nito come massimo dei valori dell'obiettivo calcolati durante
(al piµ
u) M iterazioni precedenti. Ciµ
o consente quindi che f (xk+1 ) possa essere maggiore
u facilmente accettato il passo unitario.
di f (xk ), e, di conseguenza, che possa essere piµ
Si puµ
o dimostrare che utilizzando tale criterio, si puµ
o ottenere una modi¯ca globalmente
convergente del metodo di Newton che non impone una riduzione monotona della funzione.
Un fatto notevole µe che, pur ammettendo che la funzione possa aumentare da un'iterazione
all'altra, si puµ
o continuare ad evitare la convergenza verso massimi locali.
µ possibile anche considerare criteri piµ
E
u generali di stabilizzazione in cui, sotto opportune
condizioni, non si e®ettuano controlli sul valore dell'obiettivo durante un numero pre¯ssato
di iterazioni. L'esperienza di calcolo sembra indicare che le versioni non monotone delle
modi¯che globalmente convergenti del metodo di Newton assicurano notevoli vantaggi
computazionali rispetto alle tecniche di tipo monotono, sia in termini di riduzione del
numero di valutazioni della funzione, sia in termini di riduzione del numero di iterazioni.
66
5.5
Fattorizzazione di Cholesky modi¯cata*
Una delle modi¯che piµ
u note del metodo di Newton µe la fattorizzazione di Cholesky modi¯cata proposta da Gill e Murray (1972). Accenniamo qui a una possibile realizzazione
algoritmica.
Come si µe detto in precedenza, ci si propone di de¯nire una matrice diagonale Dk in modo
tale che la matrice r2 f (xk )+Dk risulti de¯nita positiva e sia quindi fattorizzabile secondo
Cholesky nella forma
r2 f (xk ) + Dk = Lk L0k :
Per sempli¯care le notazioni, poniamo A = r2 f (xk ):
Se A µe de¯nita positiva la fattorizzazione di Cholesky di A si puµo calcolare a partire
dall'equazione A = LL0 : Possiamo costruire L per colonne, osservando che per i ¸ j
dall'equazione precedente si ottiene
j
X
aij =
lik ljk ;
k=1
da cui segue, per j = 1; : : : ; n:
ljj
v
u
j¡1
X
u
2
= tajj ¡
ljk
;
k=1
lij =
aij ¡
j¡1
X
lik ljk
k=1
ljj
;
i = j + 1; : : : ; n:
Se A non µe de¯nita positiva le formule precedenti non sono applicabili, in quanto potrebbe
essere richiesto calcolare la radice quadrata di un numero negativo oppure dividere per
un elemento nullo. Possiamo allora introdurre una perturbazione sulla diagonale di A in
modo da poter ancora e®ettuare la fattorizzazione, il che implica che la matrice perturbata
sarµa de¯nita positiva.
A tale scopo introduciamo due scalari positivi ¹1 e ¹2 tali che 0 < ¹1 < ¹2 : Calcoliamo
allora la prima colonna di L ponendo
½p
a11 se a11 > ¹1
l11 = p
¹2 altrimenti.
Successivamente poniamo:
ai1
; i = 2; : : : ; n:
l11
Note le colonne 1; 2; : : : ; j ¡ 1 di L, possiamo calcolare la colonna j¡ma ponendo:
li1 =
ljj =
8v
u
>
>
u
>
<t
>
>
>
:p
ajj ¡
j¡1
X
2
ljk
k=1
lij =
aij ¡
k=1
ljj
2
ljk
> ¹1
k=1
altrimenti,
¹2
j¡1
X
se ajj ¡
j¡1
X
ljk lik
;
67
i = j + 1; : : : ; n:
Per quanto riguarda la scelta dei parametri ¹1 e ¹2 uno de criteri proposti
quello di calcolare la quantitµ
a
3
puµ
o essere
¯
¯
¯)
(¯
¯ @ 2 f (x ) ¯
¯ @ 2 f (x ) ¯
¯
¯
k ¯
k ¯
wk = max ¯
¯;:::¯
¯ ;
¯ (@xn )2 ¯
¯ (@x1 )2 ¯
e assumere ¹1 = r1 wk ; ¹2 = r2 wk ; in cui r1 e r2 sono scelti inizialmente a valori molto
piccoli e possono essere modi¯cati in base ai risultati dell'iterazione precedente.
µ da notare che valori eccessivamente piccoli di ¹2 tendono a rendere la matrice LL0 quasi
E
singolare e quindi possono far ottenere direzioni con kdk k molto grande. Il valore di ¹1 si
puµ
o far dipendere dalla norma del gradiente, in modo da evitare che il metodo di Newton
sia alterato in prossimitµ
a di un punto di minimo. Si puµ
o assumere, ad esempio,
r1 · r minf1; krf (xk )kg
con r su±cientemente piccolo, ad esempio r = 10¡6 .
Una volta determinata la fattorizzazione di A = r2 f (xk ) si tratta di risolvere il sistema:
LL0 d = ¡rf (xk ):
Ciµo puµ
o essere ottenuto risolvendo prima (per eliminazione in avanti) il sistema triangolare
Ls = ¡rf (xk );
e successivamente (per eliminazione all'indietro) il sistema triangolare
L0 d = s:
Le soluzioni dei due sistemi si ottengono assumendo:
s1 = ¡
@f (xk )=@xi +
si = ¡
lii
@f (xk )=@x1
;
l11
Pi¡1
k=1 lik sk
e calcolando successivamente
di =
n
X
lki dk
k=i+1
lii
i = 2; : : : ; n;
sn
;
lnn
dn =
si ¡
;
;
i = n ¡ 1; : : : ; 1:
La direzione di ricerca ottenuta attraverso la fattorizzazione di Cholesky modi¯cata sarµ
a
una direzione di discesa, lungo la quale µe possibile e®ettuare, ad esempio, una ricerca
unidimensionale con il metodo di Armijo non monotono.
3
Cfr. D.P. Bertsekas Nonlinear Programming, Athena Scienti¯c, 1995
68