In [1]:
using LinearAlgebra
using Random
In [2]:
##
Random.seed!(0)
A = randn(5,5)
A = A+A'
Ainit = copy(A)
Out[2]:
5×5 Array{Float64,2}:
  1.35821    1.12575   -1.04191   -0.322427  0.46948 
  1.12575    0.129895  -0.871821  -2.12147   0.973077
 -1.04191   -0.871821   0.794965  -1.66916   0.795922
 -0.322427  -2.12147   -1.66916    4.55247   0.131077
  0.46948    0.973077   0.795922   0.131077  0.558933
In [3]:
## Example of the deflation method
λs,X = eigen(A)
λi = rand(1:size(A,2))
x = X[:,λi] # just an eigenvector
# build the householder matrix H
xe1 = zeros(size(A,1))
xe1[1] = norm(x)
v = xe1-x
H = I - 2*v*v'/(v'*v)
# apply it
display(H*A*H')
display(λs[λi])

5×5 Array{Float64,2}:
 -0.257679     -1.10042e-15   2.84239e-15  -2.59207e-15  -6.46266e-17
 -1.23125e-15   0.929135     -1.59581      -2.20476       1.33736    
  2.96355e-15  -1.59581       1.45064      -1.59513       0.465633   
 -2.47452e-15  -2.20476      -1.59513       4.54801       0.0903187  
 -8.54693e-17   1.33736       0.465633      0.0903187     0.724375   
-0.25767946156596877
In [4]:
##
eigvals((H*A*H')[2:end,2:end])
Out[4]:
4-element Array{Float64,1}:
 -2.177936876216415
  1.123268282931314
  2.904571788931462
  5.802249203273017
In [5]:
## The subspace method coverges
function subspace(A,k)
  X = Matrix(1.0I,size(A)...)
  D = Diagonal(sort(diag(X'*A*X)))
  for i = 1:k
    Y = A*X
    D = Diagonal(sort(diag(X'*Y)))
    X,R = qr(Y)
  end
  return X,D
end
X,D = subspace(Ainit,101)
display([diag(D) eigvals(Ainit)])

5×2 Array{Float64,2}:
 -2.17794   -2.17794 
 -0.257679  -0.257679
  1.12327    1.12327 
  2.90457    2.90457 
  5.80225    5.80225
In [6]:
##
X = Matrix(1.0I,size(A)...)
A = copy(Ainit)
for i = 1:101
  Q,R = qr(A)
  global A = R*Q
  global D = Diagonal(sort(diag(A)))
  #global X = X*Q
end
display([diag(D) sort(eigvals(A))])

5×2 Array{Float64,2}:
 -2.17794   -2.17794 
 -0.257679  -0.257679
  1.12327    1.12327 
  2.90457    2.90457 
  5.80225    5.80225
In [7]:
## Let's test triadiagonal steps
Ainit = Tridiagonal(randn(8,8))
Ainit = Ainit + Ainit'
Q,R = qr(Matrix(Ainit))
R*Q
Out[7]:
8×8 Array{Float64,2}:
  1.29797  -2.14369   8.88178e-16  …  -1.70296e-18  -1.80848e-18
 -2.14369  -1.2712    2.35306          7.44467e-18  -5.60355e-18
  0.0       2.35306   0.421318         1.88353e-18   2.58376e-19
  0.0       0.0      -0.167794         8.78622e-18  -7.21831e-18
  0.0       0.0       0.0             -1.94289e-16   1.2081e-16 
  0.0       0.0       0.0          …   0.734347      2.22045e-16
  0.0       0.0       0.0             -0.0122276     0.58023    
  0.0       0.0       0.0              0.58023       0.599098
In [8]:
##
A = copy(Ainit)
for i = 1:10
  Q,R = qr(A)
  A = R*Q
  D = diagm(sort(diag(A)))
  X = X*Q
  display(A)
end

MethodError: no method matching qr!(::SparseArrays.SparseMatrixCSC{Float64,Int64})
Closest candidates are:
  qr!(!Matched::Union{DenseArray{#s617,2}, Base.ReinterpretArray{#s617,2,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray}, Base.ReshapedArray{#s617,2,A,MI} where MI<:Tuple{Vararg{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray}, SubArray{#s617,2,A,I,L} where L where I<:Tuple{Vararg{Union{Int64, AbstractRange{Int64}, Base.AbstractCartesianIndex},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, Base.ReshapedArray{T,N,A,MI} where MI<:Tuple{Vararg{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, DenseArray}} where #s617<:Union{Complex{Float32}, Complex{Float64}, Float32, Float64}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.2/LinearAlgebra/src/qr.jl:244
  qr!(!Matched::Union{DenseArray{#s617,2}, Base.ReinterpretArray{#s617,2,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray}, Base.ReshapedArray{#s617,2,A,MI} where MI<:Tuple{Vararg{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray}, SubArray{#s617,2,A,I,L} where L where I<:Tuple{Vararg{Union{Int64, AbstractRange{Int64}, Base.AbstractCartesianIndex},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, Base.ReshapedArray{T,N,A,MI} where MI<:Tuple{Vararg{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, DenseArray}} where #s617<:Union{Complex{Float32}, Complex{Float64}, Float32, Float64}, !Matched::Val{false}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.2/LinearAlgebra/src/qr.jl:242
  qr!(!Matched::Union{DenseArray{#s617,2}, Base.ReinterpretArray{#s617,2,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray}, Base.ReshapedArray{#s617,2,A,MI} where MI<:Tuple{Vararg{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray}, SubArray{#s617,2,A,I,L} where L where I<:Tuple{Vararg{Union{Int64, AbstractRange{Int64}, Base.AbstractCartesianIndex},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, Base.ReshapedArray{T,N,A,MI} where MI<:Tuple{Vararg{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},N} where N} where A<:Union{Base.ReinterpretArray{T,N,S,A} where S where A<:Union{SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, SubArray{T,N,A,I,true} where I<:Union{Tuple{Vararg{Real,N} where N}, Tuple{AbstractUnitRange,Vararg{Any,N} where N}} where A<:DenseArray where N where T, DenseArray} where N where T, DenseArray}} where #s617<:Union{Complex{Float32}, Complex{Float64}, Float32, Float64}, !Matched::Val{true}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.2/LinearAlgebra/src/qr.jl:243
  ...

Stacktrace:
 [1] qr(::Tridiagonal{Float64,Array{Float64,1}}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.2/LinearAlgebra/src/qr.jl:372
 [2] top-level scope at ./In[8]:4