%% Time Varying Systems and Computation
% Klaus Diepold, Sven Gronauer, Simon Schumann
%
%% 3. Numerical Linear Algebra
%% 3.1 Near Singularity
% An n \: \text{x} \: n matrix \mathbf{A} is nonsingular if and only if its determinant is nonzero. For instance,

A=[1 2;2 4]
det(A)

% One might therefore expect that a small value for \text{det} (\mathbf{A}) indicates a matrix that is nearly singular. However, the size of \text{det} (\mathbf{A}) tells us nothing about near singularity. Indeed, since \text{det} (\alpha \mathbf{A}) =\alpha ^n \:\text{det} ( \mathbf{A}) , we can achieve any value for the determinant by multiplying by a scalar \alpha, yet \alpha \mathbf{A} is no more or less nearly singular than \mathbf{A} for  \alpha \neq 0.
% Another limitation of the determinant is shown by the two matrices

B=triu(ones(10))
C=B-diag(diag(B))+diag(-ones(10,1))

% Both matrices have unit diagonal and off-diagonal elements bounded in modulus by 1. So \text{det}(\mathbf{B}) = \text{det}(\mathbf{C}) =1:

det(B)
det(C)

% Yet, if we check their inverse matrices

inv(B)
inv(C)

% we realize that \mathbf{C} is ill conditioned. Thus, we can draw the conclusion that  determinant cannot distinguish between the ill-conditioned and well-conditioned matrices. To addres this issue, we can use condition number.

cond(B)
cond(C)

% The condition number of a matrix is a measure of how close the matrix is to being singular. The larger the condition number, the closer the matrix is to singularity.
%
%% 3.2 Evaluating Matrix Products in an Inefficient Order
% The cost of evaluating a matrix product depends on the order in which the product is evaluated (assuming the matrices are not all n \: \text{x} \: n). More precisely, matrix multiplication is associative, so (\mathbf{A} \mathbf{B} ) \mathbf{C} = \mathbf{A} (\mathbf{B} \mathbf{C}), and in general the cost of the evaluation of a product depends on where one puts the parentheses. One order may be much superior to others, so one should not simply evaluate the product in a fixed left-right or right-left order. For example, if \mathbf{x}, \mathbf{y}, \mathbf{z} and  are 10000-vectors then  \mathbf{x} \mathbf{y} ^\top \mathbf{z}  can be evaluated in two different ways:

x=100*ones(10000,1);
y=x;
z=y;
timeit(@()(x*y')*z)
timeit(@()x*(y'*z))

%
%% 3.3 Exploiting Structure in the Matrix
% One of the fundamental tenets of numerical linear algebra is that one should try to exploit any matrix structure that might be present. Sparsity  is particularly important to exploit, since algorithms intended for dense matrices may be impractical for sparse matrices because of extensive fill-in (zeros becoming nonzero).
% Ideally, linear algebra software would detect structure in a matrix and call an algorithm that exploits that structure. A notable example of such a meta-algorithm is the MATLAB backslash function:

D=ones(1000)+eye(1000)+[linspace(1,1000,1000)' zeros(1000,999)]; %1000 x 1000 dense matrix
S=eye(1000);                                                     %1000 x 1000 sparse matrix
x=linspace(1,1000,1000)';
v_1=D*x;
v_2=S*x;
timeit(@()(D\v_1))
timeit(@()(S\v_2))