cvxopt.cholmod is an interface to the Cholesky factorization
routines of the CHOLMOD package.
It includes functions for Cholesky factorization of sparse positive
definite matrices, and for solving sparse sets of linear equations with
positive definite matrices.
The routines can also be used for computing LDL
(or LDL
) factorizations of symmetric indefinite matrices
(with L unit lower-triangular and D diagonal and nonsingular) if
such a factorization exists.
A, B[, p=None [,
uplo='L']]) |
ArithmeticError
if the factorization does not exist.
As an example, we solve
>>> from cvxopt.base import matrix, spmatrix >>> from cvxopt import cholmod >>> A = spmatrix([10,3, 5,-2, 5, 2], [0,2, 1,3, 2, 3], [0,0, 1,1, 2, 3]) >>> X = matrix(range(8), (4,2), 'd') >>> cholmod.linsolve(A,X) >>> print X -1.4634e-01 4.8780e-02 1.3333e+00 4.0000e+00 4.8780e-01 1.1707e+00 2.8333e+00 7.5000e+00
A, B[, p=None [,
uplo='L']]) |
The following code computes the inverse of the coefficient matrix in (7.2) as a sparse matrix.
>>> X = cholmod.splinsolve(A, spmatrix(1.0,range(4),range(4))) >>> print X SIZE: (4,4) (0, 0) 1.2195e-01 (2, 0) -7.3171e-02 (1, 1) 3.3333e-01 (3, 1) 3.3333e-01 (0, 2) -7.3171e-02 (2, 2) 2.4390e-01 (1, 3) 3.3333e-01 (3, 3) 8.3333e-01
The functions linsolve() and splinsolve() are equivalent to symbolic() and numeric() called in sequence, followed by solve(), respectively, spsolve().
A[, p=None [, uplo='L']]) |
options['supernodal']
.
If uplo is 'L'
, only the lower triangular part of A
is accessed and the upper triangular part is ignored.
If uplo is 'U'
, only the upper triangular part of A
is accessed and the lower triangular part is ignored.
The symbolic factorization is returned as an opaque C object that can be passed to cholmod.numeric().
A, F) |
If F was created by a cholmod.symbolic with
uplo equal to 'L'
, then only the lower triangular part
of A is accessed and the upper triangular part is ignored.
If it was created with uplo is 'U'
, then only the upper
triangular part of A is accessed and the lower triangular part
is ignored.
On successful exit, the factorization is stored in F.
Raises an ArithmeticError
if the factorization does not
exist.
F, B[, sys=0]) |
sys | equation |
0 | ![]() |
1 | ![]() |
2 | ![]() |
3 | ![]() |
4 | ![]() |
5 | ![]() |
6 | ![]() |
7 | ![]() |
8 | ![]() |
(If F is a Cholesky factorization of the form (7.3),
D is an identity matrix in this table.
If A is complex, should be replaced by
.)
The matrix B is a dense 'd'
or 'z'
matrix, with the same type
as A. On exit it is overwritten by the solution.
F, B[, sys=0]) |
For the same example as above:
>>> X = matrix(range(8), (4,2), 'd') >>> F = cholmod.symbolic(A) >>> cholmod.numeric(A,F) >>> cholmod.solve(F,X) >>> print X -1.4634e-01 4.8780e-02 1.3333e+00 4.0000e+00 4.8780e-01 1.1707e+00 2.8333e+00 7.5000e+00
F) |
In the functions listed above, the default values of the control parameters described in the CHOLMOD user guide are used, except for Common->print which is set to 0 instead of 3 and Common->supernodal which is set to 2 instead of 1. These parameters (and a few others) can be modified by making an entry in the dictionary cholmod.options. The meaning of these parameters is as follows.
options['supernodal']
options['print']
As an example that illustrates diag() and the use of cholmod.options, we compute the logarithm of the determinant of the coefficient matrix in (7.2) by two methods.
>>> import math >>> from cvxopt.cholmod import options >>> from cvxopt.base import log >>> options['supernodal'] = 2 >>> F = cholmod.symbolic(A) >>> cholmod.numeric(A,F) >>> print 2.0 * sum(log(cholmod.diag(F))) 5.50533153593
>>> options['supernodal'] = 0 >>> F = cholmod.symbolic(A) >>> cholmod.numeric(A,F) >>> Di = matrix(1.0, (4,1)) >>> cholmod.solve(F,Di,sys=6) >>> print -sum(log(Di)) 5.50533153593