function [t,Y]=lindys(A,tvali,y0)
% HA 26.11.08
% Ratkaistaan lineaarinen DYS matriisiexp-funktiolla
% Perustuu Matlab:n expm-funktioon.
% Syötteet:
%    matriisi A (n x n)
%    aikaväli tvali [a,b]
%    alkuarvovektori  y0 (hetkellä t=0)
% Tulos: t:aikapisteet 100 kpl.
%        Y (100,n) : ratkaisusarakkeet
% Esim: A=[1 0 0;1 3 0;0 1 1]; y0=[1;2;3];
%       [T,Y]=lindys(A,[-2 5],y0);
%       plot(T,Y)  % aikaparvi
%       plot(Y(:,1),Y(:,2)  % Faasitaso (1 2)
%       plot3(Y(:,1),Y(:,2),Y(:,3)); grid on % Faasiavaruus

t=linspace(tvali(1),tvali(2));
[n,n]=size(A);
Y=zeros(n,100);
for j=1:100
    Y(:,j)=expm(A*t(j))*y0;
end;    
Y=Y';