function [t,y]=expeuler_rendszer(diffegy,t0,y0,h,N) 
% t0: kezdeti idõpillanat
% y0: kezdeti vektor
% N: idolepesek szama a numerikus megoldas soran
t=zeros(N+1,1); 
y=zeros(2,N+1);
t(1)=t0; 
y(:,1)=y0; 
for i=1:N
    t(i+1) = t(i) + h;
    y(:,i+1) = y(:,i) + h * diffegy(t(i),y(:,i)); 
end
plot(t,y)