%M2 lokalni extremy
syms x y                               %deklarace symbolickych promennych
% ZADANE FUNKCE
%1
 f=(x^2+4*x)*y+y^2;
% %2
% f= x^4-2*x^2+y^4+2*y^2;
% %3
% f= -x^4-2*x^2+y^4-2*y^2;
% %4
% f= x^4+2*x^2+y^4-2*y^2;
% %5
% f=1/3*x^3 +1/2*x^2 - 6*x +1/3* y^3 +1/2* y^2 -20 *y;
% %6
% f= x^5-5*x+y^3-3*y;
% %7
% f= 2*x^3+2*y^3+x^2-4*x-24*y;
% %8
% f= 1/3*x^3-x*y^2+2*y+4 ;
% %9
% f= exp(-x^2-y^2)*(2*y^2+x^2);

% Derivace 1. radu
fx=diff(f,x)
fy=diff(f,y)

% Stacionarni body
[X,Y]=solve(fx, fy)                  % reseni fx=0, fy=0
                                     % X,Y vektory - stac. body, symbolické
                                     % X(1),X(2)... x-ove souradnice
                                     % Y(1),Y(2)... y-ove souradnice
% vyhodime komplexni cisla
indkomplex=find(imag(X)~=0);
X(indkomplex)=[]; Y(indkomplex)=[];
indkomplex=find(imag(Y)~=0);
X(indkomplex)=[]; Y(indkomplex)=[];
%%%%%%%%%%%%%%%%%%%%%%%%%%

% Derivace 2.radu
H=hessian(f,[x y])               % H je matice z derivaci 2. radu (Hessova)
detH=det(H)                      % detH determinant H
                                 % H, detH jsou v symbolickych promennych

% Delta2 a Delta1 ve stacionarnich bodech

detHX=subs(detH,{x,y},{X,Y})     % vypocte determinant v kazdem stac. bode
                                 % detHX je vektor symbolickych vyrazu

Delta2 = double(detHX)           %Delta2 uz jsou realna cisla

Delta1 = subs(H(1,1),{x,y},{X,Y}) % v symbolickych promennych
Delta1 = double(Delta1)           % vektor realnych cisel

% Zaver:  "usvedceni z extremu"
for kk = 1:length(X)             % pro kazdy stacionarni bod zvlast
    fprintf('\n v bodě  ');
    if Delta2(kk)<0
        fprintf('[%5.3f, %5.3f]$ NENÍ extrém, ', double(X(kk)), double(Y(kk)));
        fprintf('protože Delta_2 = %5.3f \n'   , Delta2(kk));
    end
    if Delta2(kk)>0
    fprintf('[%5.3f, %5.3f] JE extrém, ', double(X(kk)), double(Y(kk)));
        if Delta1(kk)>0
            fprintf('MINIMUM, hodnota %5.3f, \n', double(subs(f,{x,y},{X(kk),Y(kk)})));
            fprintf('protože Delta_2 = %5.3f a Delta_1 = %5.3f \n',Delta2(kk), Delta1(kk));
        else
            fprintf('MAXIMUM, hodnota %5.3f, \n', double(subs(f,{x,y},{X(kk),Y(kk)})));
            fprintf('protože Delta_2 = %5.3f a Delta_1 = %5.3f \n',Delta2(kk), Delta1(kk));
        end
    end

if abs(Delta2(kk)) <eps
    fprintf('NEVÍME\n');
end
end %for kk

%%%%%%% OBRAZEK %%%%%%%%%%%%%
figure();
hold on
grid on                                      % mrizka
fsurf(f)                                     % vykresli f
plot3(X,Y,double(subs(f,{x,y}, {X,Y})),'r*') % stacionarno body
title(texlabel(f))                           % titulek
hold off

% tvar lze upravit pomoci otaceni (tlacitko v grafickem okne cca uprostred)
% ulozeni do souboru aktualniho okna

% print('jmenosouboru', '-djpeg')