function priklad01
%% Vstupni data
a = 1; b = 1; % [m] rozmery konstrukce
k = 1;        % [Nm-1] tuhost pruziny
F = .25;      % [N] pusobici sila
L = sqrt(a^2+b^2); % puvodni delka pruziny


ezplot(@Pi, [-4,1] );
xlabel('posun' );
ylabel('energie' );

%% Newtonova metoda
v = -2; % Pocatecni odhad
g = gfun(v); % Pocatecni hodnota gradientu
i = 0; % Citac iteraci 
clc; % Vymazej obrazovku

while( abs(g) > 1e-12 )
    disp( ['iterace = ' num2str(i) ' gradient = ' num2str(g) ' posun=' num2str(v)] );
    i = i + 1; % Zvysim citac iteraci 
    h = hfun(v); % Spocti druhe derivace
    v = v - g/h; % Aktualizuj puvodni reseni
    g = gfun(v); % Aktualizuj gradient
end

disp( ['iterace = ' num2str(i) ' gradient = ' num2str(g) ' posun=' num2str(v)] );

%% Cilova funkce a analyza citlivosti
    function val=Pi(v) % Potencialni energie systemu pro posun v
        val = .5*k*D(v).^2 + F*v;
    end

    function val=gfun(v) % Hodnota gradientu objektivni funkce
        val = k*D(v)*(b+v)/l(v) + F;
    end

    function val=hfun(v) % Hodnota druhe derivace funkce
       val=k*( ((b+v)/(l(v)))^2*(1-D(v)/l(v)) + D(v)/l(v));
    end

    function val=D(v) % Pro dany v urci prodlouzeni/zkraceni pruziny
        val = l(v) - L;
    end

    function val=l(v) % Pro dane v urciti novou delku pruziny
       val = sqrt( a^2 + (b+v)^2 ); 
    end
end
