%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                            Michael Pokojovy                             % 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[X, Y] = meshgrid(-2:0.1:2);

figure(1);

set(gcf, 'PaperUnits', 'centimeters');
xSize = 14; ySize = 10;
xLeft = (21 - xSize)/2; yTop = (30 - ySize)/2;
set(gcf,'PaperPosition', [xLeft yTop xSize ySize]);
set(gcf,'Position',[0 0 xSize*50 ySize*50]);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
hold on;
t = title('Ellipsoid: $(x - c)^{T} A^{-1} (x - c) \leq 1$', 'interpreter', 'latex');
set(t, 'FontSize', 16);

A = [ 5 -2 0;
     -2 10 4;
      0  4 6];
 
c = [-1; 3; 2];
 
[S, D] = eig(A);

N = 20;

theta = linspace(0, pi, N);
phi   = linspace(0, 2*pi, N);

X = zeros(N, N);
Y = zeros(N, N);
Z = zeros(N, N);

for i = 1:N
for j = 1:N
    X(i, j) = sin(theta(i))*cos(phi(j));
    Y(i, j) = sin(theta(i))*sin(phi(j));
    Z(i, j) = cos(theta(i));
    
    P = [X(i, j); Y(i, j); Z(i, j)];
    P = S*sqrt(D)*S'*P + c;
    X(i, j) = P(1);
    Y(i, j) = P(2);
    Z(i, j) = P(3);
end
end

colormap gray;
grid on;

surf(X, Y,  Z);
view([-51 38]);

h = text(0, 0, 6, '$A = \left(\begin{array}{ccc} 5 & -2 & 0 \\ -2 & 10 & 4 \\ 0 & 4 & 6 \end{array}\right)$');
set(h, 'FontSize', 16,'interpreter', 'latex');
h = text(0, -1, 3, '$c = \left(\begin{array}{c} -1 \\ 3 \\ 2 \end{array}\right)$');
set(h, 'FontSize', 16,'interpreter', 'latex');

xlabel('x');
ylabel('y');
zlabel('z');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%