
vn = [100:100:1000];
%vn = [5:15];
nn = length(vn);
tempo1 = zeros(nn,1);
condiz1 = zeros(nn,1);
err1 = zeros(nn,1);
tempo2 = zeros(nn,1);
err2 = zeros(nn,1);
tempo3 = zeros(nn,1);
err3 = zeros(nn,1);
resid = zeros(nn,1);
alpha = 1;
ds = 1e-8;

for i = 1:nn
	n = vn(i)
	m = 2*n;
	A = rand(m,n);
	%A = hilb(m);
	%A = A(:,1:n);
	sol = ones(n,1);
	b = A*sol;
	[Q R] = qr(A);
	b = b+alpha*Q(:,n+1);
	b = b + ds*randn(m,1);
	condiz1(i) = cond(A);
	resid(i) = norm(A*sol-b);
	% backslash
	tic
	x = A\b;
	tempo1(i) = toc;
	err1(i) = norm(x-sol);
	% equazioni normali
	tic
	G = A'*A;
	R = chol(G);
	c = A'*b;
	y = R'\c;
	x = R\y;
	tempo2(i) = toc;
	err2(i) = norm(x-sol);
	% fattorizzazione QR
	tic
	[Q R] = qr(A);
	c = Q'*b;
	R1 = R(1:n,:);
	c1 = c(1:n);
	x = R1\c1;
	tempo3(i) = toc;
	err3(i) = norm(x-sol);
end
figure(1)
semilogy(vn',[tempo1 tempo2 tempo3])
title('tempi di calcolo')
legend('backslash','eq. normali','fatt. QR','location','southeast')
figure(2)
semilogy(vn',[err1 err2 err3])
title('errori')
legend('backslash','eq. normali','fatt. QR','location','southeast')


