TraceEstimators

TraceEstimators.ChebyHutchSpaceMethod
ChebyHutchSpace(A::AbstractMatrix, a::Number, b::Number; fn::Function=invfun, dfn::Function=rademacherDistribution!, m = 4, n = 6)

Create Chebyshev-Hutchinson Workspace for chebyhutch or chebydiagonal methods.

Arguments

  • A : Symmetric Positive Definite Matrix
  • a : Bound for minimum eigenvalue of A
  • b : Bound for maximum eigenvalue of A
  • fn : Function to appy. By default uses inverse function.
  • dfn : Distribution function that returns a vector v. By default uses rademacherDistribution!
  • m : Iteration number, increase this for precision. By default m = 4
  • n : Polynomial degree, increase this for accuracy. By default n = 6
TraceEstimators.HutchWorkspaceMethod
HutchWorkspace(A; N = 30, skipverify = false)
HutchWorkspace(A, randfunc::Function; N = 30, skipverify = false)

Arguments

  • A : Symmetric Positive Definite Matrix with Low Condtion Number (k < 500)
  • randfunc : Function to generate random values for x distributed uniformly (Base: rand(-1.0:2.0:1.0, size(A,1)) Example: f(n) = rand(-1.0:2.0:1.0, n)
  • N : Number of iterations (Default: 30)
  • skipverify : If false, it will check isposdef(A) (Default: false)
TraceEstimators.SLQWorkspaceMethod
SLQWorkspace(A::AbstractMatrix; fn::Function, dfn::Function, ctol, m, nv)

Create an SLQWorkspace for supplied SPD Matrix A. Use it to calculate tr(fn(A)).

Arguments

  • A : Symmetric Positive Definite Matrix
  • fn : Function to apply. By default uses inverse function
  • dfn : Distribution function for v (random dist. with norm(v) = 1). By default uses rademacherDistribution!(v::Vector, t::Type)
  • ctol : SLQ Convergence Tolerance value. By default ctol = 0.1
  • m : Specify value for lanczos steps. By default m = 15
  • nv : Specify value for SLQ iterations. By default nb = 10
TraceEstimators.chebydiagonalMethod
chebydiagonal(w::ChebyHutchSpace)
chebydiagonal(A::AbstractMatrix, m::Integer, n::Integer; fn::Function=invfun, dfn::Function=rademacherDistribution!)

Chebyshev-Hutchinson to estimate diagonal elements of fn(A), for given matrix A and an analytic function fn.

Arguments

  • A : Symmetric Positive Definite Matrix
  • m : Iteration number, increase this for precision. By default m = 4
  • n : Polynomial degree, increase this for accuracy. By default n = 6
  • fn : Function to appy. By default uses inverse function.
  • dfn : Distribution function that returns a vector v. By default uses rademacherDistribution!
TraceEstimators.chebyhutchMethod
chebyhutch(w::ChebyHutchSpace)
chebyhutch(A::AbstractMatrix, m::Integer, n::Integer; a::Float64, b::Float64, fn::Function=invfun, dfn::Function=rademacherDistribution!)
chebyhutch(A::AbstractMatrix; fn::Function=invfun, dfn::Function=rademacherDistribution!)

Chebyshev-Hutchinson to estimate tr(fn(A)), for given matrix A and an analytic function fn.

Arguments

  • A : Symmetric Positive Definite Matrix
  • m : Iteration number, increase this for precision. By default m = 4
  • n : Polynomial degree, increase this for accuracy. By default n = 6
  • fn : Function to appy. By default uses inverse function.
  • dfn : Distribution function that returns a vector v. By default uses rademacherDistribution!
TraceEstimators.diag_invMethod
diag_inv(A::AbstractMatrix, n::Int64, p::Int64; pc = "chol", model = "linear")

Diagonal Approximation algorithm for the matrix inverse.

Arguments

  • A : Symmetric Positive Definite Matrix
  • n : Probing vector count for initial approximation.
  • p : Sample Points count for interpolation
  • pc : Preconditioner used for initial approximation and cg ("chol" - Incomplete Cholesky, "ilu" - IncompleteLU, "amg" - AlgebraicMultigrid, "cheby" - Chebyshev Approximation). Default = "chol"
  • model : Fitting model used for calculation of trace ("linear" - Linear Regression, "pchip" - PCHIP interpolation). Default = "linear".
TraceEstimators.diagappMethod
diagapp(A::AbstractMatrix; randfunc=Base.rand)

Diagonal Approximation algorithm for inverse of matrix.

Arguments

  • A : Symmetric Positive Definite Matrix with Low Condtion Number (k < 500)
  • randfunc : Random function for the particular type of matrix A (An example can be found in test/diagapp.jl)
TraceEstimators.hutch!Method
hutch!(w::HutchWorkspace; aitken = false)
hutch!(A::AbstractArray; N = 30, skipverify = false, aitken = false)

Take a HutchWorkspace object as input, apply hutchinson estimation algorithm and solve for trace of inverse of the matrix. (in-place version of hutch)

TraceEstimators.slqMethod
slq(w::SLQWorkspace; skipverify = false)
slq(A::AbstractMatrix; skipverify = false, fn::Function = invfun,
dfn::Function = Base.rand, ctol = 0.1, eps = ϵ, mtol = tol)

SLQ method to calculate tr(fn(A)) for a Symmetric Positive Definite matrix A and an analytic function fn.

Arguments

  • A : Symmetric Positive Definite Matrix
  • skipverify : Skip isposdef(A) verification. By default skipverify = false
  • fn : Function to apply on A before trace calculation. fn must be analytic λₘ and λ₁ of A. By default fn = inv
  • dfn : Distribution function for v (random dist. with norm(v) = 1). By default uses rademacherDistribution!(v::Vector, t::Type)
  • ctol : SLQ Convergence Tolerance value. Decrease this for higher precision. By default ctol = 0.1
  • eps : Error bound for lanczos steps calculation. Decrease this for higher accuracy. By default eps = 0.05
  • mtol : Tolerance for eigenvalue Convergence. Decrease this for precision. By default mtol = 0.01
TraceEstimators.tr_invMethod
tr_inv(A::AbstractMatrix, n::Int64, p::Int64; pc = "chol", model = "linear")

Diagonal Approximation algorithm for calculating trace of the matrix inverse.

Arguments

  • A : Symmetric Positive Definite Matrix
  • n : Probing vector count for initial approximation.
  • p : Sample Points count for interpolation
  • pc : Preconditioner used for initial approximation and cg ("chol" - Incomplete Cholesky, "ilu" - IncompleteLU, "amg" - AlgebraicMultigrid, "cheby" - Chebyshev Approximation). Default = "chol"
  • model : Fitting model used for calculation of trace ("linear" - Linear Regression, "pchip" - PCHIP interpolation). Default = "linear".