The Python interface extends the Python CVXOPT package with functions from the CHOMPACK library. It is based on the CVXOPT dense and sparse matrix types, and two new Python objects: a chordal matrix object for storing symmetric chordal matrices, and a chordal factor object for storing Cholesky factors of positive definite chordal matrices. These two data types can be created and manipulated, and converted to CVXOPT matrices. using the functions described below. They can be thought of as opaque objects that contain the values of the matrix, plus the sparsity pattern and a perfect elimination ordering for the sparsity pattern.
When the following documentation states a requirement that two chordal matrix objects and/orfactors have the same sparsity pattern, we mean by this that they were created from the same chordal matrix object via a series of calls to functions that create a new chordal matrix or factor from an existing one (such as cholesky, completion, copy, llt, partial_inv, or project). Chordal matrices that were created from CVXOPT matrices via the function embed are not recognized to have the same sparsity pattern, even though their sparsity patterns may be equal (mathematically).
CVXOPT integer matrices are used to represent permutation matrices. The
relation between the CVXOPT integer matrix p (a permutation of the
column matrix with entries 0, 1, ..., n-1) and the permutation matrix
it represents is as follows: if the CVXOPT matrix X has
value
, then the CVXOPT matrix X[p, p] has value
.
The functions in the package can be divided in different groups.
Computes a chordal embedding of a sparse matrix.
Y, nfill = embed(X, p = None)
Returns a chordal embedding of the sparsity pattern of X, projects X on the embedding, and returns the result as a chordal matrix object. The argument p is a permutation with as default value the natural ordering matrix([0, 1, ..., n-1]). The embedding is computed via a symbolic Cholesky factorization of X[p, p].
Parameters: |
|
---|---|
Returns: | tuple (Y, nfill) with Y a chordal matrix embedding and nfill the number of nonzero entries added in the embedding |
Projects a CVXOPT sparse matrix on a chordal sparsity pattern.
C = project(X, Y)
Projects the CVXOPT sparse matrix Y on the sparsity pattern of the chordal matrix X, and returns the result as a chordal matrix. Only the lower triangular part of X is referenced.
Parameters: |
|
---|---|
Returns: | chordal matrix with the same sparsity pattern as X |
Converts a chordal matrix or factor to a CVXOPT sparse matrix.
L = sparse(X)
If X is a chordal matrix, the function returns the lower triangular part of X[p, p] as a CVXOPT sparse matrix. If X is a chordal factor for a Cholesky factorization (1) the function returns the lower triangular sparse matrix L.
Parameters: | X – chordal matrix or factor |
---|---|
Returns: | lower triangular CVXOPT sparse square matrix of doubles |
Cholesky factorization.
L = cholesky(X)
Computes a zero fill-in Cholesky factorization
(1)
of a positive definite chordal matrix X.
Parameters: | X – chordal matrix |
---|---|
Returns: | chordal factor if the factorization was successful |
Raises : | an ArithmethicError if the matrix is not positive definite |
Solves a factored set of equations, or multiplies with Cholesky factors
solve(L, X, mode = 0)
L contains the factors of a factorization (1) of a positive definite sparse chordal matrix. X is a CVXOPT dense matrix of doubles with the same number of rows as L. On exit, X is overwritten with one of the four matrices in the table.
action | mode |
---|---|
![]() |
0 |
![]() |
1 |
![]() |
2 |
![]() |
3 |
Parameters: |
|
---|
Maximum-determinant positive definite completion.
L = completion(X)
Returns the Cholesky factor of the inverse of the maximum-determinant positive definite completion of a symmetric chordal matrix X, ie, the Cholesky factor of the inverse of the solution of
The inverse has the same sparsity pattern as
and satisfies the nonlinear equation
completion returns the factor in the factorization
.
Parameters: | X – chordal matrix |
---|---|
Returns: | chordal factor with the same sparsity pattern as X |
Raises : | an ArithmethicError if the matrix does not have a positive definite completion |
Evaluates the projection of the inverse of the matrix
on the sparsity pattern of
.
Y = partial_inv(L)
Computes
where is a positive definite chordal matrix specified by
its Cholesky factor
.
Parameters: | L – chordal factor |
---|---|
Returns: | chordal matrix with the same sparsity pattern as L |
The mapping
is the Hessian of the log-det barrier at a positive definite chordal
matrix , applied to a symmetric chordal matrix
.
The Hessian operator can be factored as
where the mappings on the right hand side are adjoint mappings that map chordal symmetric matrices to chordal symmetric matrices.
hessian(L, Y, U, adj = False, inv = False)
evaluates these mappings or their inverses for a list of symmetric chordal matrices U = [ U[0], ..., U[N-1] ], and overwrites the matrices with the results. The following table lists the possible actions.
Action | adj | inv |
---|---|---|
![]() |
False | False |
![]() |
False | True |
![]() |
True | False |
![]() |
True | True |
The input argument L is the Cholesky factor of ,
as computed by the command L = cholesky(X).
The input argument Y is the partial inverse of the inverse of
, as computed by the command Y = partial_inv(L).
The input argument U is a list of chordal matrices with the
same sparsity pattern as L and Y.
The matrices can be computed by two calls
hessian(L, Y, U, adj = False, inv = False)
hessian(L, Y, U, adj = True, inv = False)
The matrices can be computed as
hessian(L, Y, U, adj = True, inv = True)
hessian(L, Y, U, adj = False, inv = True)
Parameters: |
|
---|
Returns a copy of a chordal matrix.
Parameters: | X – chordal matrix |
---|---|
Returns: | chordal matrix with the same sparsity pattern and numerical values as X |
Evaluates
Parameters: |
|
---|
Evaluates
Parameters: |
|
---|
Inner product of symmetric chordal sparse matrices.
Returns the inner product
of two symmetric sparse matrices with the same chordal sparsity pattern.
Parameters: |
|
---|---|
Return type: | float |
Computes projected rank 2 update of a chordal matrix
where is of order n, and
and
are n-vectors.
Parameters: |
|
---|---|
Return type: | float |
Computes a symmetric matrix from its Cholesky factorization
X = llt(L)
Computes X from its Cholesky factorization .
Parameters: | L – chordal factor |
---|---|
Returns: | chordal matrix with the same sparsity pattern as L |
Returns the logarithm of the determinant of a Cholesky factor L.
Parameters: | L – chordal factor |
---|---|
Returns: | float |
Returns a dictionary with information about a chordal sparsity pattern.
Parameters: | X – chordal matrix or factor |
---|---|
Returns: | python dictionary representation of the sparsity pattern |
Maximum cardinality search ordering of a sparse chordal matrix.
Returns the maximum cardinality search ordering of a symmetric chordal matrix X. The maximum cardinality search ordering is a perfect elimination ordering for the Cholesky factorization.
Parameters: | X – CVXOPT sparse square matrix of doubles. Only the sparsity pattern of the lower triangular part of the matrix is accessed |
---|---|
Returns: | CVXOPT dense integer matrix of length n, if n is the order of X |
Checks whether an ordering is a perfect elmimination order.
Returns True if the permutation p is a perfect elimination order for a Cholesky factorization of X.
Parameters: |
|
---|
Performs a symmetric permutation of a square sparse matrix.
Y = perm(X, p)
This is equivalent to but more efficient than Y = X[p, p].
Parameters: |
|
---|---|
Returns: | CVXOPT sparse square matrix of doubles |
Symmetrizes a lower triangular matrix. Returns
where is a lower triangular matrix.
Parameters: | X – CVXOPT sparse square matrix of doubles. Must be lower triangular |
---|---|
Returns: | CVXOPT sparse square matrix of doubles |
Returns the lower triangular part of a sparse matrix X.
Parameters: | X – CVXOPT sparse square matrix of doubles |
---|---|
Returns: | CVXOPT sparse square matrix of doubles |
Computes maximal chordal subgraph of sparsity graph and returns the projection of X on the chordal subgraph as a chordal matrix. The optional parameter k determines the last vertex in a perfect elimination ordering of the maximal chordal subgraph. A node of maximum degree is chosen if k is not specified.
param X: CVXOPT sparse square matrix of doubles. Only the sparsity pattern of the lower triangular part of the matrix is accessed. param k: integer between and
if
is the order of X
returns: chordal matrix