Where My program Wrong ?!

when i use Newton iterative to solve a small ODEs , i met with some
problems.
the equation is :
--------------------------------------------
function dydt = test_Euler(t,y)
dydt = [1+y(1)^2*y(2)+y(2);
(1-y(1)^2)*y(2)-y(1)];
---------------------------------------------
Jacobian is :
---------------------
function dydt = Jac_test_Euler(t,y)
dydt = [ 2*y(1)*y(2), y(1)^2+1;
-2*y(1)-1, 1-y(1)^2];
----------------------
My main program is :
-----------------------------
function y=EX_BDF(h,y0,delta,epsilon,maxNum)

y = zeros(20/h,2);
y(1,:) = y0;

%%
单步BDF方法,即向后的Eu
ler方法

P = y(1,:);
%t1=0;t2=0;t3=0;
F = feval(@test_Euler,[],y0);

for i = 1:20/h-1;
J = eye(2)-h*feval(@Jac_test_Euler,[],y(i,:));
%F = feval(@test_Euler,[],y(i,:));
%P = y(i,:);
for k = 1:maxNum
Q = P + (-P + h*F'+y(i,:))'\J;
err = norm(Q-P);
relerr = err/(norm(Q)+epsilon);
P = Q;
F = feval(@test_Euler,[],P);
if k >= maxNum
k
end
if (err<delta)|(relerr<epsilon)|(abs(F)<epsilon)
break
end
end
inter = k ;
y(i+1,:) = P;
end

% plot and compare with the result of using ode45 (red)
%hold on
%plot(0:h:20-h,y(:,1));
%plot(0:h:20-h,y(:,2));
%[t,y]=ode45(@test_Euler,[0 20],y0);
%plot(t,y(:,1),'r');
%hold off
------------------------------
use the function EX_BDF like this:
EX_BDF(0.01,[0 1],1.0e-4,1.0e-6,100000);
------------------------------

the function has some problems,but i don't know it....

please help me .

thanks a lot !!!!~~~~~~~~
.