Harj. 6 LV

pe 8.3.02 HA

Alustukset

> restart:

Warning, the name changecoords has been redefined

> with(LinearAlgebra): with(linalg):

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

> alias(Inv=MatrixInverse,IdM=IdentityMatrix,Diag=DiagonalMatrix,Tr=Transpose):

> with(plots): with(plottools):

Warning, the name arrow has been redefined

> V2L:=vek->convert(vek,list):
v2l:=vek->convert(vek,list):
L2V:=lis->convert(lis,Vector):
l2v:=lis->convert(lis,vector):

> jana3dpv:=(p,v)->spacecurve([p,V2L(Vector(p)+Vector(v))],thickness=2,color=blue):

1.

> f:=(x,y)->arctan(y/x);

> x0:=2:y0:=2;

> T:=f(x0,y0)+D[1](f)(x0,y0)*(x-x0)+D[2](f)(x0,y0)*(y-y0);

> subs(x=x0,y=y0,T);
f(2,2);

> hpmax:=1.5:htmax:=0.8:

> pinta:=plot3d(f(x,y),x=x0-hpmax..x0+hpmax,y=y0-hpmax..y0+hpmax,style=patchcontour):

> taso:=plot3d(T,x=x0-htmax..x0+htmax,y=y0-htmax..y0+htmax,style=PATCH,color=red,grid=[7,7]):

> display(pinta,taso);

> N:=[D[1](f)(x0,y0),D[2](f)(x0,y0),-1];

> display(pinta,taso,jana3dpv([x0,y0,f(x0,y0)],N),jana3dpv([x0,y0,f(x0,y0)],-N));

Pinta, tangenttitaso ja normaali (pyörittele vähän, niin näkyy paremmin).

2 ja 3.

> A:=<<1/2,1/4,1/4>|<1/4,1/2,1/4>|<1/4,1/4,1/2>>;

> A-Tr(A);

Olishan tuon nähnyt päältäpäinkin.

Koska matriisi on symmetrinen , niin peruslause => ominaisvekroreista voidaan muodostaa ortonormaali kanta.

> p:=det(A-lambda*IdM(3));factor(p);# Karakteristinen polynomi

> lam:=solve(%=0,lambda);

> M1:=A-lam[1]*IdM(3);
M2:=A-lam[2]*IdM(3);

> gaussjord(M1),gaussjord(M2);

x3 voidaan valita vapaasti, esim. x3=1. x2=x1 (=1), x3=x1 (=1). Siis ominaisvektori [1,1,1] ja normeerattuna:

> X1:=1/sqrt(3)*<1,1,1>;

. x3 ja x2 voidaan valita vapaasti: x1 = -x2 - x3.

Saadaan 2 LRV ominaisvektoria valitsemalla a) x2=1, x3 =0 => x1=-1 ja b) x2=0, x3=1 => x1=-1, eli:

> X2:=1/sqrt(2)*<-1,1,0>;X3:=1/sqrt(2)*<-1,0,1>;

Tehtävä 2 on valmis (ja vähän ylikin, ei siinä pyydetty edes normeerausta).

> Tr(X1).X2;Tr(X1).X3;

Teorian mukaan eri ominaisarvoihin liittyvät ominaisvektorit ovat ortogonaaliset, joten näin käy aina (jos oikein laskettu).

Sensijaan kaksinkertaiseen ominaisarvoon liittyvät ominaisvektorit eivät ole automaattisesti ortogonaaliset.

Ortogonaalinen kanta saadaan ajattelemalla geometrisen havainnollisesti: Projisioidaan toinen vektori toiselle ja vähennetää tämä projektiovektori. Normeerataan ensin (sehän jo tehtiin), jolloin kaava tulee yksinkertaiseen muotoon:

X2 - (X2,X3)X3

> V2:=X2-(Tr(X2).X3).X3;

> Tr(V2).X3;

> X2:=V2/Norm(V2,2);

> X:=<X1|X2|X3>;

> Tr(X).X;

Käsin laskettuna saattoi olla hieman vaikea, kun ortogonalsiointitekniikkaa ei ole kunnolla esitelty. Sinänsä, kun tämän on ajatellut havainnollisesti, asia on hyvin luonteva.

Sovitaanpa, että j os tällainen tehtävä tulee kokeeseen, niin tuo projektiokaava ohjeineen annetaan.

Maplella automaattisemmin

> omsys:=Eigenvectors(A);

> V:=omsys[2];

Kaksi ensimmäistä vektoria eivät taaskaan ole ortogonaalisia. Help:stä löytyy GramScmidt (edellyttää tämän yhdistelmän assosiontia ortogonalisointiin).

Kaiken lisäksi GramScmidt näkyy toimivan versiossa 7, mutta ei 6:ssa. (Helppiesimerkki (cut/paste) menee virheseen.)

Antaapa sen nyt siis olla ja todetaan, että loppu tehdään, kuten yllä.

> #?LinearAlgebra

>

> GramSchmidt([<-1,1,0>,<-1,0,1>]); # Eipä toimi, linalg:n vastaava varmasti toimii, mutta olkoonpa nyt.

Error, (in GramSchmidt) 1st argument must be a non-empty list or set of vectors

Error, (in GramSchmidt) 1st argument must be a non-empty list or set of vectors

4.

a)

> q:=2*x1^2+4*x2^2+x1*x2;

> A:=<<2,1/2>|<1/2,4>>;

> x:=<x1,x2>;

> Tr(x).A.x;expand(%);

Oikein on. Katsotaan vielä Hessian ja pannaan vierelle 2*A

> hessian(q,[x1,x2]),2*A; # hessian asustaa vanhassa linalg:ssa

> Eigenvalues(A);

> omsys:=Eigenvectors(A);

> eigenvectors(A);

> DD:=Diag(omsys[1]);

> X1:=omsys[2][1..2,1]: X1:=X1/Norm(X1,2):X1:=map(simplify,%);

> Norm(X1,2);simplify(%);

> X2:=omsys[2][1..2,2]: X2:=X2/Norm(X2,2):X2:=map(simplify,%);

> X:=<X1|X2>;

> XT:=Tr(X);

> A-map(simplify,X.DD.XT);map(simplify,%);

> Q:=<y1|y2>.DD.<y1,y2>;

> implicitplot(q=1,x1=-2..2,x2=-2..2,grid=[50,50]);

> implicitplot(Q=1,y1=-2..2,y2=-2..2,grid=[50,50]);

>

>

b)

> A:=<<0,1,0>|<1,0,1>|<0,1,1>>;

> <x1|x2|x3>.A.<x1,x2,x3>;q:=expand(%);

> p:=det(A-lambda*IdM(3));

> lam:=fsolve(p=0,lambda);

> Q:=lam[1]*y1^2+lam[2]*y2^2+lam[3]*y3^2;

> implicitplot3d(Q=1,y1=-2..2,y2=.2..2,y3=-2..2);

>

>

> A:=<<0,1,0>|<1,0,1>|<0,1,1.>>;

> lam:=Eigenvalues(A);

> map(Im,lam);lam:=map(Re,lam);

> Q:=<y1|y2|y3>.Diag(lam).<y1,y2,y3>;

c)

> A:=<<1,0,0,0>|<0,0,-2,0>|<0,-2,0,0>|<0,0,0,2>>;

> <x1|x2|x3|x4>.A.<x1,x2,x3,x4>;

> lam:=Eigenvalues(A);

> Q:=add(lam[k]*y[k]^2,k=1..4);

> sum(lam[k]*y[k]^2,k=1..4); # mystistä!

Error, bad index into Vector

> lam[4]; # Vektorin indekseissä ei ole oikeasti vikaa, toimiihan myös add.

> ll:=convert(lam,list);

> sum(ll[i]*y[i]^2,i=1..4);

Listaksi konversion jälkeen toimii. "Bad index into vector" on tosi outo juttu, mutta olkoonpa nyt.

5.

> A:=<<1,12>|<12,-6>>;

> q:=expand(<x1|x2>.A.<x1,x2>);

> (lam,X):=Eigenvectors(A);

> normit:=Norm(X[1..2,1],2),Norm(X[1..2,2],2);

> X:=X/normit[1]; # koska kerran samoja

> Inv(X), Tr(X); # Ovat samoja (hyvä aina tarkistaa).

> y1:=(1/sqrt(3))*sinh(t);y2p:=(1/sqrt(2))*cosh(t);y2m:=-y2p;

> pos:=plot([y1,y2p,t=-2..2],view=[-4..4,-4..4],scaling=constrained):

> neg:=plot([y1,y2m,t=-2..2],view=[-4..4,-4..4],scaling=constrained):

> display(pos,neg);

> x1x2p:=X.<y1,y2p>;

> x1x2m:=X.<y1,y2m>;

> phaara:=plot([op(V2L(x1x2p)),t=-3..3],x1=-8..8,x2=-8..8):

> nhaara:=plot([op(V2L(x1x2m)),t=-3..3]):

> display(phaara,nhaara);

> implicitplot(q=5,x1=-8..8,x2=-8..8,grid=[50,50]);

>

>

> V2L(x1x2p);

>

>