Usage
To use the plotting capabilities, one should first explicitly import either Plots
or PyPlot
. (An interface to Makie
is in progress.)
A simple interface for modestly-sized matrices is available through the spectralportrait
function.
Typical "advanced" use of the package is as follows:
using Pseudospectra
using Plots
A = your_matrix_generating_function()
ps_data = new_matrix(A)
driver!(ps_data)
options = Dict{Symbol,Any}()
# modify `options` to concentrate on a region of interest, increase resolution, etc.
driver!(ps_data,options)
This should show contour plots of $\log_{10}(\epsilon)$ in the vicinity of the spectrum - the standard display of a spectral portrait. Eigenvalues of A
are displayed as black points, if they are available.
The new_matrix
function constructs a stateful object with information about A
conducive to pseudospectral analysis; driver!
then manages the appropriate computations and plotting.
Requirements
The integrated plotting capabilities require that the Plots
and/or PyPlot
packages be installed. These are not formal package requirements because much of the Pseudospectra
package is useful without plotting.
Use without the integrated plotters
One can use the psa_compute
function to "simply" evaluate resolvent norms on a grid, as follows:
using LinearAlgebra, Pseudospectra
# matrix must be complex to get a true Schur decomposition
A = your_matrix_generator() .+ 0im
Tschur, U, eigA = schur(A)
ax = [-3,3,-3,3] # domain of results
npts = 100 # mesh size
opts = Dict{Symbox,Any}()
Z, x, y, levels = psa_compute(Tschur, npts, ax, eigA, opts)
Use with extended precision
The dense-matrix methods work with matrix element types of extended precision such as Float128
(from the Quadmath
package) and BigFloat
; they do require linear algebra methods which are implemented in the GenericLinearAlgebra
and GenericSchur
packages.
using Pseudospectra, Plots, GenericLinearAlgebra, GenericSchur, Quadmath
Aq = Float128.(A_matrix)
spectralportrait(A)
Note that linear algebra with non-BLAS types is expensive, so computation takes many minutes even for moderate-sized matrices.
Parallelism
The computation of spectral portraits lends itself to parallel processing. Multiprocessing is implemented here for some algorithm variants.
The core dense-matrix algorithms can be set to distribute the Z
-grid over Julia threads with the threaded
entry in options
. This is especially useful for extended precision, and for working on fine grids with BLAS types. (One may need to adjust thread usage in the BLAS to balance if the library is aggressive).
Unfortunately BigFloat
operations currently make heavy use of the heap, invoking the thread-unfriendly garbage collector, so for extended precision Float128
may be more practical.
Large matrices
Currently, the methods for handling large matrices are those ported from Eigtool. They can be invoked via options to the user-facing computational functions.
Direct-sparse
This is a Lanczos scheme which does a sparse LU decomposition for each point in the Z
grid.
Regional subspace projection
This uses one Schur factorization of the matrix and is not always reliable. [More documentation needed]
Arnoldi projection
This uses a specialized wrapper of the ARPACK library, and may also be used for linear maps not explicitly represented as matrices. [More documentation needed]