PA081: Programování numerických
výpočtů
PA081: Programování numerických výpočtů
12. Automatické vyhodnocení derivací
jaro 2011
1/16
Automatické derivování
Obrácené vyhodnocení
Motivace
► kromě funkce dokážeme vyhodnotit i její derivace
► první, druhé, ...
► parciální
► mnohé numerické problémy jsou lépe řešitelné metodami s dostupnou derivací
► metody s derivací zpravidla rychleji konvergují
► Newton vs. metoda sečen
Motivace
► kromě funkce dokážeme vyhodnotit i její derivace
► první, druhé, ...
► parciální
► mnohé numerické problémy jsou lépe řešitelné metodami s dostupnou derivací
► metody s derivací zpravidla rychleji konvergují
► Newton vs. metoda sečen
aproximace konečnou diferencí nestačí
► zpracovávaná funkce nebývá triviální
► analytické vyjádření derivace je náročné
► zdroj nebezpečných chyb
► komplikace při změnách
Automatické derivování
► funkce je implementována jako program
► hledáme program, který vedle funkční hodnoty spočte i derivace
► http://www.autodiff.org
Automatické derivování
► funkce je implementována jako program
► hledáme program, který vedle funkční hodnoty spočte i derivace
► http://www.autodiff.org
► příklad s majákem
	II!                       & / ď/ ^^Jriti = 7Ifl
*--v—	
► bod dopadu paprsku na nábřeží v tan cot
yi = —;-r    j2 = yyi
y - tan cot
Automatické derivování
► vyjádřeno jako jednoduchý program
a := tan cot;      yx:= va/(y - a); y2-=yyi\
► chceme spočítat i rychlost pohybu, tj. derivaci v čase
► výpočet současně s původní funkcí
t := 1
a := tanai a.2 := y - a yi := va/a.2 y2 ■■= yyi
a\ := cot
ái := cot
ä := ái/(cosai)2
Č12 '■= —á.
ýi := v(äa.2 - aä2)Ia\ y2 ■= yýi
► pouze mechanická aplikace základních pravidel
Automatické derivování
► transformace zdrojového kódu (např. ADIC)
► přetížení operátorů a funkcí (ADOL-C)
► zjednodušený příklad
Automatické derivování
Obrácené vyhodnocení
class adouble { double x; double dx;
}
adouble operátor * (adouble a, adouble b) { adouble r; r.x = a.x * b.x; r.dx = a.x * b.dx + a.dx * b.x; return r;
}
adouble vl, v2, v3; vl = v2 * v3:
Automatické derivování
PA081: Programování numerických
výpočtů
Automatické derivování
► aktivní proměnná
Obrácené vyhodnocení
► nezávislé proměnné, podle kterých chceme derivovat
► všechny další proměnné na nich během výpočtu závisející
► vyhrazený typ adoubl e
► přiřazení adoubl e ->• doubí e je chybné
► ztrácíme informaci o derivacích tam, kde je potřebná
► explicitní přístup metodou value()
► více nezávislých proměnných
► s každou aktivní proměnnou udržujeme vektor parciálních
derivací
6/16
Typický program
#define NUMBER_DIRECTIONS 3 #define ADOLC_TAPELESS #inc"lude <adol c/adoubl e. h>
func(adouble[] , const adouble[]);
adoubl e x [3] ; adoubl e y [4] ;
for (i=0; i<3; i++)
for (j=0; i<3; i++)
x [i] .setADValue(j ,i==j ? 1.0 : 0.0);
func(y,x);
for (i=0; i<4; i++) {
cout « y.getValueO « ": "; for (j=0; j<3; j++) cout « y.getADValue(j) « cout « endl;
}
Automatické derivovaní
Obrácené vyhodnocení
Není to tak jednoduché
► nediferencovatelné funkce (fmod)
► nelze použít
► absolutní hodnota (f abs)
► v 0 pravá nebo levá derivace ± 1
► případně položíme uměle rovnu 0
► minimum a maximum (f mi n, f max)
► lze nahradit (a + b ± \a- b\)/2
► nedefinované na celém R. (sqrt, pow, tan, ...)
► derivaci položíme uměle NaN, resp. Inf
PA081: Programování numerických
výpočtů
Automatické derivování
Obrácené vyhodnocení
8/16
Obrácené vyhodnocení
► množství potřebných operací 0(mk) pro /: [R
► m počet nezávislých proměnných
► n počet závislých proměnných
► k počet kroků výpočtu
► nepříjemné s rostoucím m
► velká třída problémů velkým m ale malým n
► optimalizace funkce více proměnných (n = 1)
Obrácené vyhodnocení
► množství potřebných operací 0(mk) pro /: Rm ->• Rn
► m počet nezávislých proměnných
► n počet závislých proměnných
► k počet kroků výpočtu
► nepříjemné s rostoucím m
► velká třída problémů velkým m ale malým n
► optimalizace funkce více proměnných (n = 1)
► počítali jsme citlivost pomocných proměnných na změny vstupu
► obrácený postup, citlivost výstupu (závislých proměnných) na změny pomocných
► výpočet od konce
PA081: Programování numerických
výpočtů
Automatické derivování
Obrácené vyhodnocení
9/16
Obrácené vyhodnocení
PA081: Programování numerických
výpočtů
► např. funkce y = sin(%i
► označujeme ä = d y Ida
/x2)eX3
Automatické derivování
Obrácené vyhodnocení
a.\ := x\\x2 a2 '■= sinai
a3 := e*3
y := a2a3 y:= 1 ä2 := a3y a-i := a2y x3 := eX3ä3 ä\ := ä2 cosai x2 := -ä\lx2 x\ := ä2/x2
► v O(nk) operacích máme všechny dy/dxt
10/16
Stopa výpočtu
derivace se vyhodnocují od konce
► tj. až poté, co proběhl výpočet funkce
► proto „obrácené" (reverse) vyhodnocení
► výpočet funkce je nutné nějak zaznamenat
► ADOL-C: datová struktura (i soubor) páska (tape)
► záznam posloupnosti operací, nikoli konkrétních hodnot
► lze recyklovat pro jiný vstup
► případné rozdílnosti vyhodnocení (a>b ? c : d) jsou detekovány
► použije se pro vyhodnocení funkce, gradientu/Jakobiánu, Hessiánu, ...
Instrumentace pro zpětné vyhodnocení
		
#include <adol c/adoubl e. h>		
void func(		
double px[4],       /* vstup */		
double py[4],       /* vystup */		
) {		
adouble x[4],y[4],aux;		
trace_on(0);	/* paska 0	V
for (i=0; i<4;          x[i] «= px[i];		
aux = x[l] * x[4] + exp(x[2]);	/* vlastni	vypočet */
y[4] = si n(aux);		
for (i=0; i<4; i++) y[i] »= py[i];		
trace off();		
}		
Automatické derivování
Obrácené vyhodnocení
12/16
Drivery
► ADOL-C poskytuje celou řadu „driver" funkcí
► vstupem je vždy páska a konkrétní vstupní hodnoty
► driver provede výpočet na základě záznamu na pásce
► funkční hodnota funkce Rn -* Rm
► function(short tag, int m, int n, double x[n],double y[m])
► gradient funkce Rn ->■ R.
► gradient(short tag, int n, double x[n],double g[n])
► Jakobián funkce Rn - Rm
► jacobian(short tag, int m, int n, double x[n],double J[m][n])
► Hessián (matice druhých derivací) funkce Rn ->■ R.
► hessian(short tag, int m, int n, double x[n],double H[n][n])
Drivery
► optimalizované varianty pro Newtonovu metodu
► řešení různých typů diferenciálních rovnic speciální varianty pro řídké Jakobiány a Hessiány
► tenzory vyšších derivací
► implicitně definované funkce
G(x,y) = 0
► obecné - plná kontrola nad zpracováním pásky
Rekapitulace obráceného vyhodnocení
► aktivní proměnné deklarovat jako adoubl e
► výpočet označit trace_on () ... trace_off ()
► zavolat instrumentovanou funkci jednou
► s typickými hodnotami vstupů
► vyprodukuje záznam výpočtu - pásku
► opakovaně volat potřebné drivery
► případně na různé vstupy
► pro jinou posloupnost výpočtu zavolat instrumentovanou funkci s jinou páskou
Další materiál
► rozsáhlá problematika, řada implementací
► včetně obsáhlé teorie
► http://www.autodiff.org
► přehled nástrojů, FAQ, ...
► https://projects.coi n-or.org/ADOL-C
► implementac a dokumentace ADOL-C
► včetně řady příkladů
► Griewank, Walter. Evaluating Derivatives. 2008