Jyrkimmän laskeuman menetelmä

The method of Steepest Descent, gradienttimenetelmä

3.4.2002

..L/SD.mws

> restart:

Warning, the name changecoords has been redefined

> with(LinearAlgebra):with(linalg):with(plottools): with(plots):setoptions3d(axes=boxed,orientation=[-30,50]):

Warning, the previous binding of the name GramSchmidt has been removed and it now has an assigned value

Warning, the protected names norm and trace have been redefined and unprotected

Warning, the name arrow has been redefined

Warning, the name arrow has been redefined

> V2L:=vek->convert(vek,list):

Steepest Descent

Pieni kiusa: Joudutaan konvertoimaan Vectorin ja listan välillä. Sitä varten on yllä V2L .

Eräs mahdollisuus, jota päädyin suosittelemaan myös harj8 tehtävissä on määritellä erikseen vektorimuuttujan funktio ja muuttujajonon funktio. (Matemaattisesti nämä ovat samoja asioita, mutta Maplen kannalta edellisellä on vain yksi (vektori)argumentti ja jälkimmäisellä useita (skalaari)argumentteja.)

Tämä näkyy toimivan oikein hienosti.

> f:=(x,y)->4*x^2-4*x*y+2*y^2;

> fv:=unapply(f(X[1],X[2]),X);

> fv(X);
fv:=X->4*X[1]^ 2-4*X[1]*X[2]+2*X[2]^2; # X saa olla vektori tai lista.

Gradientti:

> g:=[D[1](f),D[2](f)];
gv:=v->g(v[1],v[2]); # Tässä vektorargumenttiimuoto.

> g(x,y); g(1,2);gv([1,2]); gv(<1,2>); # Toimii, kuten pitää.

Ryhdytään iteroimaan. Kerätään pisteet indeksoituun (taulukko)muuttujaan X.

>

>

> X[0]:=<1,-2>;

> u:=-Normalize(Vector(gv(X[0])));

> evalm(X[0]+t*u);

Tässä on suoran parametriesitys etsintäsuunnalla u. Parametri t > 0.

> phi:=t->simplify(fv(evalm(X[0]+t*u))); # Nyt on kätevää, kun meillä on vektorimuoto!

> phi(t);

> tmin:=solve(diff(phi(t)=0,t));

> X[1]:=X[0]+tmin*u;

Näin saimme ensimmäisen SD-askeleen valmiiksi. Nyt voimme tehdä pienen for-silmukan.

Huom! Ennekuin testaat, talleta ws. Kokeile aina aluksi pienellä N:n arvolla.

> X:='X':X[0]:=<1,-2>:

> N:=5:for i to N do u:=-Normalize(Vector(gv(X[i-1])));
phi:=t->simplify(fv(evalm(X[i-1]+t*u)));
tmin:=solve(diff(phi(t)=0,t));
X[i]:=(evalm(X[i-1]+tmin*u));
od:

> map(print,[seq(X[i],i=0..5)]);

> XX:=Matrix(map(V2L,[seq(X[i],i=0..5)]));

> reitti:=convert(XX,listlist);

¨Tai suoraan ilman Matrix-vaihetta:

> #map(V2L,[seq(X[i],i=0..5)]); # Matrix välivaihe ehkä auttaa näkemään paremmin.

> reitinkuva:=plot(reitti,color=blue,scaling=constrained):

> korkeusparvi:=contourplot(f(x,y),x=-1..2,y=-1..3):

> reittikorkeudet:=contourplot(f(x,y),x=-1..2,y=-1..1.6,contours=map(fv,reitti),grid=[40,40]):

> display(reitinkuva,korkeusparvi,reittikorkeudet);

>

> display(%,view=[-0.6..0.4,-1..0.3]);

>

>

> map(V2L,[seq(X[i],i=0..5)]);

>

Kokoamme proseduureiksi

Maple-ohjelmat kannattaa yleensä kirjoittaa tekstitiedostoon, joka luetaan Mapleen:

Muuta hakemistopolku tarpeen mukaan/vaihda kommentti.

> #read("c:\\opetus\\maple\\v202.mpl"):
read("/home/www/public_html/opetus/v/maple/v202.mpl");
#read("/p/edu/mat-1.414/maple/v202.mpl"); # Tämä toimii ATK-keskuksen UNIX:ssa.

Yllä olevat interaktiiviset komennot on koottu proc:ksi. Lisäpiirteenä on interpolointi 2. asteen polynomilla kussakin

1-ulotteisessa viivahaussa.

> #op(SD2); # Poista alkukommentti, niin näet koodin.

Hyödyllisiä havainnollistukseen ovat nämä, jälkimmäinen myös numeroi, mikä on tarpeen usein erit. Newtonin menetelmässä.

> #op(piirrareitti);

> #op(piirrareittin);

> with(LinearAlgebra):with(plots):

>

Warning, the assigned name GramSchmidt now has a global binding

> f:=(x,y)->exp(x)*(4*x^2+2*y^2+4*x*y+2*y+1);

> fv:=xx->f(xx[1],xx[2]);

> reitti:=SD2(f,[0,0],10);

> farvot:=map(fv,reitti);

> piirrareitti(reitti,-1..1,-1..1,scaling=constrained);

> contourplot(f(x,y),x=-1..1,y=-1..1,contours=farvot);

> display([%%,%]);

>

>