clear all


f=inline('x*(3-x-2*y)','x', 'y');
g=inline('y*(2-x-y)','x','y');





    
for k   = 1:9 

    k
    x(1)   = 0.1+k*0.2;
    y(1)   = 0.5 
    t      = 0;
    Tfin   = 10;
    DeltaT = 0.1;%10^(-k);
    N      = round(Tfin/DeltaT);

    
        
                   
                        % Runge-Kutta method: errore globale O(h^4)

                   for i=1:N
                          
                          k1x = DeltaT*f(x(i),y(i));
                          k1y = DeltaT*g(x(i),y(i));
                          
                          k2x = DeltaT*f(x(i)+k1x/2, y(i)+k1y/2);
                          k2y = DeltaT*g(x(i)+k1x/2, y(i)+k1y/2);
                          
                          k3x = DeltaT*f(x(i)+k2x/2, y(i)+k2y/2);
                          k3y = DeltaT*g(x(i)+k2x/2, y(i)+k2y/2);                     
                     
                          k4x = DeltaT*f(x(i)+k3x,y(i)+k3y);
                          k4y = DeltaT*g(x(i)+k3x,y(i)+k3y);
                          
                  
                          
                          x(i+1)=x(i)+(k1x+2*k2x+2*k3x+k4x)/6;
                          y(i+1)=y(i)+(k1y+2*k2y+2*k3y+k4y)/6;
                          
                          t=[t i*DeltaT];
    
                   
                   
           
             
                   end    
       
                   
            figure(1)
            plot(x,y,'g')
            hold on 
            plot(1,1,'ro')
            plot(3,0,'o')
            plot(0,2,'o')
            plot(0,0,'ro')
pause
%            figure(2)
%            plot(t,x)
            
%h = [h DeltaT];




end