Harj 8 LV
pe 8.11.02
Alustukset
> | restart:with(LinearAlgebra):with(linalg): |
Warning, the name changecoords has been redefined
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(Id=IdentityMatrix): |
> | alias(rref=ReducedRowEchelonForm): alias(Diag=DiagonalMatrix):alias(ref=GaussianElimination): |
3.
> | A:=<<alpha,-2>|<2,0>>; |
> | p:=det(A-lambda*Id(2)); |
> | solve(p=0,lambda); |
1) => kompleksijuuret, lähdefokus (epästabiili), jos , nielufokus (stabiili), jos (Jos , saadaan keskus, joka on heikosti stabiili)
2)
Molemmat ominaisarvot samanmerkkiset (tulo =4). Jos , niin positiiviset, jos , negatiiviset.
Siten alpha > 4 ==> lähde (epästabiili), alpha < 0 ==> nielu (stabiili).
Mielenkiintoisimmat ovat ehkä dgeneraatiotapaukset 4 ja -4.
Tällöin käytetään Cayley-Hamiltonia:
> | alpha:=4: lambda:=alpha/2; |
> | At:='lambda*t*Id(2)'+'t*(A-lambda*Id(2))'; |
> | expAt:=exp(lambda*t)*(Id(2)+t*(A-lambda*Id(2))); |
> | evalm(%); |
> | exponential(A,t); |
4.
a)
diag(d1,d2, .. dn)^k = diag(d1^k,d2^k,...,dn^k). (Diagonaalimatriisien kertolaskussa diagonaalialkiot kerrotaan keskenään.)
> | At:=<<t*d1,0>|<0,t*d2>>; |
> | add((At^k)/k!,k=0..10); |
> | evalm(%); |
Tässä näkyy selkeästi, miten matriisisarja muodostuu. Diagonaalialkiosarjat ovat :n ja :n sarjoja. Eli:
> | exponential(At); |
Tämä tietysti pätee yleisellä n:llä.
Tässä siis:
> | A:=<<1,0>|<0,2>>; exponential(A,t); |
b.
> | A:=<<0,0>|<1,0>>; At:=<<0,0>|<t*1,0>>; # LinearAlgebra on huono sieventelemään matriiislausekkeita. Siksi tämä At. |
> |
> | At^2; |
> | expAt:=Id(2)+At; # Tästä eteenpäin 0:aa |
> | exponential(A,t); # Tarkistus |
c.)
> | A:=<<0,-1>|<1,0>>; |
> | seq((A^k),k=0..10); |
4:n pituinen jakso. Siitä nähdään kerroinjonot.
> | At:=<<0,-t>|<t,0>>; |
> | add((At^k)/k!,k=0..10);evalm(%); |
Kun vielä kertomat otetaan mukaan, saadaan tällaisia osasummia alkiosarjoille, ts.
> | expAt:=<<cos(t),-sin(t)>|<sin(t),cos(t)>>; |
Toki voidaan myös diagonalisoida ja käyttää Eulerin kaavaa.
5.
> | A:=<<1,-1>|<4,-3>>; |
> | g:=t-><t,1>; |
> | (lambda,ov):=Eigenvectors(A); |
Lasketaan taiteen sääntöjen mukaan, vaikka tehtävässä ei vaadittu. (Olis kannattanut):
> | 'A'*t=(-t)*'Id'+('A'+'Id')*t; |
> | eAt:=exp(-t)*(Id(2)+(A+Id(2))*t); |
> | (A+Id(2))^2; # Toden totta, Cayley-Hamilton toimii! |
> | eAt;eAt:=evalm(%); |
> | exponential(A,t); |
Saatiin sama!
> | eAt:=Matrix(eAt); |
> | yh:=eAt.<C1,C2>; |
Sitten se kaunis EHY-kaava:
> | expA:=unapply(eAt,t): |
> | y(t)=e^('A'*t)*C+Int(e^'A(t-s)'*'g(s)',s=0..t); #tekstinkäsittelyä |
> | yp:=map(int,expA(t-s).g(s),s=0..t); |
> | y:=yh+yp; |
Tarkistetaan, että toteuttaa diffyhtälön:
> | vas:=map(diff,y,t); oik:=A.y+g(t); |
> | map(simplify,vas); |
> | map(simplify,oik); |
> | %-%%; |
Hyvin kävi!
Sitten alkuehto:
> | subs(t=0,y); |
Voiko olla totta, no voi näinkin käydä. Siispä:
> | yAA:=subs(C1=1,C2=-1,y); |
> | subs(t=0,yAA); |
Niin se vaan on!
6.
> | A:=<<1,1,0,0>|<0,1,0,0>|<1,2,2,1>|<2,1,0,1>>; |
> | (lambda,ov):=Eigenvectors(A); |
> | v1:=ov[1..-1,1]; |
Etsitään ominaiarvoa 1 vastaavat yleistetyt ominaisvektorit:
> | u1:=ov[1..-1,2]; # Ensin se oikea ominaisvektori |
> | M2:=(A-Id(4))^2; |
> | rref(M2); |
> | u2:=<1,0,0,0>; |
> | M3:=(A-Id(4))^3; |
> | rref(M3); |
> | u3:=<0,0,0,1>; |
> | v1,u1,u2,u3; |
Tästä saadaan tulos, mutta nyt en enää ehdi kirjoitella.