# Svdvals is inefficient for tall matrices

**URL:** https://discourse.julialang.org/t/svdvals-is-inefficient-for-tall-matrices/50330
**Category:** Performance
**Created:** [November 17, 2020, 8:54pm UTC](https://discourse.julialang.org/t/svdvals-is-inefficient-for-tall-matrices/50330 "2020-11-17T20:54:36Z")
**Posts on this page:** 3
**Page:** 1

<div class="post-metadata">

### Author: ![JeffFessler](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jefffessler/32/6650_2.png) [@JeffFessler](https://discourse.julialang.org/u/JeffFessler)
#### Post date: [November 17, 2020, 8:54pm UTC](https://discourse.julialang.org/t/svdvals-is-inefficient-for-tall-matrices/50330/1 "2020-11-17T20:54:36Z")

</div>

For a tall matrix X, calling `svdvals(X)` directly is noticeably slower than calling it on `svdvals(X'X)` despite the overhead of the extra matrix-matrix multiply:

```julia
using BenchmarkTools: @btime
using LinearAlgebra: svdvals
X = rand(ComplexF32, 1000,100)
f1(x) = svdvals(x) # tall
f2(x) = sqrt.(svdvals(x'x)) # small and square
@assert f2(X) ≈ f1(X) # yes, they match
@btime f1($X)
@btime f2($X)
versioninfo()

```

If this is consistent for tall matrices in general, should we build this behavior into `svdvals` directly? Or is the “squaring the condition number” concern the reason for keeping it as is?

Output:

```julia
  5.257 ms (9 allocations: 842.94 KiB)
  3.316 ms (12 allocations: 218.50 KiB)
Julia Version 1.5.3
Commit 788b2c77c1 (2020-11-09 13:37 UTC)
Platform Info:
  OS: Linux (x86_64-pc-linux-gnu)
  CPU: Intel(R) Xeon(R) CPU E5-2698 v4 @ 2.20GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-9.0.1 (ORCJIT, broadwell)
Environment:
  JULIA_NUM_THREADS = 80

```

Perhaps this question is somewhat related to

> [@svdvals is alarmingly slow](https://discourse.julialang.org/t/svdvals-is-alarmingly-slow/9259/21):
>
> @simonbyrne I had a chance to try it, and the results are not good: julia\> @btime svdvals($B); 1.371 ms (11 allocations: 138.20 KiB) julia\> versioninfo() Julia Version 0.6.2 Commit d386e40c17\* (2017-12-13 18:08 UTC) Platform Info: OS: Windows (x86\_64-w64-mingw32) CPU: Intel(R) Core(TM) i7-7820HQ CPU @ 2.90GHz WORD\_SIZE: 64 BLAS: libopenblas (USE64BITINT DYNAMIC\_ARCH NO\_AFFINITY Prescott) LAPACK: libopenblas64\_ LIBM: libopenlibm LLVM: libLLVM-3.9.1 (ORCJIT, broadwell) This is t…

---

<div class="post-metadata">

### Author: ![moeddel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moeddel/32/18641_2.png) [@moeddel](https://discourse.julialang.org/u/moeddel)
#### Post date: [November 18, 2020, 7:00am UTC](https://discourse.julialang.org/t/svdvals-is-inefficient-for-tall-matrices/50330/2 "2020-11-18T07:00:18Z")

</div>

What you found is a well known trick to speed up the calculation of the SVD for substantially wide or tall matrices as described [here](https://www.rdocumentation.org/packages/corpcor/versions/1.6.9/topics/fast.svd). We also made use of this in this topic:

> [@Very slow execution time in comparison even to Python](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/19):
>
> For data with only 3 dimensions I think it will be hard to make it faster, otherwise truncated svd or RandomizedLinAlg might be worth a try. Using MKL and Float32 instead of Float64 can speed up svd quite a bit.

Maybe it is worth it to implement a `fastsvd` function into `LinearAlgebra` as well.

---

<div class="post-metadata">

### Author: ![juthohaegeman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juthohaegeman/32/8620_2.png) [@juthohaegeman](https://discourse.julialang.org/u/juthohaegeman)
#### Post date: [November 19, 2020, 2:46pm UTC](https://discourse.julialang.org/t/svdvals-is-inefficient-for-tall-matrices/50330/3 "2020-11-19T14:46:08Z")

</div>

If you make `X'*X` you want to compute a symmetric eigenvalue decomposition (`eig(X'*X)`) instead of an SVD of `X'*X`, that should be even faster because it knows that in that case there is only one instead of two unitary matrices involved. In principle this approach is less accurate.

A better approach is to to first do a thin QR decomposition of X and then doing an SVD of the square and small matrix R, such that you have  
`X = QR = (Q U) S V'`

e.g., if you only want the singular values: `f3(x) = svdvals(qr(X).R)`.

This is faster than directly `svdvals(X)` (which I find strange, I would think Lapack should be smart enough to do this itself), but slower than the appraoch based on `X'*X`. The QR decomposition is not that more costly as matrix multiplication `X'*X`, but I assume that the latter is just much more optimized.
