(*Definition of System*)

mu[x_]:=0.954+(x-0.5)*0.5

h=0.05

A={{-2/mu[h],1/mu[h],0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{1/mu[2*h],-2/mu[2*h],1/mu[2*h],0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,1/mu[3*h],-2/mu[3*h],1/mu[3*h],0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,1/mu[4*h],-2/mu[4*h],1/mu[4*h],0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,1/mu[5*h],-2/mu[5*h],1/mu[5*h],0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,1/mu[6*h],-2/mu[6*h],1/mu[6*h],0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,1/mu[7*h],-2/mu[7*h],1/mu[7*h],0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,1/mu[8*h],-2/mu[8*h],1/mu[8*h],0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,1/mu[9*h],-2/mu[9*h],1/mu[9*h],0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,1/mu[10*h],-2/mu[10*h],1/mu[10*h],0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,1/mu[11*h],-2/mu[11*h],1/mu[11*h],0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,1/mu[12*h],-2/mu[12*h],1/mu[12*h],0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,1/mu[13*h],-2/mu[13*h],1/mu[13*h],0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,1/mu[14*h],-2/mu[14*h],1/mu[14*h],0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,1/mu[15*h],-2/mu[15*h],1/mu[15*h],0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,1/mu[16*h],-2/mu[16*h],1/mu[16*h],0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1/mu[17*h],-2/mu[17*h],1/mu[17*h],0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1/mu[18*h],-2/mu[18*h],1/mu[18*h]},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1/mu[19*h],-2/mu[19*h]}}

Eigenvalues[A]

(*A=Transpose[Q].R -- Note: different from definition in class with Q->Q^T*)

{qq, rr}=QRDecomposition[A]

AA=rr.Transpose[qq]

(*Eigenvalue check using pre-defined function in Mathematica*)
(*Eigenvalues were not altered by one iterative step of QR + switch*)

Eigenvalues[AA]

(*To get eigenvalues using QR, repeat the above two steps multiple of times*)

(* Initialization*)
AA=A;err=10;counter=1;tol=10^-6;

(* A simple iterative loop to do the QR + Switch steps*)
While[ err>tol,
oldval=AA[[1]][[1]];
{qq,rr}=QRDecomposition[AA];
AA=rr.Transpose[qq];
err=Abs[oldval-AA[[1]][[1]]];
counter++;
Print[counter]
]

(*it took less than 100 iterations to attain the convergence fo the largest  diagonal element (in magnitude) in A to the desired tol*)

(* Assigning diagonal elements as the eigenvalues*)
For[i=1,i≤19,i=i+1,
eigen[i]=AA[[i]][[i]];
Print[eigen[i]]]