Module poptools.linalg

Functions

def posdef_linesearch(x: numpy.ndarray | poptools.linalg._blockmat.BlockMatArray,
dx: numpy.ndarray | poptools.linalg._blockmat.BlockMatArray,
precomputed: Literal['nothing', 'chol', 'cholinv'] = 'nothing',
min_step: float = 1e-08,
step_frac: float = 1.0) ‑> float

Given a symmetric positive definite matrix x and a symmetric matrix direction dx, this function computes the maximum step size alpha such that x + alpha * dx remains positive definite. It will need to use the inverse of the lower Cholesky factor of x. Depending on the value of precomputed, it will:

  • 'nothing': compute the Cholesky factor of x internally and invert it,
  • 'chol': assume x is already the lower Cholesky factor of x but still needs to be inverted, or
  • 'cholinv': assume x is already the inverse of the lower Cholesky factor of x.

If the maximum possible step size is less than min_step, it will raise numpy.linalg.LinAlgError. The found step size is optionally multiplied by step_frac (for making sure to stay in the interior of the PSD cone).

def frobenius(a: poptools.linalg._blockmat.BlockMatArray,
b: poptools.linalg._blockmat.BlockMatArray | None = None) ‑> numpy.ndarray

Computes the Frobenius inner product of two block matrix arrays. If b is None, it computes the squared Frobenius norms of each matrix in a and returns them in a 1D array. If b is provided, it computes the Frobenius inner products between all matrices of a and b and returns them in a 2D array.

def cho_factor(x: poptools.linalg._blockmat.BlockMatArray) ‑> poptools.linalg._blockmat.BlockMatArray
def cho_solve(chol: poptools.linalg._blockmat.BlockMatArray,
b: poptools.linalg._blockmat.BlockMatArray) ‑> poptools.linalg._blockmat.BlockMatArray
def maxeigsh(x: poptools.linalg._blockmat.BlockMatArray) ‑> float
def symmetric_part(x: poptools.linalg._blockmat.BlockMatArray) ‑> poptools.linalg._blockmat.BlockMatArray

Returns the average of the block matrix and its transpose.

Classes

class VecSymDomain (n: int)

Handles vectorization of symmetric matrices by storing only the upper triangular part. Thus, it reduces the number of parameters needed to represent a symmetric matrix of size n from n*n to n*(n+1)/2. The entries corresponding to the diagonal are scaled by sqrt(2) so that the Frobenius inner product after vectorization becomes twice the Euclidean inner product.

>>> n = 10
>>> vecsym = VecSymDomain(n)
>>> vecsym.dim
55
>>> a = np.random.normal(size=(n, n))
>>> vec = vecsym.project(a)
>>> vec.shape
(55,)
>>> assert np.isclose(vecsym.frobenius(vec), np.sum((0.5 * (a + a.T))**2))
>>> assert np.allclose(vecsym.unvectorize(vec), 0.5 * (a + a.T))

Methods

def vectorize(self, a: numpy.ndarray) ‑> numpy.ndarray

Vectorizes symmetric matrices stored in the last two dimensions of a into a vector containing only the upper (resp. lower) triangular part. The diagonal entries are scaled by 1/sqrt(2) for faster computation of the Frobenius inner product.

def project(self, a: numpy.ndarray) ‑> numpy.ndarray

Symmetrizes the input matrix (or stack of matrices) a by averaging it with its transpose, and returns it vectorized.

def unvectorize(self, vec: numpy.ndarray) ‑> numpy.ndarray

Exactly reverts the action of vectorize, reconstructing the symmetric matrices from vectors.

def frobenius(self, veca: numpy.ndarray, vecb: numpy.ndarray | None = None) ‑> numpy.ndarray

Computes the Frobenius inner product of two stacks of vectorized symmetric matrices. If vecb is None, it computes the squared Frobenius norm of veca. If vecb is provided, it computes the Frobenius inner product between all vectorized matrices in veca and vecb stored in the respective last dimension.

def matmul_project(self, *vecsymmats_or_mats: numpy.ndarray) ‑> numpy.ndarray

Multiplies the matrices using standard matrix multiplication, and then projects the result to the vectorized symmetric matrix space. The matrices can be either vectorized symmetric matrices or regular matrices, and the decision is made based on their shapes. This is somewhat similar to the Jordan product 0.5 * (a @ b + b @ a) but for multiple matrices.

class BlockMatArray (structure: list[tuple[typing.Literal['dense', 'diagonal'], int]],
blocks: list[numpy.ndarray])

Represents an array of block-diagonal matrices with a fixed structure.

Static methods

def zeros(stack_size: int, s: list[tuple[typing.Literal['dense', 'diagonal'], int]]) ‑> poptools.linalg._blockmat.BlockMatArray

Creates a block diagonal matrix with zero matrices in the blocks.

def identity(stack_size: int, s: list[tuple[typing.Literal['dense', 'diagonal'], int]]) ‑> poptools.linalg._blockmat.BlockMatArray

Creates a block diagonal matrix with identity matrices in the blocks.

def diagonal(stack_size: int,
s: list[tuple[typing.Literal['dense', 'diagonal'], int]],
values: numpy.ndarray) ‑> poptools.linalg._blockmat.BlockMatArray

Creates a block diagonal matrix with given values in the diagonal blocks. The values should be a 1D array with the same length as the number of blocks.

def stack(*bmas: BlockMatArray) ‑> poptools.linalg._blockmat.BlockMatArray

Instance variables

prop shape : tuple[int, int, int]

Returns the shape of the block matrix array. The shape is determined by the sum of the sizes of each block in the structure.

prop TBlockMatArray

Transposes the block matrix array. The structure remains the same, but the blocks are transposed.

Methods

def copy(self) ‑> poptools.linalg._blockmat.BlockMatArray

Returns a copy of the block matrix array.

def dense(self) ‑> numpy.ndarray

Converts the block matrix array to an array of dense matrices.