TraceEstimators
TraceEstimators.ChebyHutchSpace
TraceEstimators.HutchWorkspace
TraceEstimators.SLQWorkspace
TraceEstimators.chebydiagonal
TraceEstimators.chebyhutch
TraceEstimators.diag_inv
TraceEstimators.diagapp
TraceEstimators.hutch!
TraceEstimators.slq
TraceEstimators.tr_inv
TraceEstimators.ChebyHutchSpace
— MethodChebyHutchSpace(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 Matrixa
: Bound for minimum eigenvalue of Ab
: Bound for maximum eigenvalue of Afn
: 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 = 4n
: Polynomial degree, increase this for accuracy. By default n = 6
TraceEstimators.HutchWorkspace
— MethodHutchWorkspace(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.SLQWorkspace
— MethodSLQWorkspace(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 Matrixfn
: Function to apply. By default uses inverse functiondfn
: 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.1m
: Specify value for lanczos steps. By default m = 15nv
: Specify value for SLQ iterations. By default nb = 10
TraceEstimators.chebydiagonal
— Methodchebydiagonal(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 Matrixm
: Iteration number, increase this for precision. By default m = 4n
: Polynomial degree, increase this for accuracy. By default n = 6fn
: Function to appy. By default uses inverse function.dfn
: Distribution function that returns a vector v. By default uses rademacherDistribution!
TraceEstimators.chebyhutch
— Methodchebyhutch(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 Matrixm
: Iteration number, increase this for precision. By default m = 4n
: Polynomial degree, increase this for accuracy. By default n = 6fn
: Function to appy. By default uses inverse function.dfn
: Distribution function that returns a vector v. By default uses rademacherDistribution!
TraceEstimators.diag_inv
— Methoddiag_inv(A::AbstractMatrix, n::Int64, p::Int64; pc = "chol", model = "linear")
Diagonal Approximation algorithm for the matrix inverse.
Arguments
A
: Symmetric Positive Definite Matrixn
: Probing vector count for initial approximation.p
: Sample Points count for interpolationpc
: 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.diagapp
— Methoddiagapp(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!
— Methodhutch!(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.slq
— Methodslq(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 Matrixskipverify
: Skip isposdef(A) verification. By default skipverify = falsefn
: Function to apply on A before trace calculation. fn must be analytic λₘ and λ₁ of A. By default fn = invdfn
: 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.1eps
: Error bound for lanczos steps calculation. Decrease this for higher accuracy. By default eps = 0.05mtol
: Tolerance for eigenvalue Convergence. Decrease this for precision. By default mtol = 0.01
TraceEstimators.tr_inv
— Methodtr_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 Matrixn
: Probing vector count for initial approximation.p
: Sample Points count for interpolationpc
: 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".