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);
>
>