function [time,solution] = myEuler(t0,te,n,y0,fhandle)

  h = (te-t0)/n;
  solution = zeros(1,n+1);
  time     = zeros(1,n+1);
  functionindex   = 3;

  solution(1) = y0;
  time(1)     = t0;
  y = y0;

  for i=1:n

      t = t0 + i*h;
      y = y + h*feval(fhandle,t,y,functionindex);
      time(i+1) = t;
      solution(i+1) = y;

  end


