+++ /dev/null
-"""
-Utility functions for working with CVXOPT matrices (instances of the
-class:`cvxopt.base.matrix` class).
-"""
-
-from math import sqrt
-from cvxopt import matrix
-from cvxopt.lapack import gees, syevr
-
-from . import options
-
-
-def append_col(left, right):
- """
- Append two matrices side-by-side.
-
- Parameters
- ----------
-
- left, right : matrix
- The two matrices to append to one another.
-
- Returns
- -------
-
- matrix
- A new matrix consisting of ``right`` appended to the right
- of ``left``.
-
- Examples
- --------
-
- >>> A = matrix([1,2,3,4], (2,2))
- >>> B = matrix([5,6,7,8,9,10], (2,3))
- >>> print(A)
- [ 1 3]
- [ 2 4]
- <BLANKLINE>
- >>> print(B)
- [ 5 7 9]
- [ 6 8 10]
- <BLANKLINE>
- >>> print(append_col(A,B))
- [ 1 3 5 7 9]
- [ 2 4 6 8 10]
- <BLANKLINE>
-
- """
- return matrix([left.trans(), right.trans()]).trans()
-
-
-def append_row(top, bottom):
- """
- Append two matrices top-to-bottom.
-
- Parameters
- ----------
-
- top, bottom : matrix
- The two matrices to append to one another.
-
- Returns
- -------
-
- matrix
- A new matrix consisting of ``bottom`` appended below ``top``.
-
- Examples
- --------
-
- >>> A = matrix([1,2,3,4], (2,2))
- >>> B = matrix([5,6,7,8,9,10], (3,2))
- >>> print(A)
- [ 1 3]
- [ 2 4]
- <BLANKLINE>
- >>> print(B)
- [ 5 8]
- [ 6 9]
- [ 7 10]
- <BLANKLINE>
- >>> print(append_row(A,B))
- [ 1 3]
- [ 2 4]
- [ 5 8]
- [ 6 9]
- [ 7 10]
- <BLANKLINE>
-
- """
- return matrix([top, bottom])
-
-
-def eigenvalues(symmat):
- """
- Return the eigenvalues of the given symmetric real matrix.
-
- On the surface, this appears redundant to the :func:`eigenvalues_re`
- function. However, if we know in advance that our input is
- symmetric, a better algorithm can be used.
-
- Parameters
- ----------
-
- symmat : matrix
- The real symmetric matrix whose eigenvalues you want.
-
- Returns
- -------
-
- list of float
- A list of the eigenvalues (in no particular order) of ``symmat``.
-
- Raises
- ------
-
- TypeError
- If the input matrix is not symmetric.
-
- Examples
- --------
-
- >>> A = matrix([[2,1],[1,2]], tc='d')
- >>> eigenvalues(A)
- [1.0, 3.0]
-
- If the input matrix is not symmetric, it may not have real
- eigenvalues, and we don't know what to do::
-
- >>> A = matrix([[1,2],[3,4]])
- >>> eigenvalues(A)
- Traceback (most recent call last):
- ...
- TypeError: input must be a symmetric real matrix
-
- """
- if not norm(symmat.trans() - symmat) < options.ABS_TOL:
- # Ensure that ``symmat`` is symmetric (and thus square).
- raise TypeError('input must be a symmetric real matrix')
-
- domain_dim = symmat.size[0]
- eigs = matrix(0, (domain_dim, 1), tc='d')
- syevr(symmat, eigs)
- return list(eigs)
-
-
-def eigenvalues_re(anymat):
- """
- Return the real parts of the eigenvalues of the given square matrix.
-
- Parameters
- ----------
-
- anymat : matrix
- The square matrix whose eigenvalues you want.
-
- Returns
- -------
-
- list of float
- A list of the real parts (in no particular order) of the
- eigenvalues of ``anymat``.
-
- Raises
- ------
-
- TypeError
- If the input matrix is not square.
-
- Examples
- --------
-
- This is symmetric and has two real eigenvalues:
-
- >>> A = matrix([[2,1],[1,2]], tc='d')
- >>> sorted(eigenvalues_re(A))
- [1.0, 3.0]
-
- But this rotation matrix has eigenvalues `i` and `-i`, both of whose
- real parts are zero:
-
- >>> A = matrix([[0,-1],[1,0]])
- >>> eigenvalues_re(A)
- [0.0, 0.0]
-
- If the input matrix is not square, it doesn't have eigenvalues::
-
- >>> A = matrix([[1,2],[3,4],[5,6]])
- >>> eigenvalues_re(A)
- Traceback (most recent call last):
- ...
- TypeError: input matrix must be square
-
- """
- if not anymat.size[0] == anymat.size[1]:
- raise TypeError('input matrix must be square')
-
- domain_dim = anymat.size[0]
- eigs = matrix(0, (domain_dim, 1), tc='z')
-
- # Create a copy of ``anymat`` here for two reasons:
- #
- # 1. ``gees`` clobbers its input.
- # 2. We need to ensure that the type code of ``dummy`` is 'd' or 'z'.
- #
- dummy = matrix(anymat, anymat.size, tc='d')
-
- gees(dummy, eigs)
- return [eig.real for eig in eigs]
-
-
-def identity(domain_dim):
- """
- Create an identity matrix of the given dimensions.
-
- Parameters
- ----------
-
- domain_dim : int
- The dimension of the vector space on which the identity will act.
-
- Returns
- -------
-
- matrix
- A ``domain_dim``-by-``domain_dim`` dense integer identity matrix.
-
- Raises
- ------
-
- ValueError
- If you ask for the identity on zero or fewer dimensions.
-
- Examples
- --------
-
- >>> print(identity(3))
- [ 1 0 0]
- [ 0 1 0]
- [ 0 0 1]
- <BLANKLINE>
-
- """
- if domain_dim <= 0:
- raise ValueError('domain dimension must be positive')
-
- entries = [int(i == j)
- for i in range(domain_dim)
- for j in range(domain_dim)]
- return matrix(entries, (domain_dim, domain_dim))
-
-
-def inner_product(vec1, vec2):
- """
- Compute the Euclidean inner product of two vectors.
-
- Parameters
- ----------
-
- vec1, vec2 : matrix
- The two vectors whose inner product you want.
-
- Returns
- -------
-
- float
- The inner product of ``vec1`` and ``vec2``.
-
- Raises
- ------
-
- TypeError
- If the lengths of ``vec1`` and ``vec2`` differ.
-
- Examples
- --------
-
- >>> x = [1,2,3]
- >>> y = [3,4,1]
- >>> inner_product(x,y)
- 14
-
- >>> x = matrix([1,1,1])
- >>> y = matrix([2,3,4], (1,3))
- >>> inner_product(x,y)
- 9
-
- >>> x = [1,2,3]
- >>> y = [1,1]
- >>> inner_product(x,y)
- Traceback (most recent call last):
- ...
- TypeError: the lengths of vec1 and vec2 must match
-
- """
- if not len(vec1) == len(vec2):
- raise TypeError('the lengths of vec1 and vec2 must match')
-
- return sum([x*y for (x, y) in zip(vec1, vec2)])
-
-
-def norm(matrix_or_vector):
- """
- Return the Frobenius norm of a matrix or vector.
-
- When the input is a vector, its matrix-Frobenius norm is the same
- thing as its vector-Euclidean norm.
-
- Parameters
- ----------
-
- matrix_or_vector : matrix
- The matrix or vector whose norm you want.
-
- Returns
- -------
-
- float
- The norm of ``matrix_or_vector``.
-
- Examples
- --------
-
- >>> v = matrix([1,1])
- >>> print('{:.5f}'.format(norm(v)))
- 1.41421
-
- >>> A = matrix([1,1,1,1], (2,2))
- >>> norm(A)
- 2.0
-
- """
- return sqrt(inner_product(matrix_or_vector, matrix_or_vector))
-
-
-def vec(mat):
- """
- Create a long vector in column-major order from ``mat``.
-
- Parameters
- ----------
-
- mat : matrix
- Any sort of real matrix that you want written as a long vector.
-
- Returns
- -------
-
- matrix
- An ``len(mat)``-by-``1`` long column vector containign the
- entries of ``mat`` in column major order.
-
- Examples
- --------
-
- >>> A = matrix([[1,2],[3,4]])
- >>> print(A)
- [ 1 3]
- [ 2 4]
- <BLANKLINE>
-
- >>> print(vec(A))
- [ 1]
- [ 2]
- [ 3]
- [ 4]
- <BLANKLINE>
-
- Note that if ``mat`` is a vector, this function is a no-op:
-
- >>> v = matrix([1,2,3,4], (4,1))
- >>> print(v)
- [ 1]
- [ 2]
- [ 3]
- [ 4]
- <BLANKLINE>
- >>> print(vec(v))
- [ 1]
- [ 2]
- [ 3]
- [ 4]
- <BLANKLINE>
-
- """
- return matrix(mat, (len(mat), 1))
projects / dunshire.git / blobdiff