%------------------------------------------------------------------------
%%%               ENGINEERING MATHEMATICS III - MATH 2Z03
%%%                               LAB #5
%%%                 LINEAR ALGEBRA: EIGENVALUE PROBLEMS
%------------------------------------------------------------------------


%------------------------------------------------------------------------
% Covers:
% - "Numerical Mathematics" by M. Grasselli and D. Pelinovsky,
%    Sections 4.1-4.3
%------------------------------------------------------------------------


%------------------------------------------------------------------------
% Instructions:
% - Execute the commands given below in MATLAB. You can do this in
%   one of the following three ways:
%   (i) retype the commands one by one at the MATLAB command line (>>),
%   (ii) copy and paste the commands one by one to the MATLAB
%        command line (>>),
%   (iii) run all the commands at once by calling the present file
%         as a script, i.e. >> m2z03lab
% - Analyze the output and compare it with the comments embedded in
%   the file below.
%------------------------------------------------------------------------

clc; close all; clear all;

%------------------------------------------------------------------------
% Section 4.1: Matrix Eigenvalue Problems
disp(' -----------------------------------------------------------')
disp(' Section 4.1: Matrix Eigenvalue Problems ')

%%% Problem: Find all the eigenvalues and eigenvectors of an n by n
%%%          complex matrix A

% Eigenvalues and eigenvectors of square matrices can be computed with
% the MATLAB function 'eig'

% Use command 'help' to obtain information about MATLAB function 'eig'
help eig

%%% Example
% Square matrix A
A = [2 4 2; 6 0 8; 3 2 -8]

% Computing the vector of eigenvalues (lambda)
lambda = eig(A)

% Computing the vector of eigenvalues
% (represented as a diagonal matrix lambda)
% and the matrix of eigenvectors V
[V,lambda] = eig(A)

disp('  Press any key to continue ... ')
pause
%------------------------------------------------------------------------


%------------------------------------------------------------------------
% Section 4.2: Properties of Eigenvalues
disp(' -----------------------------------------------------------')
disp(' Section 4.2: Properties of Eigenvalues ')

%%% Algebraic multiplicity
A = [1 2 2 -2; 0 2 0 -3; 3 -6 -4 6; -1 2 2 -3]
lambda = eig(A)'

% algebraic multiplicity (1) for eigenvalue lambda = -7
null((A-lambda(1)*eye(4))^4)

% algebraic multiplicity (2) for eigenvalue lambda = 2
null((A-lambda(2)*eye(4))^4)

% algebraic multiplicity (1) for eigenvalue lambda = -1
null((A-lambda(3)*eye(4))^4)

%%% Block-triangular matrices
% Eigenvalues of matrices (diagonal blocks) C22 and C33 coincide with
% some eigenvalues of the matrix C (block triangular):
C = [4 3 -2 5 7; 0 5 3 6 9; 0 -5 2 1 -7;0 0 0 4 3; 0 0 0 -2 -1]
C22 = [5 3; -5 2]
C33 = [4 3; -2 -1]

% Compare:
eig(C)

eig(C22)

eig(C33)

disp('  Press any key to continue ... ')
pause
%------------------------------------------------------------------------


%------------------------------------------------------------------------
% Section 4.3: Properties of Eigenvectors
disp(' -----------------------------------------------------------')
disp(' Section 4.3: Properties of Eigenvectors ')

%%% If matrix A has distinct eigenvalues, then the corresponding
% non-zero eigenvectors are linearly independent.

% Matrix A has three distinct eigenvalues and linearly independent eigenvectors
A = [2 4 2; 6 0 8; 3 2 -8]

% The corresponding matrix of eigenvectors has full rank three
[V,lambda] = eig(A)

diag(lambda)'

rank(V)

%%% Geometric multiplicity
A = [2 6 -15; 1 1 -5; 1 2 -6]

eig(A)'

% lambda = -1 is the only eigenvalue of A with algebraic multiplicity
% equal to 3 and geometric multiplicity equal to 2, meaning that A has
% only two linearly independent eigenvectors.
null((A+eye(3))^3)

null(A+eye(3))


%%% Diagonalization
% The diagonal matrix of eigenvalues in the diagonalization formula
% is precisely the second output of the MATLAB function 'eig':
A = [2 4 5; 3 6 7; -1 4 -2]

[V,lambda] = eig(A)

% The diagonal matrix of eigenvalues
inv(V)*A*V

disp('  This is the end of the lab tutorial #5. ')