% Copyright (c) 2017, Domenico L. Gatti
% All rights reserved.
% 
% Redistribution and use in source and binary forms, with or without 
% modification, are permitted provided that the following conditions are 
% met:
% 
%     * Redistributions of source code must retain the above copyright 
%       notice, this list of conditions and the following disclaimer.
%     * Redistributions in binary form must reproduce the above copyright 
%       notice, this list of conditions and the following disclaimer in 
%       the documentation and/or other materials provided with the 
%       distribution
%       
% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
% IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
% THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
% PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
% CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
% EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
% PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
% PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
% LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
% NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
% SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
%
%% General dependencies
% We always start from the CODE directory and we add to the path
% subdirectories containing various tools described in the book chapters.
addpath(genpath('../GENERAL_SCRIPTS_FUNCTIONS'));
% addpath(genpath('../DATABASE'));
% addpath(genpath('../TOOLBOXES'));

%% CHAPTER 9: Eigenvalues and eigenvectors
clear, clc
M = randi(10,4)
[S,D] = eig(M)
S*D*inv(S), S*D/S
trace(M)
sum(diag(D))

[~,A] = gauss_elim_step_by_step(M,[1 2 3 4]')
prod(diag(A))
prod(diag(D))
det(M)

[L,U] = lu(M)
s = det(L)
s*prod(diag(U))

M
M2 = M + 0.1*eye(4)
[S2,D2] = eig(M2)
S
S2
diag(D)
diag(D2)
diag(D2)-diag(D)

M_ut = triu(M)
[S_ut,D_ut] = eig(M_ut)
M_lt = tril(M)
[S_lt,D_lt] = eig(M_lt)
M_diag = diag(diag(M))
[S_diag,D_diag] = eig(M_diag)

eig(M)
eig(M')
eig(inv(M))
1./eig(M)

S = rand(4)
D1 = diag(rand(4,1))
M1 = S*D1*inv(S)
D2 = diag(rand(4,1))
M2 = S*D2*inv(S)
M1*M2
M2*M1
eig(M1*M2)
eig(M2*M1)
D1+D2
eig(M1+M2)

M1 = randi(10,4)
M1 = M1'*M1
[S,D1] = eig(M1)
D2 = diag(rand(4,1))
M2 = S*D2*inv(S)
M1*M2
M2*M1
eig(M1*M2)
eig(M2*M1)
D1+D2
eig(M1+M2)

M1 = rand(4)
M1 = M1'*M1
[S,D] = eig(M1)
M1^4
S*D^4/S
expm(M1)
S*expm(D)/S
S*diag(exp(diag(D)))/S

(exp(1))^M1
5^M1
S*diag(5.^diag(D))/S

M3 = rand(4)
M3 = M3*M3'
M3 = M3'*M3
[S,D] = eig(M3)
S'*S
S*S'
S*D*S'

M3 = rand(4)
M3 = triu(M3)-triu(M3)'+diag(diag(M3))
[S,D] = eig(M3)
S'*S
S*S'

Q1 = [cos(pi/2) sin(pi/2) 0;...
    -sin(pi/2) cos(pi/2) 0;...
    0 0 1]
Q1 = [1 0 0;...
      0 cos(pi/2) sin(pi/2) ;...
      0 -sin(pi/2) cos(pi/2)]
[S,D] = eig(Q1)

S'*S
S*S'

A = rand(4)
A = A'*A
[S,D] = eig(A)

