Homotopy Test

n=1200;
A=RandomReal[{-1,1},{n,n}];
evals=Eigenvalues[A];
evals2D=Table[z=evals[[i]];{Re[z],Im[z]},{i,n}];
ListPlot[
evals2D,
AspectRatio->Automatic,
PlotRange->All
]

ListPlot[
evals/.{z_Complex:>{Re[z],Im[z]},
AspectRatio->Automatic,
PlotRange->All
]z_

SVD: Two ways

{n,m}={523,524};
SqrtSV=RandomVariate[NormalDistribution[],Min[n,m]];
Histogram[Data]

U=QRDecomposition[RandomReal[{-1,1},{n,n}]][[1]];
V=QRDecomposition[RandomReal[{-1,1},{m,m}]][[1]];
SV=SparseArray[Band[{1,1}]->(SqrtSV^2),{n,m}];
A=Transpose[U].SV.V;
Dimensions[A]
MatrixPlot[A]

SVD1=SingularValueList[A,Tolerance->(10^(-32))];
SVD2=(Sqrt[Eigenvalues[Transpose[A].A]]);
ListLogPlot[ {SVD1,SVD2, Reverse[Sort[(SqrtSV^2)]]},
PlotRange->All,
PlotStyle->{{PointSize[0.03],Red},{PointSize[0.015], Blue}, Green}]

SVD: Bi Diagonalization

Remember our formula for reflectors that zer vectors.

Householder[x_List,i_]:= Module[{n=Length[x],u,H, τ=Norm[x[[i;;-1]]]},
H=SparseArray[Band[{1,1}]->1,{n,n}];
τ=Sign[x[[i]]] τ; 
u=Join[Table[0,{i-1}],{1},((x[[i+1;;-1]])/(τ+x[[i]]))];
H-(2/((Norm[u])^2))Outer[Times,u,u]
]

Here it is working the way we had planned it to zero columns.

{n,m}={12,13};
AOriginal=RandomReal[{-1,1},{n,m}];
A=AOriginal;
x=A[[All,1]];
H=Householder[x,1];
MatrixPlot[Chop[H.A]]

Here it is zeroing rows! The integer says leave the stuff with indices higher than this alone.

{n,m}={12,13};
AOriginal=RandomReal[{-1,1},{n,m}];
A=AOriginal;
x=A[[1]];
H=Householder[x,2];
MatrixPlot[Chop[A.H]]

Putting them together and checking it does not change the SVs.

{n,m}={12,13};
AOriginal=RandomReal[{-1,1},{n,m}];
A=AOriginal;
x=A[[All,1]];
H=Householder[x,1];
A=Chop[H.A];
x=A[[1]];
H=Householder[x,2];
A=Chop[A.H];
MatrixPlot[A]
Norm[SingularValueList[A]-SingularValueList[AOriginal]]

Putting them together in a loop.

{n,m}={12,13};
AOriginal=RandomReal[{-1,1},{n,m}];
A=AOriginal;
Table[
x=A[[All,i]];
H=Householder[x,i];
A=Chop[H.A];
x=A[[i]];
H=Householder[x,i+1];
A=Chop[A.H];
MatrixPlot[A,PlotLabel->i],
{i,1,Min[n,m]-1}]
Norm[SingularValueList[A]-SingularValueList[AOriginal]]

{n,m}={4,23};
AOriginal=RandomReal[{-1,1},{n,m}];
A=AOriginal;
Table[
x=A[[All,i]];
H=Householder[x,i];
A=Chop[H.A];
x=A[[i]];
H=Householder[x,i+1];
A=Chop[A.H];
MatrixPlot[A,PlotLabel->i],
{i,1,Min[n,m]-1}]
Norm[SingularValueList[A]-SingularValueList[AOriginal]]