% P1
randn('seed',0)
rand('seed',0)
% Example 1
clc
disp('Example 1.')
% You should get good results here even without reorthogonalization.
n = 10;
d = randn(n,1); d = sort(d);
q = rand(n,1); q = q/norm(q);
[alpha,beta] = MakeTridiag(d,q);
T = diag(alpha) + diag(beta,-1) + diag(beta,1);
[Q,D] = eig(T);
format long
clc
TheVecs = abs([Q(n,:)' q])
TheEigs = [diag(D) d]

disp(' ')
disp('Example 2.')
% You'll get poor results here without reorthogonalization...
n = 100;
d = randn(n,1); d = sort(d);
q = rand(n,1); q = q/norm(q);
[alpha,beta] = MakeTridiag(d,q);
T = diag(alpha) + diag(beta,-1) + diag(beta,1);
[Q,D] = eig(T);
format long
CosVecs = abs(Q(n,:))*abs(q)
EigDiff = norm(diag(D)-d)
TheEigs = [diag(D) d];
FirstTenEigs = TheEigs(1:10,:)