I have a simple Matlab script which aims to compute $k$ singular values of a matrix $A$. $A$ is a random dense square matrix of size $5000\times5000$, with 100 of its singular values constrained to be 0 (though that last detail does not seem to matter for my question).
I'm doing this in Matlab via
[Uk, Sk, Vk] = svds(A, k);. According to the documentation,
svds uses Lanczos bidiagonalization to compute these values. I looked at the function definition (
edit svds) and do not see any relevant branching, e.g. using different algorithms under the hood based on different conditions. However, when I increase $k$, I get very curious scaling/performance:
The docs mention
Increasing k can sometimes improve performance, especially when the matrix has repeated singular values.
But I interpret this to mean performance would be improved per $k$, rather than some huge reduction in total overall runtime.
Is this a known behavior of Lanczos bidiagonalization (an algorithm I'm not very familiar with)? Or does anyone have any speculations as to why the performance of
svds is like this?
Edit: Here is a minimal version of my script so others can try to reproduce:
results = ; A = rand(5000, 5000); [U, S, V] = svd(A); dS = diag(S); dS(4900:5000) = 0; A = U*diag(dS)*V; b = rand(5000, 1); for k = 100 : 100 : 4500 tic [Uk, Sk, Vk] = svds(A,k); Ahat = Vk*diag(1./diag(Sk))*Uk'; test = Ahat * b; time_k = toc results = [results; k time_k]; end plot(results(:,1), results(:,2))