f = inline('y*(y - 1)','y','t')
   
    t0 = 0;
    y0 = 0.5;
	h = 0.5;
	n = 5;
    t1 = 0;
    
    if (n <=0)
        n = abs(t1 - t0)/h;
    end
    
    Yr =[];
    Tr =[];
    
     y = y0;
	 t = t0;
     
     for j = 1:n
	 k1 = feval(f,y,t);
	 k2 = feval(f,y+ h*k1, t+h);
     y = y + (h/2)*(k1+k2);
	 t = t+h;
	 Yr(j) = y;
	 Tr(j) = t;
	end
        
	R = [Tr ; Yr];
    R'