ODE with Newton method

I am trying to solve the ODE y'=x*cos(x^2)*y^2 using newton's method but my code keeps running forever.

Please, I really appreciate some help with this.

clear all
close all

n=1001
y0=1;
xmin=0; xmax=10; 
x=linspace(xmin,xmax,n); %grid
step=(xmax-xmin)/n; %Step Size
yy=zeros(1,n); %Array for results for EEE

yy(1) = y0; % Initial value


%Define the Function and derivative
f=@(x,y) x*cos(x^2)*y^2;
df=@(x,y) 2*x*cos(x^2)*y ;

F=@(x,y,yn) y-yn-step*f(x,y); %Function to set to zero
J=@(x,y) 1 - step*df(x,y); %Jacobian

tol=1.e-5; % Tolerance

for i=xmin+1:n-1
    yy(i+1)=yy(i); %initial guess for Newton's method

    res=-J(yy(i+1))\F(yy(i+1),yy(i)); 
    while (norm(res,inf)>1.e-10)
        yy(i+1)=yy(i+1) + res;
        res=-J(yy(1,1))\F(yy(i+1),yy(i));
    end
    yy(i+1)= yy(i+1) + res;
 end

plot(x,yy,'k--')
xlabel('t')
ylabel('y')

Thanks

[Torsten]

For comparison:

fun = @(x,y)x*cos(x^2)*y^2;
n = 1001;
xmin=0; xmax=10; 
x=linspace(xmin,xmax,n); %grid
y0 = 1;
[X,Y] = ode45(fun,x,y0);
plot(X,Y)

Answer 1

[VBBV]

clear all
close all

n=1001

n = 1001

y0=1;
xmin=0; xmax=10; 
x=linspace(xmin,xmax,n); %grid
step=(xmax-xmin)/n; %Step Size
yy=zeros(1,n); %Array for results for EEE

yy(1) = y0; % Initial value


%Define the Function and derivative
f=@(x,y) x*cos(x^2)*y^2;
df=@(x,y) 2*x*cos(x^2)*y ;

F=@(x,y,yn) y-yn-step*f(x,y) %Function to set to zero

F = function_handle with value:
    @(x,y,yn)y-yn-step*f(x,y)

J=@(x,y) 1 - step*df(x,y) %Jacobian

J = function_handle with value:
    @(x,y)1-step*df(x,y)

tol=1.e-5; % Tolerance

for i=xmin+1:n-1
    yy(i+1)=yy(i); %initial guess for Newton's method

    res=-J(yy(i+1),x(i))\F(x(i),yy(i+1),yy(i));  % why missing input arguments ???
    while (norm(res,inf)>1.e-10)
        yy(i+1)=yy(i+1) + res;
        res=-J(yy(1,1),x(i))\F(x(i),yy(i+1),yy(i));
    end
    yy(i+1)= yy(i+1) + res;
 end

plot(x,yy,'k--')
xlabel('t')
ylabel('y')