%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Module V: Quadratic Form, Characteristic Roots, Matrix Decomposition
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% NOTE: to see results in the Command Window, 
% simply delete the ";" at the end of a given command.

%reconsider the following square, symmetric matrix:
M=[1 4 7;4 2 6;7 6 3];

% get its eigenvalues and eigenvectors
[C L]=eig(M); % L has some negative elements, so M is not pos. def.

%Verify spectral decomposition rule:
test=C*L*C'; %this returns M (with perhaps some minute rounding errors)

% An attempt to perform a Cholesky transformation on M will fail:
% Mc=chol(M); % you will need to comment-out this line to run the rest of your script

% However, if M is full rank, we know that M'M will be pos. def. Let's
% first check on rank:
rM=rank(M); % OK - full rank

M2=M'*M;
[C L]=eig(M2); % all char. roots are positive
L=chol(M2);

% Verify Cholesky decomposition rule:
test=L'*L; %this yields the original M2

% Practice:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
A=[20 4 3 1;4 8 6 2;3 6 5 7;1 2 11 4];
A1=A'*A;
[C L]=eig(A1);

%check A-80:
test1=sum(sum(A1*C - C*L)); 
%Instead of inspecting if all elements of this difference are zero, 
% we can take the row sum of the column sum & check if this scalar is zero.
% It should be close to zero, short of rounding errors. You will get a
% number such as 3.2663e-013, which means it is zero to the 13th decimal -
% good enough!

%check A-84:
test2=sum(sum(C'*A1*C-L));

%check A-85:
test3=sum(sum(C*L*C'-A1));

d=det(A1);
r=rank(A1);