FastCholesky

This package exports fastcholesky function, which works exactly like the cholesky from the Julia, but faster!

FastCholesky.cholinvMethod
cholinv(input)

Calculate the inverse of the input matrix input using Cholesky factorization. This function is an alias for inv(fastcholesky(input)).

julia> A = [4.0 2.0; 2.0 5.0];

julia> A_inv = cholinv(A);

julia> A_inv ≈ inv(A)
true
source
FastCholesky.cholinv_logdetMethod
cholinv_logdet(input)

Calculate the inverse and the natural logarithm of the determinant of the input matrix input simultaneously using Cholesky factorization.

julia> A = [4.0 2.0; 2.0 5.0];

julia> A_inv, logdet_A = cholinv_logdet(A);

julia> isapprox(A_inv * A, I)
true

julia> isapprox(logdet_A, log(det(A)))
true
source
FastCholesky.chollogdetMethod
chollogdet(input)

Calculate the log-determinant of the input matrix input using Cholesky factorization. This function is an alias for logdet(fastcholesky(input)).

julia> A = [4.0 2.0; 2.0 5.0];

julia> logdet_A = chollogdet(A);

julia> isapprox(logdet_A, log(det(A)))
true
source
FastCholesky.cholsqrtMethod
cholsqrt(input)

Calculate the Cholesky square root of the input matrix input. This function is an alias for fastcholesky(input).L. NOTE: This is not equal to the standard matrix square root used in literature, which requires the result to be symmetric.

julia> A = [4.0 2.0; 2.0 5.0];

julia> A_sqrt = cholsqrt(A);

julia> isapprox(A_sqrt * A_sqrt', A)
true
source
FastCholesky.fastcholesky!Function
fastcholesky!(input; fallback_gmw81=true, symmetrize_input=true, gmw81_tol=PositiveFactorizations.default_δ(input), symmetric_tol=1e-8)

Calculate the Cholesky factorization of the input matrix input in-place. This function is an in-place version of fastcholesky. It first checks if the input matrix is symmetric, and if it is, it will use the built-in Cholesky factorization by wrapping the input in a Hermitian matrix. If the input matrix is not symmetric and symmetrize_input=true, it will symmetrize the input matrix and try again recursively. If the input matrix is not symmetrics, not positive definite, and fallback_gmw81=true, it will use the GMW81 algorithm as a fallback. In other cases, it will throw an error.

Keyword arguments

  • fallback_gmw81::Bool=true: If true, the function will use the GMW81 algorithm as a fallback if the input matrix is not positive definite.
  • symmetrize_input::Bool=true: If true, the function will symmetrize the input matrix before retrying the Cholesky factorization in case if the first attempt failed.
  • gmw81_tol::Real=PositiveFactorizations.default_δ(input): The tolerance for the positive-definiteness of the input matrix for the GMW81 algorithm.
  • symmetric_tol::Real=1e-8: The tolerance for the symmetry of the input matrix.

Environment Variables

The behavior of this function regarding non-symmetric matrices can be controlled through environment variables:

  • JULIA_FASTCHOLESKY_NO_WARN_NON_SYMMETRIC=1: Set this to suppress warnings about non-symmetric input matrices
  • JULIA_FASTCHOLESKY_THROW_ERROR_NON_SYMMETRIC=1: Set this to make the function error instead of warn when encountering non-symmetric matrices
Note

Use this function only when you expect the input matrix to be nearly positive definite.

julia> C = fastcholesky!([ 1.0 0.5; 0.5 1.0 ]);

julia> C.L * C.L' ≈ [ 1.0 0.5; 0.5 1.0 ]
true
source
FastCholesky.fastcholeskyMethod
fastcholesky(input)

Calculate the Cholesky factorization of the input matrix input. This function provides a more efficient implementation for certain input matrices compared to the standard LinearAlgebra.cholesky function. By default, it falls back to using LinearAlgebra.cholesky(PositiveFactorizations.Positive, input), which means it does not require the input matrix to be perfectly symmetric.

Note

This function assumes that the input matrix is nearly positive definite, and it will attempt to make the smallest possible adjustments to the matrix to ensure it becomes positive definite. Note that the magnitude of these adjustments may not necessarily be small, so it's important to use this function only when you expect the input matrix to be nearly positive definite.

Environment Variables

The behavior of this function regarding non-symmetric matrices can be controlled through environment variables:

  • JULIA_FASTCHOLESKY_NO_WARN_NON_SYMMETRIC=1: Set this to suppress warnings about non-symmetric input matrices
  • JULIA_FASTCHOLESKY_THROW_ERROR_NON_SYMMETRIC=1: Set this to make the function error instead of warn when encountering non-symmetric matrices
julia> C = fastcholesky([ 1.0 0.5; 0.5 1.0 ]);

julia> C.L * C.L' ≈ [ 1.0 0.5; 0.5 1.0 ]
true
source