function [T,Y]=eulerV(fnimi,Tspan,ya,n)
% fnimi : funktiomääritys, inline tai m-tiedosto, määritellään kuten
%         ode45:lle ym.  (help ode45)
%         Tspan on muotoa [a b]  ya:alkuarvovektori (y(a)=ya)
%          n: välin (a,b) jakovälien lukumäärä
% Tulos: T: sarake, aikapisteet.
%        Y: Matriisi, jonka j:s sarake ilmaisee ratkaisufunktion
%           y_j arvot T-aikapisteissä.
%           i:s rivi: T(i)  Y(i,1), Y(i,2), ..., Y(i,N)
% Siis jos on 2x2-systeemi, niin Y:ssä on 2 saraketta. 
%
% Esim: Fsys=inline('[y(1)-y(1)*y(2)/2;-2*y(2)+y(1)*y(2)]','t','y'); 
%       y0=[20;10];
%       [T,Y]=eulerV(Fsys,[0,1],y0,20); plot(T,Y)
     %  Huom! Jos Fsys määritellään m-funktioksi, on kutsussa oltava 
     %  'Fsys' tai @Fsys.
%
a=Tspan(1);b=Tspan(2);
h=(b-a)/n;
N=length(ya);
Y=zeros(n+1,N);  % i:s rivi:  Y(i,1), Y(i,2), ..., Y(i,N)
T=a:h:b;
Y(1,:)=ya';

for i=1:n
   Y(i+1,:)=Y(i,:)+h*(feval(fnimi,T(i),Y(i,:)))'; 
   % Ffunktio 'fnimi' palauttaa pystyvektorin, siksi pitää transponoida.
end;
