function result = mvnewton(f,g,aprox,tol)

%Argumenty fce:
% f,g - matlab function_handle
% aprox - pocatecni iterace (bod o dvou souradnicich)
% tol - volitelny argument, tolerance pri jejiz dosazeni iterovani konci
%Return fce:
% result - priblizna hodnota hledaneho reseni (bod o dvou souradnicich)

if nargin < 3
    error('Nedostatecny pocet zadanych argumentu')
end

if(~isa(f,'function_handle') || ~isa(g,'function_handle'))
    error('f nebo g nejsou matlabovsky function_handle')
end

if (~exist('tol', 'var'))
        tol = 0.01;
end

iter = 1;
maxiter = 1000; % maxiter - maximalni pocet iteraci pri kterem iterace konci
syms x;
syms y;
f = sym(f);
g = sym(g);
fx = diff(f,x);
fy = diff(f,y);
gx = diff(g,x);
gy = diff(g,y);

F = [f;g];
J = [fx, fy ; gx, gy];

x = aprox(1); y = aprox(2);

xy = [x;y] - inv(eval(J))*eval(F);
x = xy(1); y = xy(2);
if isnan(xy)
    error('x nebo y nabylo v prubehu iterace hodnoty NaN')
end
result = xy;


while (abs(eval(f))+ abs(eval(g))) > tol
    if(iter > maxiter)
        error('1 000 iteraci nedostacovalo k dosahnuti pozadovane presnosti')
    end

    if(~isreal(aprox) || isnan(x) || isnan(y) || isinf(x) || isinf(y))
        error('V prubehu iterace souradnice x,y nebyly konecne realne c.')
    end

    result = xy;

    xy = [x;y] - inv(eval(J))*eval(F);
    x = xy(1); y = xy(2);
    iter = iter+1;
end
end