linalg (Linear Algebra Routines)

Module that defines solvers for matrices that are frequently encountered in RBF applications. Most functions in this module can take either scipy sparse matrices or numpy arrays as input

class rbf.linalg.GMRESSolver(A, drop_tol=0.005, fill_factor=2.0, normalize_inplace=False)

Solves the system of equations Ax = b for x iteratively with GMRES and an incomplete LU decomposition.

Parameters:

A ((n, n) CSC sparse matrix) –
drop_tol (float, optional) – Passed to scipy.sparse.linalg.spilu. This controls the sparsity of the ILU decomposition used for the preconditioner. It should be between 0 and 1. smaller values make the decomposition denser but better approximates the LU decomposition. If the value is too large then you may get a “Factor is exactly singular” error.
fill_factor (float, optional) – Passed to scipy.sparse.linalg.spilu. I believe this controls the memory allocated for the ILU decomposition. If this value is too small then memory will be allocated dynamically for the decomposition. If this is too large then you may get a memory error.
normalize_inplace (bool) – If True and A is a csc matrix, then A is normalized in place.

solve(b, tol=1e-10)

Solve Ax = b for x

Parameters:

b ((n,) array) –
tol (float, optional) –

Returns:

(n,) array

class rbf.linalg.PartitionedPosDefSolver(A, B, build_inverse=False, factor_inplace=False)

Solves the system of equations

\[\begin{split}\left[ \begin{array}{cc} A & B \\ B^T & 0 \\ \end{array} \right] \left[ \begin{array}{c} x \\ y \\ \end{array} \right] = \left[ \begin{array}{c} a \\ b \\ \end{array} \right]\end{split}\]

for x and y, where A is a positive definite matrix. Rather than naively building and solving the system, this class partitions the inverse as

\[\begin{split}\left[ \begin{array}{cc} C & D \\ D^T & E \\ \end{array} \right]\end{split}\]

where

\[C = A^{-1} - (A^{-1} B) (B^T A^{-1} B)^{-1} (A^{-1} B)^T\]

\[D = (A^{-1} B) (B^T A^{-1} B)^{-1}\]

\[E = - (B^T A^{-1} B)^{-1}\]

The inverse of A is not computed, but instead its action is performed by solving the Cholesky decomposition of A. A can be a scipy sparse matrix or a numpy array. B can also be either a scipy sparse matrix or a numpy array but it will be converted to a numpy array.

Parameters:

A ((n, n) array or sparse matrix) –
B ((n, p) array or sparse matrix) –
build_inverse (bool, optional) – If True, the inverse of the partitioned matrix is built, which makes subsequent calls to solve faster.
factor_inplace (bool, optional) – If True and if A is fortran contiguous, then the Cholesky factorization of A is done inplace.

Note

This class stores the factorization of A, which may be sparse, the dense matrix A^-1 B, and the dense factorization of B^T A^-1 B. If the number of columns in B is large then this may take up too much memory.

solve(a, b=None)

Solves for x and y given a and b.

Parameters:

a ((n, ...) array or sparse matrix) –
b ((p, ...) array or sparse matrix) –

Returns:

(n, …) array
(p, …) array

class rbf.linalg.PartitionedSolver(A, B, build_inverse=False, check_cond=False)

Solves the system of equations

\[\begin{split}\left[ \begin{array}{cc} A & B \\ B^T & 0 \\ \end{array} \right] \left[ \begin{array}{c} x \\ y \\ \end{array} \right] = \left[ \begin{array}{c} a \\ b \\ \end{array} \right]\end{split}\]

for x and y. This class builds the system and then factors it with an LU decomposition. As opposed to PartitionedPosDefSolver, A is not assumed to be positive definite. A can be a sparse matrix or an array. B can also be a sparse matrix but it will be converted to an array.

Parameters:

A ((n, n) array or sparse matrix) –
B ((n, p) array or sparse matrix) –
build_inverse (bool, optional) – If True, the inverse of the partitioned matrix is built, which makes subsequent calls to solve faster.
check_cond (bool, optional) – If True, a warning is raised if the partitioned matrix is ill-conditioned. Ignored if A is sparse.

solve(a, b=None)

Solves for x and y given a and b.

Parameters:

a ((n, ...) array or sparse matrix) –
b ((p, ...) array or sparse matrix, optional) – If not given, it is assumed to be zeros

Returns:

(n, …) array
(p, …) array

class rbf.linalg.PosDefSolver(A, build_inverse=False, factor_inplace=False)

Sparse or dense positive definite matrix solver.

Factors the positive definite matrix A as LL^T = A and provides an efficient method for solving Ax = b for x. Additionally provides a method to solve Lx = b, get the log determinant of A, and get L.

Parameters:

A ((n, n) array or sparse matrix) – Positive definite matrix
build_inverse (bool, optional) – If True, the inverse of A is built, which makes subsequent calls to solve faster.
factor_inplace (bool, optional) – If True and if A is fortran contiguous, then the Cholesky factorization of A is done inplace.

L()

Returns the factorization L.

Returns:: (n, n) array or sparse matrix

log_det()

Returns the log determinant of A.

Returns:: float

solve(b)

solves Ax = b for x.

Parameters:: b ((n, ...) array or sparse matrix) –
Returns:: (n, …) array

solve_L(b)

solves Lx = b for x.

Parameters:: b ((n, ...) array or sparse matrix) –
Returns:: (n, …) array

class rbf.linalg.Solver(A, build_inverse=False, check_cond=False, factor_inplace=False)

Sparse or dense matrix solver.

Parameters:

A ((n, n) array or sparse matrix) –
build_inverse (bool, optional) – If True, the inverse of A is built, which makes subsequent calls to solve faster.
check_cond (bool, optional) – If True, a warning is raised if A is ill-conditioned. Ignored if A is sparse.
factor_inplace (bool, optional) – If True and if A is fortran contiguous, then the LU factorization of A is done inplace.

solve(b)

solves Ax = b for x.

Parameters:: b ((n, ...) array or sparse matrix) –
Returns:: (n, …) array

rbf.linalg.as_array(A, dtype=None, copy=False): Return A as an array if it is not already. This properly handles when A is sparse.

rbf.linalg.as_sparse_or_array(A, dtype=None, copy=False): If A is a scipy sparse matrix then return it as a csc matrix. Otherwise, return it as an array.

rbf.linalg.is_positive_definite(A): Tests if A is positive definite. This is done by testing whether the Cholesky decomposition finishes successfully. A can be a sparse matrix or array.