
N=100;
A=randn(N);A=A'*A;
tic;B=inv(sqrtm(A));toc
norm(B*B-inv(A),1)
tic;B2=sqrtm(inv(A));toc
norm(B2*B2-inv(A),1)
tic;[Q,D]=eig(A);B1=Q*diag(1./sqrt(diag(D)))*Q';toc
norm(B1*B1-inv(A),1)


#problem from chapter 5 of lecture notes (animal breeding)
clear
load rozd5data.txt


%1st solution create the normal system and solve it...
tic;
G=[X'*X,X'*Z;Z'*X,Z'*Z+k*inv(A)];
r=[X'*y;Z'*y];
s1=G\r
t1=toc

%LLSP insteat of normal system
tic;
S2=sqrt(k)*inv(sqrtm(A));
%norm(S2'*S2-k*inv(A),1)
B=[X,Z;zeros(rows(A),2), S2];
s2=B\[y;zeros(rows(A),1)]
t2=toc


%avoidinig sqrt of a matrix
tic;
S3=sqrt(k)*inv(chol(A))';
%norm(S3'*S3-k*inv(A),1)
B=[X,Z;zeros(rows(A),2), S3];
s3=B\[y;zeros(rows(A),1)]
t3=toc

%avoiding inverse of L (but we have t ocompute Z*L and 
tic;
L=chol(A)';
S4=sqrt(k)*eye(size(A));
B=[X,Z*L;zeros(rows(A),2),S4];
s4=B\[y;zeros(rows(A),1)];
s4(3:end)=L*s4(3:end)
t4=toc
differences=norm((ones(4,length(s1))*diag(s1))'-[s1,s2,s3,s4],1)
time=[t1;t2;t3;t4]
disp("time are not important as size is veeery small");
%try t ogenerate random large X,Z,A,y and apply all solutions