From 4c964d34ba81b4cfe191e490fcd9ac5ea25a79a3 Mon Sep 17 00:00:00 2001 From: Pablo Zubieta Date: Mon, 19 Dec 2016 17:53:04 -0600 Subject: [PATCH] eigvecs(::Bidiagonal) speedup --- base/linalg/bidiag.jl | 15 +++++++++++---- 1 file changed, 11 insertions(+), 4 deletions(-) diff --git a/base/linalg/bidiag.jl b/base/linalg/bidiag.jl index 72aa768bda274..a3b4d72d1685e 100644 --- a/base/linalg/bidiag.jl +++ b/base/linalg/bidiag.jl @@ -559,23 +559,30 @@ function eigvecs{T}(M::Bidiagonal{T}) n = length(M.dv) Q = Array{T}(n, n) blks = [0; find(x -> x == 0, M.ev); n] + v = zeros(T, n) if M.isupper for idx_block = 1:length(blks) - 1, i = blks[idx_block] + 1:blks[idx_block + 1] #index of eigenvector - v=zeros(T, n) + fill!(v, zero(T)) v[blks[idx_block] + 1] = one(T) for j = blks[idx_block] + 1:i - 1 #Starting from j=i, eigenvector elements will be 0 v[j+1] = (M.dv[i] - M.dv[j])/M.ev[j] * v[j] end - Q[:, i] = v/norm(v) + c = norm(v) + for j = 1:n + Q[j, i] = v[j] / c + end end else for idx_block = 1:length(blks) - 1, i = blks[idx_block + 1]:-1:blks[idx_block] + 1 #index of eigenvector - v = zeros(T, n) + fill!(v, zero(T)) v[blks[idx_block+1]] = one(T) for j = (blks[idx_block+1] - 1):-1:max(1, (i - 1)) #Starting from j=i, eigenvector elements will be 0 v[j] = (M.dv[i] - M.dv[j+1])/M.ev[j] * v[j+1] end - Q[:,i] = v/norm(v) + c = norm(v) + for j = 1:n + Q[j, i] = v[j] / c + end end end Q #Actually Triangular