% Zgled za resevanje nelienarnega sistema z Newtonovo metodo
% Iscemo tisto vrednost parametra r in seveda tudi tocko (x,y),
% kjer se krivulji x^2+y^2=r^2 in x^2-2xy+y=1 dotikata

% Bor Plestenjak 2004

% najprej narisemo, kako zgledajo presecisca za razlicne r

[X,Y]=meshgrid(-1.8:0.01:1.8,-1.8:0.01:1.8);

Z1=sqrt(X.^2+Y.^2);
Z2=X.^2-2*X.*Y+Y-1;

[C,h]=contour(X,Y,Z1,[1.2 1 0.8 0.6]);
clabel(C,h)
hold on
contour(X,Y,Z2,[0 0],'b');
hold off

disp('Pritisni enter')
pause

% iscemo radij, pri katerem se x^2+y^2=r^2 dotakne krivulje x^2-2xy+y-1
% takrat dobimo se tretjo enacbo - determinanta Jacobijeve matrike, ki mora biti singularna
% resimo nelinearni sistem 3x3 in dobimo resitev

F=inline(strcat(...
   '[ x(1)^2 + x(2)^2 - x(3)^2,',...
   '  x(1)^2 - 2*x(1)*x(2) + x(2) - 1,',...
   '  2*x(1) - 4*x(1)^2 - 4*x(1)*x(2) + 4*x(2)^2]'));

JF=inline(strcat(...
   '[ 2*x(1), 2*x(2), -2*x(3);',...    
   '  2*x(1)-2*x(2), -2*x(1)+1, 0;',... 
   '  2-8*x(1)-4*x(2), -4*x(1)+8*x(2), 0]'));

newton(F,JF,[-0.5;0.5;0.7])
