Luento
ke 30.10.02 (katsottiin jo viime to) Lähteet: KRE aivan erityisesti

Lineaariset diffyhtälösysteemit

Lähdetään opiskelemaan diffyhtälöitä hieman nurinkurisesti. Aloitetaan lineaarisista 2. kertaluvun yhtälöistä.

Syy: Haluamme jatkaa hyvään vauhtiin päässyttä lineaarialgebratoimintaa ilman pitempiä keskeytyksiä.

2. kertaluvun lineaariset diffyhtälöt

KRE CH 2 s. 64 ->

Yhteenvetoa 1. kertaluvusta

Palataan ...
Lineaariset Epälineaariset
   (EHY)  y' + p(x)y = r(x)
  1. (HY):n yleinen + (EHY):n erityinen
  2. Integroiva tekijä johtaa aina ratkaisuun
  3. E1 helppo, perustuu ratkaisukaavaan, luonteeltaan globaali
  1. Analyyttisia menetelmiä vain erikoistapauksissa:
    • Separoituva
    • Bernoulli (ohimennen)
    • Eksakti+integroiva tekijä (yleissivistäen)
  2. E1 matemaattisesti vaativa, tulos "vain" lokaali (silti hyvin merkittävä tulos).
  3. Turvauduttava useimmiten numeerisiin menetelmiin.

2. kertaluvun lineaaristen motivaatiota


2.1 Lineaarinen homogeeninen

  (EHY)  y'' + p(x)y' + q(x)y = r(x)
  (HY)   y'' + p(x)y' + q(x)y = 0
Kaikki muut ovat epälineaarisia.


Differentiaaliyhtälösysteemit

Merkintöjä:

(KRE-kirjassa käytetään y:tä x:n sijasta, no ai haitanne, vaikka vaihdellaan.)

Yleinen tapaus

    x1' = f1(t,x1, ..., xn)
    x2' = f2(t,x1., ..., xn)
     .
     .
     .

    xn' = fn(t,x1, ..., xn)  

Alkuehdot: x1(a)=b1, ... ,xn(a)=bn

Vektorimuodossa:

(DYS)     X'=F(t,X)

(AE)      X(a)=B

Lineaariset diffyhtälösysteemit

KRE 4.1-> s. 158 ->
GRE 3.9-> s. 156 ->
BdB Ch 7 s. 335 ->

Tyypillisiä ilmiöitä mallinnukseen

Yleinen lineaarinen systeemi


    x1'=a11(t)x1 + a12(t)x2 + ... + a1n(t)xn + f1(t)
    x2'=a21(t)x1 + a22(t)x2 + ... + a2n(t)xn + f2(t)
    .
    .
    .
    xn'=an1(t)x1 + an2(t)x2 + ... + ann(t)xn + fn(t)

Yllä oleva voidaan esittää matriisimuodossa:
 (LSYS)   X'=A(t)X + F(t)
Jos F(t)=0 , systeemiä sanotaan homogeeniseksi.

Käytämme lyhenteitä HY (homogeeniyhtälö) ja EHY (epähomogeeniyhtälö).

Alkuarvotehtävä (AA-teht.)

Olkoon annettu alkuehdot x1(a)=b1,x2(a)=b2 ... ,xn(a)=bn

Yllä oleva diffyhtälösysteemi yhdessä alkuehtojen kanssa muodostaa (AA)-tehtävän. Matriisimuodossa:

  (LAA)     X'=A(t)X + F(t), X(a)=B
missä B=(b1, ... ,bn)T

Basic Concepts and Theory

KRE 3.2 s. 159 ->

Peruslauseita

Lineaaristen systeemien peruslauseet

Lause 1
  1. Jos x1 ja x2 toteuttavat (HY):n, niin c1x1 + c2x2 toteuttaa (HY):n
  2. (EHY):n kaikki ratkaisut saadaan muodossa xh+xp, missä xh käy läpi kaikki (HY):n ratkaisut ja xp on jokin (EHY):n erityisratkaisu.
Tod: Täsmälleen sama kuin vastaava differenssiyhtälöille (kyse on vain ja ainoastaan lineaarisuudesta).
Lause 2, olemassaolo ja yksikäsitteisyys

Olkoot kerroinfunktiot ai,j(t) ja oikeat puolet gi(t), i,j=1 ... n jatkuvia välillä [a,b] ja olkoon t0 \in [a,b]

Tällöin jokaista alkuarvovektoria d=[d1, d2,...,dn] kohti on olemassa yksikäsitteinen (EHY):n ratkaisu x(t) siten, että

x(t0) = d.

Tod: Tätä ei todisteta. Viitataan siihen adjektiiviattribuutilla "syvällinen"

Korkeamman kertaluvun yhtälö(systeemi)n muuntaminen 1. kertal. systeemiksi

Kts. KRE s. ..., GRE s. ..., myös harj. teht   
Tämä on tärkeä seikka ja siitä seuraa:

Maple-avusteisia ratkaisutapoja

Matriisieksponenttifunktio ja yleistetyt ominaisvektorit

Skannattuja luentokalvoja

eAt Maplella

eAt:n hyödyntäminen Matlabilla

Tässä Mathworksin koodi, joka laskee suoraan sarjakehitelmällä. Opetus- ja havainnollistustarkoituksiin.
function E = expm2(A)
%EXPM2  Matrix exponential via Taylor series.
%   E = expm2(A) illustrates the classic definition for the
%   matrix exponential.  As a practical numerical method,
%   this is often slow and inaccurate.
%
%   See also EXPM, EXPM1, EXPM3.

%   Copyright (c) 1984-98 by The MathWorks, Inc.
%   $Revision: 5.4 $  $Date: 1997/11/21 23:38:26 $

% Taylor series for exp(A)
E = zeros(size(A));
F = eye(size(A));
k = 1;
while norm(E+F-E,1) > 0
   E = E + F;
   F = A*F/k;
   k = k+1;
end

expm-funktio, tositarkoitukseen

>> help expm 

 EXPM   Matrix exponential.
    EXPM(X) is the matrix exponential of X.  EXPM is computed using
    a scaling and squaring algorithm with a Pade approximation.
 
    Although it is not computed this way, if X has a full set
    of eigenvectors V with corresponding eigenvalues D, then
    [V,D] = EIG(X) and EXPM(X) = V*diag(exp(diag(D)))/V.
 
    EXPM1, EXPM2 and EXPM3 are alternative methods.
 
    EXP(X) (that's without the M) does it element-by-element.
 
    See also EXPM1, EXPM2, EXPM3, LOGM, SQRTM, FUNM.


expm-käyttöinen lineaarisen systeemin ratkaisija

Seuraavaksi oma koodi, jonka "tapahtumarivi" on tämä:
y(:,i) = expm(A * t(i)) * x0; 

function [t,y] = linsys(A,x0,T,siz)
% Input: Matriisi A, alkuarvovektori x0, loppuaika T (aikaväli: 0 .. T)
% Valinnainen siz: kuinka moneen osaan aikaväli jaetaan, oletus 100.
%
% Output: t - vektori:  diskretoitu aika-akseli (oletus 100-pituinen)
%         y - matriisi: 100 x n, kukin sarake edustaa ratkaisufunktion
%                       arvoja t-aikapisteissä.
% 
% Esim: A=[1 0 0;1 3 0;0 1 1];x0=[1;2;3];
%       [t,y]=linsys(A,x0,2);
%       plot(t,y)                       % aikariippuvuusparvi (mieliv. n)
%       plot(y(:,1),y(:,2))             % faasitaso (1,2)  
%       plot3(y(:,1),y(:,2),y(:,3))     % faasiavaruus (3-d)
%       plot(y(:,i),y(:,j))             % projektio (i,j)-faasitasossa
%                                       % (jos n > 2)
if nargin < 4, siz = 100; end; 
t = linspace(0,T,siz);
m=size(A);m=m(1);
y = zeros(m,siz);
for i=1:siz,
y(:,i) = expm(A * t(i)) * x0;
end;
y=y';      % Transponoidaan y-matriisi,jotta yhdenmukainen ode-funkt. kanssa

2x2-systeemien faasitasoluokittelu

Tässä eräs kuva (Robert Pichet'n kurssimateriaalia)