Harj. 4 Av ratkaisuja

V3 10.10.2000

Teht. 1

> restart:with(linalg):

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

> A:=matrix([[cos(Theta),-sin(Theta)],[sin(Theta),cos(Theta)]]);

> eigenvectors(A);

Error, (in eigenvectors) eigenvects only works for a matrix of rationals, rational functions, algebraic numbers, or algebraic functions at present

> id:=diag(1,1);

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

> solve(p=0,lambda);map(simplify,[%],{cos(Theta)^2=1-sin(Theta)^2});

> lambda:=cos(Theta)+I*sin(Theta);

> Cmat:=evalm(A-lambda*id);gausselim(%);

> v1:=vector([I,1]);v2:=map(conjugate,v1);

Teht. 2

> restart: with(linalg):

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

a)

> A:=matrix([[2,-3],[0,0]]);

> Cmat:=A-lambda*diag(1,1); p:=det(Cmat);

Kaksi erillistä reaalista ominaisarvoa, joten 2 LRT om. vektoria. Siten matriisi on diagva.

> Cmatit:=[subs(lambda=0,evalm(Cmat)),subs(lambda=2,evalm(Cmat))];

> map(gaussjord,Cmatit);

> v0:=vector([3/2,1]);v2:=vector([1,0]);

> T:=augment(v0,v2);DD:=diag(0,2);TI:=inverse(T);
# Huom! D ei sovellu nimeksi, se on (derivaatalle) varattu nimi.

Tarkistus: Pannaan rinnakkain: Samathan siitä tulee.

> evalm(T&*DD&*TI),op(A);

b)

> A:=matrix([[2,4],[-1,-2]]);eigenvectors(A);

Kaksinkert. ominaisarvo 0, vastaava ominaisavaruus vain 1-ulotteinen. Ts. matriisilla

ei ole kahta LRT ominaisvektoria, joten se ei diagonalisoidu.

Lasketaan vielä "käsin".

> restart: with(linalg):A:=matrix([[2,4],[-1,-2]]);

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

> Cmat:=evalm(A-lambda*(diag(1,1))); p:=det(Cmat);

Siis on kaksinkertainen ominaisarvo. Tutkitaan ominaisavaruuden domensiota, eli määrätään ominaisvektorit:

> Cmat:=subs(lambda=0,op(Cmat));gaussjord(Cmat);

Saadaan vain yksi LRT ominaisvektori [2,-1].

Teht. 3

> with(LinearAlgebra):

Warning, the assigned name GramSchmidt now has a global binding

> A:=Matrix([[2,10],[-1,4]]);

> omsys:=Eigenvectors(A);

> omvekt:=omsys[2];

> oma:=omsys[1];

> w:=Column(omvekt,1);

> lambda:=oma[1];

> alpha:=Re(lambda);beta:=Im(lambda);

> u:=map(Re,w);v:=map(Im,w);

> uv:=<u | v>;

> X:=uv.Transpose(<exp(alpha*t)*(a*cos(beta*t)+b*sin(beta*t))|exp(alpha*t)*(b*cos(beta*t)-a*sin(beta*t))>);

>

>

> AEa:={subs(t=0,X[1])=1,subs(t=0,X[2])=0};
AEb:={subs(t=0,X[1])=0,subs(t=0,X[2])=1};

> kerta:=solve(AEa,{a,b});kertb:=solve(AEb,{a,b});

> Xa:=subs(kerta,X);Xb:=subs(kertb,X);

> plot([Xa[1],Xa[2],t=-3..0.5],x=-5..2,y=-1..1);plot([Xa[1],Xa[2]],t=-3..0.5);

> plot([Xb[1],Xb[2],t=-3..2],x=-1.2..0.3,y=-0.2..1.5);plot([Xb[1],Xb[2]],t=-3..2,xy=-2..2);

>

Ominaisarvot/vektorit käsinlaskutyylillä

> restart:with(linalg):A:=matrix([[2,10],[-1,4]]);

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

Kompleksisten ominaisarvojen tapauksessa rarkaisu on

> X:=exp(alpha*t)*(a*cos(beta*t)+b*sin(beta*t))*u+exp(alpha*t)*(b*cos(beta*t)-a*sin(beta*t))*v;

> p:=det(A-lambda*diag(1,1));

> lambda:='lambda':lam:=solve(p=0,lambda);lambda:=lam[1];alpha:=Re(lambda);beta:=Im(lambda);

> Cmat:=evalm(A-lambda*diag(1,1));gaussjord(Cmat);

> w:=vector([1-3*I,1]);u:=map(Re,w);v:=map(Im,w);

> uv:=augment(u,v);

a)

> ab:=linsolve(uv,vector([1,0]));

> Xa:=subs(a=aba[1],b=aba[2],X);
evalm(Xa);

b)

> ab:=linsolve(uv,vector([0,1]));

> Xb:=subs(a=ab[1],b=ab[2],X);
evalm(Xb);map(simplify,%);

dsolve-"nappia painamalla"

> restart:dys:=diff(x(t),t)=2*x(t)+10*y(t),diff(y(t),t)=-x(t)+4*y(t);

> xa:=dsolve({dys,x(0)=1,y(0)=0},{x(t),y(t)});

> xb:=dsolve({dys,x(0)=0,y(0)=1},{x(t),y(t)});

Teht. 4

> restart:with(LinearAlgebra):

> A:=Matrix([[-2/25,1/50],[2/25,-2/25]]);

> esys:=Eigenvectors(A);

> lambda:=esys[1];V:=esys[2];

> X:=C1*exp(lambda[1]*t)*Column(V,1)+C2*exp(lambda[2]*t)*Column(V,2);

> xy:=evalm(%);

> AE:=subs(t=0,xy[1])=25,subs(t=0,xy[2])=0;

> Ceet:=solve({AE},{C1,C2});

> xy:=subs(Ceet,op(xy));

> x:=xy[1];y:=xy[2];

> plot([x,y],t=0..100);

> solve(x=y,t);evalf(%);

> plot([x,y,t=0..100]);

>

Katsotaan, miten pannaan pakettiin dsolvella

> restart:dys:=diff(x(t),t)=-2/25*x(t)+1/50*y(t),diff(y(t),t)=2/25*x(t)-2/25*y(t);

> dsolve({dys,x(0)=25,y(0)=0},{x(t),y(t)});

>