Arpack

Arpack

This package provides bindings to ARPACK, which can be used to perform iterative solutions for eigensystems (using eigs) or singular value decompositions (using svds).

eigs calculates the eigenvalues and, optionally, eigenvectors of its input(s) using implicitly restarted Lanczos or Arnoldi iterations for real symmetric or general nonsymmetric matrices respectively.

For the single matrix version,

eigs(A; nev=6, ncv=max(20,2*nev+1), which=:LM, tol=0.0, maxiter=300, sigma=nothing, ritzvec=true, v0=zeros((0,))) -> (d,[v,],nconv,niter,nmult,resid)

the following keyword arguments are supported:

whichtype of eigenvalues
:LMeigenvalues of largest magnitude (default)
:SMeigenvalues of smallest magnitude
:LReigenvalues of largest real part
:SReigenvalues of smallest real part
:LIeigenvalues of largest imaginary part (nonsymmetric or complex A only)
:SIeigenvalues of smallest imaginary part (nonsymmetric or complex A only)
:BEcompute half of the eigenvalues from each end of the spectrum, biased in favor of the high end. (real symmetric A only)

We can see the various keywords in action in the following examples:

julia> A = Diagonal(1:4);

julia> λ, ϕ = eigs(A, nev = 2, which=:SM);

julia> λ
2-element Array{Float64,1}:
 1.0000000000000002
 2.0

julia> B = Diagonal([1., 2., -3im, 4im]);

julia> λ, ϕ = eigs(B, nev=1, which=:LI);

julia> λ
1-element Array{Complex{Float64},1}:
 1.3322676295501878e-15 + 4.0im

julia> λ, ϕ = eigs(B, nev=1, which=:SI);

julia> λ
1-element Array{Complex{Float64},1}:
 -2.498001805406602e-16 - 3.0000000000000018im

julia> λ, ϕ = eigs(B, nev=1, which=:LR);

julia> λ
1-element Array{Complex{Float64},1}:
 2.0000000000000004 + 4.0615212488780827e-17im

julia> λ, ϕ = eigs(B, nev=1, which=:SR);

julia> λ
1-element Array{Complex{Float64},1}:
 -8.881784197001252e-16 + 3.999999999999997im

julia> λ, ϕ = eigs(B, nev=1, sigma=1.5);

julia> λ
1-element Array{Complex{Float64},1}:
 1.0000000000000004 + 4.0417078924070745e-18im
Note

The sigma and which keywords interact: the description of eigenvalues searched for by which do not necessarily refer to the eigenvalues of A, but rather the linear operator constructed by the specification of the iteration mode implied by sigma.

sigmaiteration modewhich refers to eigenvalues of
nothingordinary (forward)$A$
real or complexinverse with level shift sigma$(A - \\sigma I )^{-1}$
Note

Although tol has a default value, the best choice depends strongly on the matrix A. We recommend that users always specify a value for tol which suits their specific needs.

For details of how the errors in the computed eigenvalues are estimated, see:

  • B. N. Parlett, "The Symmetric Eigenvalue Problem", SIAM: Philadelphia, 2/e (1998), Ch. 13.2, "Accessing Accuracy in Lanczos Problems", pp. 290-292 ff.
  • R. B. Lehoucq and D. C. Sorensen, "Deflation Techniques for an Implicitly Restarted Arnoldi Iteration", SIAM Journal on Matrix Analysis and Applications (1996), 17(4), 789–821. doi:10.1137/S0895479895281484

For the two-input generalized eigensolution version,

eigs(A, B; nev=6, ncv=max(20,2*nev+1), which=:LM, tol=0.0, maxiter=300, sigma=nothing, ritzvec=true, v0=zeros((0,))) -> (d,[v,],nconv,niter,nmult,resid)

the following keyword arguments are supported:

whichtype of eigenvalues
:LMeigenvalues of largest magnitude (default)
:SMeigenvalues of smallest magnitude
:LReigenvalues of largest real part
:SReigenvalues of smallest real part
:LIeigenvalues of largest imaginary part (nonsymmetric or complex A only)
:SIeigenvalues of smallest imaginary part (nonsymmetric or complex A only)
:BEcompute half of the eigenvalues from each end of the spectrum, biased in favor of the high end. (real symmetric A only)

eigs returns the nev requested eigenvalues in d, the corresponding Ritz vectors v (only if ritzvec=true), the number of converged eigenvalues nconv, the number of iterations niter and the number of matrix vector multiplications nmult, as well as the final residual vector resid.

We can see the various keywords in action in the following examples:

julia> A = sparse(1.0I, 4, 4); B = Diagonal(1:4);

julia> λ, ϕ = eigs(A, B, nev = 2);

julia> λ
2-element Array{Float64,1}:
 1.0000000000000002
 0.5

julia> A = Diagonal([1, -2im, 3, 4im]); B = sparse(1.0I, 4, 4);

julia> λ, ϕ = eigs(A, B, nev=1, which=:SI);

julia> λ
1-element Array{Complex{Float64},1}:
 -1.5720931501039814e-16 - 1.9999999999999984im

julia> λ, ϕ = eigs(A, B, nev=1, which=:LI);

julia> λ
1-element Array{Complex{Float64},1}:
 0.0 + 4.000000000000002im
Note

The sigma and which keywords interact: the description of eigenvalues searched for by which do not necessarily refer to the eigenvalue problem $Av = Bv\\lambda$, but rather the linear operator constructed by the specification of the iteration mode implied by sigma.

sigmaiteration modewhich refers to the problem
nothingordinary (forward)$Av = Bv\\lambda$
real or complexinverse with level shift sigma$(A - \\sigma B )^{-1}B = v\\nu$
Arpack.eigsMethod.
eigs(A; nev=6, ncv=max(20,2*nev+1), which=:LM, tol=0.0, maxiter=300, sigma=nothing, ritzvec=true, v0=zeros((0,))) -> (d,[v,],nconv,niter,nmult,resid)

Computes eigenvalues d of A using implicitly restarted Lanczos or Arnoldi iterations for real symmetric or general nonsymmetric matrices respectively. See the manual for more information.

eigs returns the nev requested eigenvalues in d, the corresponding Ritz vectors v (only if ritzvec=true), the number of converged eigenvalues nconv, the number of iterations niter and the number of matrix vector multiplications nmult, as well as the final residual vector resid.

Examples

julia> using Arpack

julia> A = Diagonal(1:4);

julia> λ, ϕ = eigs(A, nev = 2);

julia> λ
2-element Array{Float64,1}:
 4.0
 3.0
source
Arpack.eigsMethod.
eigs(A, B; nev=6, ncv=max(20,2*nev+1), which=:LM, tol=0.0, maxiter=300, sigma=nothing, ritzvec=true, v0=zeros((0,))) -> (d,[v,],nconv,niter,nmult,resid)

Computes generalized eigenvalues d of A and B using implicitly restarted Lanczos or Arnoldi iterations for real symmetric or general nonsymmetric matrices respectively. See the manual for more information.

source
Arpack.svdsFunction.
svds(A; nsv=6, ritzvec=true, tol=0.0, maxiter=1000, ncv=2*nsv, v0=zeros((0,))) -> (SVD([left_sv,] s, [right_sv,]), nconv, niter, nmult, resid)

Computes the largest singular values s of A using implicitly restarted Lanczos iterations derived from eigs.

Inputs

  • A: Linear operator whose singular values are desired. A may be represented as a subtype of AbstractArray, e.g., a sparse matrix, or any other type supporting the four methods size(A), eltype(A), A * vector, and A' * vector.
  • nsv: Number of singular values. Default: 6.
  • ritzvec: If true, return the left and right singular vectors left_sv and right_sv. If false, omit the singular vectors. Default: true.
  • tol: tolerance, see eigs.
  • maxiter: Maximum number of iterations, see eigs. Default: 1000.
  • ncv: Maximum size of the Krylov subspace, see eigs (there called nev). Default: 2*nsv.
  • v0: Initial guess for the first Krylov vector. It may have length min(size(A)...), or 0.

Outputs

  • svd: An SVD object containing the left singular vectors, the requested values, and the right singular vectors. If ritzvec = false, the left and right singular vectors will be empty.
  • nconv: Number of converged singular values.
  • niter: Number of iterations.
  • nmult: Number of matrix–vector products used.
  • resid: Final residual vector.

Examples

julia> A = Diagonal(1:4);

julia> s = svds(A, nsv = 2)[1];

julia> s.S
2-element Array{Float64,1}:
 4.0
 2.9999999999999996
Implementation

svds(A) is formally equivalent to calling eigs to perform implicitly restarted Lanczos tridiagonalization on the Hermitian matrix $A^\prime A$ or $AA^\prime$ such that the size is smallest.

source