# Calculating the \`k\` smallest positive eigenvalues of a matrix with ARPACK

**URL:** <https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864>\
**Category:** General Usage\
**Created:** [September 17, 2019, 7:37pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864 "2019-09-17T19:37:42Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 17, 2019, 7:37pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/1 "2019-09-17T19:37:42Z")

</div>

I have a sqaure matrix that is not necessarily positive definite, for example `rand(10000,10000)`, and I’m trying to obtain the `k` smallest positive eigenvalues with the command `eigs()` from `ARPACK`. How can I achieve that provided that the smallest eigenvalues are negative and I want to select the `k` smallest eigenvalues only among the positive eigenvalues? The options `which=:SM` and `which=:SR` don’t do the trick. Thank you!

---

<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:** [September 17, 2019, 9:02pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/2 "2019-09-17T21:02:14Z")

</div>

Interior eigenvalues are essentially impossible with Krylov based methods if you can only multiply with the matrix, you’ll want to use shift and invert for this (or in this case only invert, since you are looking for eigenvalues closest to zero, so the shift is zero). Thus, you need to be able to efficiently apply the inverse of that matrix.

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 12:04am UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/3 "2019-09-18T00:04:50Z")

</div>

Ok 🙂 Thanks. Will try to find another way then

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [September 18, 2019, 12:08am UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/4 "2019-09-18T00:08:51Z")

</div>

If your matrix is dense but small enough to store, e.g. `rand(10000,10000)`, you might as well use `eigvecs`. Really large matrices usually have some structure (e.g. sparsity) that allow you to solve `A \ b` more quickly (e.g. sparse-direct methods or iterative solvers), which allow you to use shift-and-invert or other methods (e.g. [contour-integration](http://www.ecs.umass.edu/~polizzi/feast/) methods, which also require shifted linear solves).

---

<div class="post-metadata">

**Author:** ![pablosanjose](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pablosanjose/32/7006_2.png) [@pablosanjose](https://discourse.julialang.org/u/pablosanjose)\
**Post date:** [September 18, 2019, 9:06am UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/5 "2019-09-18T09:06:09Z")

</div>

ARPACK’s `eigs` implements shift and invert AFAIU, but you have to provide a `sigma = ...` option, indicating the shift (zero in your case, I believe). Don’t use `which`.

ARPACK will not focus on one side of the shift point, however. Have a look at LOBPCG in IterativeSolvers.jl for this. Also check out the other native-Julia packages for these things, e.g KrylovKit.jl and ArnoldiMethod.jl

EDIT: example

```julia
julia> using Arpack

julia> r = rand(10000, 10000); r += transpose(r);

julia> eigs(r, sigma = 0)[1]
6-element Array{Float64,1}:
 -0.0025761027840885565
  0.021348703719922396 
 -0.02638392700654548  
  0.030827180680079686 
 -0.044944020253078275 
 -0.053163665138762084 

```

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [September 18, 2019, 9:14am UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/6 "2019-09-18T09:14:38Z")

</div>

Depends on how many negative eigenvalues you have. If you just have a few and want the positive among these, just use ARPACK or LOBPCG with a large number of sought eigenvectors. If you have too many of them so that you can’t afford to compute them all, then yeah it’s an interior eigenproblem and you’re basically screwed unless you can afford efficient inverses. There’s also a “shift-and-square” method that squares the matrix A to transform the eigenvalues close to 0 into exterior ones, but I wouldn’t try it unless you’re really desperate.

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 2:21pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/7 "2019-09-18T14:21:25Z")

</div>

I want to avoid the use of `eigvecs` because I only need a few eigenvalues and not all eigenvalues.

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 2:22pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/8 "2019-09-18T14:22:34Z")

</div>

Ok, that sounds good! Will try to explore your suggestions!!

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [September 18, 2019, 6:12pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/9 "2019-09-18T18:12:05Z")

</div>

> [@sergevic](#):
>
> I want to avoid the use of `eigvecs` because I only need a few eigenvalues and not all eigenvalues.

I understand, but if you have a dense matrix with no structure it might not be any slower to just compute all the eigenvalues than to try to use an iterative method — you can’t get around the fact that you need O(n³) operations for interior eigenvalues, so you are just playing with the constant factor. If your matrices are very large, it is advisable to look for some special structure to exploit.

(One thing you could try to do is to compute the Hessenberg factorization of your dense matrix, and then use the Hessenberg factorization for shift-and-invert solves, which have [new efficient solvers in Julia 1.3](https://github.com/JuliaLang/julia/pull/31853) — this is useful if you need a lot of different shifts for the same matrix, but doesn’t get around the O(n³) scaling.)

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 6:49pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/10 "2019-09-18T18:49:59Z")

</div>

I think the `KrylovKit.jl` will do the trick, but I’m stuck with the arguments options of the function `eigsolve().` The documentation is very complete but, as I always say when I criticize the documentation of Julia: there are no examples with the commands. For example, the simple code `eigsolve(rand(1000,1000), howmany=15, which=:SM)` doesn’t work and throw the error `got unsupported keyword argument "howmany"` even if `homany` and `which` are the arguments of the function. I even tried `eigsolve(rand(1000,1000), [howmany=15, which=:SM])` with no success. I’m sure the option `EigSorter(f; rev = false)` will do the trick for what I want, instead of `which=:SM`, but if the command doesn’t work with the basic arguments, then with this one it will be difficult to implement it…

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 6:52pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/11 "2019-09-18T18:52:20Z")

</div>

Yes, good point indeed! And thanks for the Hessenberg factorization!

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [September 18, 2019, 6:59pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/12 "2019-09-18T18:59:33Z")

</div>

Documentation getting out of date is an unfortunate fact of life (and it’s a good example of why the mantra “always document everything a lot” is not always the right thing to do…), and the only thing you can do is to open an issue (or even better submit a PR). That said, Julia packages have good tests, which often double up as examples. It tends to rot less fast, as CI yells at people when tests break; in this case [https://github.com/Jutho/KrylovKit.jl/blob/master/test/eigsolve.jl](https://github.com/Jutho/KrylovKit.jl/blob/master/test/eigsolve.jl) should have what you’re looking for.

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 7:01pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/13 "2019-09-18T19:01:21Z")

</div>

I agree totally with you! Thanks for the link 🙂 Github pages are always a better source for that!

---

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [September 18, 2019, 7:12pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/14 "2019-09-18T19:12:47Z")

</div>

Be careful to the caracter `;` howmany is not a keyward argument!! It should be

```julia
eigsolve(rand(1000,1000), 15, :SM) 

```

Also, you should be carefull that to target specific parts of the spectrum, shift-invert methods are required, `:SM` will not do the trick (completely).

> And `:SM` is not valid…

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 7:18pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/15 "2019-09-18T19:18:31Z")

</div>

Ok, perfect! 🙂 I think I’m close to the solution! I’m trying to use the option `EigSorter` to get the `k` positive eigenvalues I want and tried `essai=eigsolve(rand(1000,1000),15,EigSorter(function pos(x) x[x.>0] end))` , but still have an error! I’m still motivated to find the solution! 😛 Perhaps the problem is that the eigenvalues I want are not close to the periphery of the spectrum, and so perhaps the better is to compute all the eigenvalues, as @stevengj was saying, despite I want to avoid that …

---

<div class="post-metadata">

**Author:** ![rveltz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rveltz/32/2707_2.png) [@rveltz](https://discourse.julialang.org/u/rveltz)\
**Post date:** [September 18, 2019, 7:32pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/16 "2019-09-18T19:32:59Z")

</div>

It depends on the size of the matrix, you can use a Cayley transform or the FEAST algorithm. Roughly speaking, power methods for `M` gives you the largest EV, so the same method for `M^{-1}` will give the smallest ones.

Linear algebra is not simple as soon as the “matrix” is large, see @stevengj, @juthohaegeman, @pablosanjose answers above. Hence, you cannot expect to use a general method to work in all cases. Rather you should use the specificities of the matrix.

---

<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:** [September 18, 2019, 8:02pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/17 "2019-09-18T20:02:22Z")

</div>

Is your matrix hermitian (or, since it is likely real, symmetric)? In that case you can indeed square it and look for the `:SR` eigenvalues as suggested above. I did not think of this as you used `rand(10000,10000)` as an example in your initial post. Squaring has a negative impact on convergence speed (it squares the ratio between the gaps in the spectrum), but at least gives you a way to move eigenvalues from the interior to the peripheral. If it is not symmetric and eigenvalues will be complex, you cannot expect squaring the matrix will help.

For using `EigSorter`, the function `f` should be a scalar function, it’s really as using `sort(..., by = f, rev = ...)`, e.g. `:LM` is `EigSorter(abs, rev=true)` , `:SR` is `EigSorter(real)`.

Any improvements in the documentation are welcome.

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 8:31pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/18 "2019-09-18T20:31:13Z")

</div>

Yes I’m aware of that… You are right…

---

<div class="post-metadata">

**Author:** ![sergevic](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sergevic/32/9958_2.png) [@sergevic](https://discourse.julialang.org/u/sergevic)\
**Post date:** [September 18, 2019, 8:32pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/19 "2019-09-18T20:32:14Z")

</div>

It’s a real and symmetric matrix. The `rand(10000,10000)` was just an example to give a reproducible code.

---

<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:** [September 18, 2019, 8:41pm UTC](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864/20 "2019-09-18T20:41:16Z")

</div>

Then I would use the squaring idea, unless indeed there are only a few negative eigenvalues so that you can just query sufficient eigenvalues to also get the first positive ones.

[Next page](https://discourse.julialang.org/t/calculating-the-k-smallest-positive-eigenvalues-of-a-matrix-with-arpack/28864.md?page=2)