tensor.slinalg – Linear Algebra Ops Using Scipy

Note

This module is not imported by default. You need to import it to use it.

API

pytensor.tensor.slinalg.block_diag(*matrices)

Construct a block diagonal matrix from a sequence of input tensors.

Given the inputs A, B and C, the output will have these arrays arranged on the diagonal:

[[A, 0, 0],

[0, B, 0], [0, 0, C]]

Parameters:

• A (tensors) – Input tensors to form the block diagonal matrix. last two dimensions of the inputs will be used, and all inputs should have at least 2 dimensins.

• B (tensors) – Input tensors to form the block diagonal matrix. last two dimensions of the inputs will be used, and all inputs should have at least 2 dimensins.

• ... (C) – Input tensors to form the block diagonal matrix. last two dimensions of the inputs will be used, and all inputs should have at least 2 dimensins.

Returns:

out – The block diagonal matrix formed from the input matrices.

Return type:

tensor

Examples

Create a block diagonal matrix from two 2x2 matrices:

..code-block:: python

import numpy as np from pytensor.tensor.linalg import block_diag

A = pt.as_tensor_variable(np.array([[1, 2], [3, 4]])) B = pt.as_tensor_variable(np.array([[5, 6], [7, 8]]))

result = block_diagonal(A, B, name=’X’) print(result.eval()) Out: array([[1, 2, 0, 0],

[3, 4, 0, 0], [0, 0, 5, 6], [0, 0, 7, 8]])

pytensor.tensor.slinalg.cho_solve(c_and_lower, b, *, check_finite=True, b_ndim=None)

Solve the linear equations A x = b, given the Cholesky factorization of A.

Parameters:

• c_and_lower (tuple of (TensorLike, bool)) – Cholesky factorization of a, as given by cho_factor

• b (TensorLike) – Right-hand side

• check_finite (bool, optional) – Whether to check that the input matrices contain only finite numbers. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.

• b_ndim (int) – Whether the core case of b is a vector (1) or matrix (2). This will influence how batched dimensions are interpreted.

pytensor.tensor.slinalg.cholesky(x, lower=True, *, check_finite=False, overwrite_a=False, on_error='raise')

Return a triangular matrix square root of positive semi-definite x.

L = cholesky(X, lower=True) implies dot(L, L.T) == X.

Parameters:

• x (tensor_like) –

• lower (bool, default=True) – Whether to return the lower or upper cholesky factor

• check_finite (bool, default=False) – Whether to check that the input matrix contains only finite numbers.

• overwrite_a (bool, ignored) – Whether to use the same memory for the output as a. This argument is ignored, and is present here only for consistency with scipy.linalg.cholesky.

• on_error (['raise', 'nan']) – If on_error is set to ‘raise’, this Op will raise a scipy.linalg.LinAlgError if the matrix is not positive definite. If on_error is set to ‘nan’, it will return a matrix containing nans instead.

Returns:

Lower or upper triangular Cholesky factor of x

Return type:

TensorVariable

Example

import pytensor
import pytensor.tensor as pt
import numpy as np

x = pt.tensor('x', shape=(5, 5), dtype='float64')
L = pt.linalg.cholesky(x)

f = pytensor.function([x], L)
x_value = np.random.normal(size=(5, 5))
x_value = x_value @ x_value.T # Ensures x is positive definite
L_value = f(x_value)
assert np.allclose(L_value @ L_value.T, x_value)

pytensor.tensor.slinalg.lu(a, permute_l=False, check_finite=True, p_indices=False, overwrite_a=False)

Factorize a matrix as the product of a unit lower triangular matrix and an upper triangular matrix:

… math:

A = P L U

Where P is a permutation matrix, L is lower triangular with unit diagonal elements, and U is upper triangular.

Parameters:

• a (TensorLike) – Matrix to be factorized

• permute_l (bool) – If True, L is a product of permutation and unit lower triangular matrices. Only two values, PL and U, will be returned in this case, and PL will not be lower triangular.

• check_finite (bool) – Whether to check that the input matrix contains only finite numbers.

• p_indices (bool) – If True, return integer matrix indices for the permutation matrix. Otherwise, return the permutation matrix itself.

• overwrite_a (bool) – Ignored by Pytensor. Pytensor will always perform computation inplace if possible.

Returns:

• P (TensorVariable) – Permutation matrix, or array of integer indices for permutation matrix. Not returned if permute_l is True.

• L (TensorVariable) – Lower triangular matrix, or product of permutation and unit lower triangular matrices if permute_l is True.

• U (TensorVariable) – Upper triangular matrix

pytensor.tensor.slinalg.lu_factor(a, *, check_finite=True, overwrite_a=False)

LU factorization with partial pivoting.

Parameters:

• a (TensorLike) – Matrix to be factorized

• check_finite (bool) – Whether to check that the input matrix contains only finite numbers.

• overwrite_a (bool) – Unused by PyTensor. PyTensor will always perform the operation in-place if possible.

Returns:

• LU (TensorVariable) – LU decomposition of a

• pivots (TensorVariable) – An array of integers representin the pivot indices

pytensor.tensor.slinalg.lu_solve(LU_and_pivots, b, trans=False, b_ndim=None, check_finite=True, overwrite_b=False)

Solve a system of linear equations given the LU decomposition of the matrix.

Parameters:

• LU_and_pivots (tuple[TensorLike, TensorLike]) – LU decomposition of the matrix, as returned by lu_factor

• b (TensorLike) – Right-hand side of the equation

• trans (bool) – If True, solve A^T x = b, instead of Ax = b. Default is False

• b_ndim (int, optional) – The number of core dimensions in b. Used to distinguish between a batch of vectors (b_ndim=1) and a matrix of vectors (b_ndim=2). Default is None, which will infer the number of core dimensions from the input.

• check_finite (bool) – If True, check that the input matrices contain only finite numbers. Default is True.

• overwrite_b (bool) – Ignored by Pytensor. Pytensor will always compute inplace when possible.

pytensor.tensor.slinalg.qr(A, mode='full', overwrite_a=False, pivoting=False, lwork=None)

QR Decomposition of input matrix a.

The QR decomposition of a matrix A is a factorization of the form :math`A = QR`, where Q is an orthogonal matrix (\(Q Q^T = I\)) and R is an upper triangular matrix.

This decomposition is useful in various numerical methods, including solving linear systems and least squares problems.

Parameters:

• A (TensorLike) – Input matrix of shape (M, N) to be decomposed.

• mode (str, one of "full", "economic", "r", or "raw") –

  How the QR decomposition is computed and returned. Choosing the mode can avoid unnecessary computations, depending on which of the return matrices are needed. Given input matrix with shape Choices are:

  • ”full” (or “complete”): returns Q and R with dimensions (M, M) and (M, N).

  • ”economic” (or “reduced”): returns Q and R with dimensions (M, K) and (K, N),

    where K = min(M, N).

  • ”r”: returns only R with dimensions (K, N).

  • ”raw”: returns H and tau with dimensions (N, M) and (K,), where H is the matrix of

    Householder reflections, and tau is the vector of Householder coefficients.

• pivoting (bool, default False) – If True, also return a vector of rank-revealing permutations P such that A[:, P] = QR.

• overwrite_a (bool, ignored) – Ignored. Included only for consistency with the function signature of scipy.linalg.qr. Pytensor will always automatically overwrite the input matrix A if it is safe to do sol.

• lwork (int, ignored) – Ignored. Included only for consistency with the function signature of scipy.linalg.qr. Pytensor will automatically determine the optimal workspace size for the QR decomposition.

Returns:

• Q or H (TensorVariable, optional) – A matrix with orthonormal columns. When mode = ‘complete’, it is the result is an orthogonal/unitary matrix depending on whether a is real/complex. The determinant may be either +/- 1 in that case. If mode = ‘raw’, it is the matrix of Householder reflections. If mode = ‘r’, Q is not returned.

• R or tau (TensorVariable, optional) – Upper-triangular matrix. If mode = ‘raw’, it is the vector of Householder coefficients.

pytensor.tensor.slinalg.solve(a, b, *, lower=False, overwrite_a=False, overwrite_b=False, check_finite=True, assume_a='gen', transposed=False, b_ndim=None)

Solves the linear equation set a * x = b for the unknown x for square a matrix.

If the data matrix is known to be a particular type then supplying the corresponding string to assume_a key chooses the dedicated solver. The available options are

diagonal	‘diagonal’
tridiagonal	‘tridiagonal’
banded	‘banded’
upper triangular	‘upper triangular’
lower triangular	‘lower triangular’
symmetric	‘symmetric’ (or ‘sym’)
hermitian	‘hermitian’ (or ‘her’)
positive definite	‘positive definite’ (or ‘pos’)
general	‘general’ (or ‘gen’)

If omitted, 'general' is the default structure.

The datatype of the arrays define which solver is called regardless of the values. In other words, even when the complex array entries have precisely zero imaginary parts, the complex solver will be called based on the data type of the array.

Parameters:

• a ((..., N, N) array_like) – Square input data

• b ((..., N, NRHS) array_like) – Input data for the right hand side.

• lower (bool, default False) – Ignored unless assume_a is one of 'sym', 'her', or 'pos'. If True, the calculation uses only the data in the lower triangle of a; entries above the diagonal are ignored. If False (default), the calculation uses only the data in the upper triangle of a; entries below the diagonal are ignored.

• overwrite_a (bool) – Unused by PyTensor. PyTensor will always perform the operation in-place if possible.

• overwrite_b (bool) – Unused by PyTensor. PyTensor will always perform the operation in-place if possible.

• check_finite (bool, optional) – Whether to check that the input matrices contain only finite numbers. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.

• assume_a (str, optional) – Valid entries are explained above.

• transposed (bool, default False) – If True, solves the system A^T x = b. Default is False.

• b_ndim (int) – Whether the core case of b is a vector (1) or matrix (2). This will influence how batched dimensions are interpreted. By default, we assume b_ndim = b.ndim is 2 if b.ndim > 1, else 1.

pytensor.tensor.slinalg.solve_continuous_lyapunov(A, Q)

Solve the continuous Lyapunov equation \(A X + X A^H + Q = 0\).

Parameters:

• A (TensorLike) – Square matrix of shape N x N.

• Q (TensorLike) – Square matrix of shape N x N.

Returns:

X – Square matrix of shape N x N

Return type:

TensorVariable

pytensor.tensor.slinalg.solve_discrete_are(A, B, Q, R, enforce_Q_symmetric=False)

Solve the discrete Algebraic Riccati equation \(A^TXA - X - (A^TXB)(R + B^TXB)^{-1}(B^TXA) + Q = 0\).

Discrete-time Algebraic Riccati equations arise in the context of optimal control and filtering problems, as the solution to Linear-Quadratic Regulators (LQR), Linear-Quadratic-Guassian (LQG) control problems, and as the steady-state covariance of the Kalman Filter.

Such problems typically have many solutions, but we are generally only interested in the unique stabilizing solution. This stable solution, if it exists, will be returned by this function.

Parameters:

• A (TensorLike) – Square matrix of shape M x M

• B (TensorLike) – Square matrix of shape M x M

• Q (TensorLike) – Symmetric square matrix of shape M x M

• R (TensorLike) – Square matrix of shape N x N

• enforce_Q_symmetric (bool) – If True, the provided Q matrix is transformed to 0.5 * (Q + Q.T) to ensure symmetry

Returns:

X – Square matrix of shape M x M, representing the solution to the DARE

Return type:

TensorVariable

pytensor.tensor.slinalg.solve_discrete_lyapunov(A, Q, method='bilinear')

Solve the discrete Lyapunov equation \(A X A^H - X = Q\).

Parameters:

• A (TensorLike) – Square matrix of shape N x N

• Q (TensorLike) – Square matrix of shape N x N

• method (str, one of "direct" or "bilinear") –

  Solver method used, . "direct" solves the problem directly via matrix inversion. This has a pure PyTensor implementation and can thus be cross-compiled to supported backends, and should be preferred when

  N is not large. The direct method scales poorly with the size of N, and the bilinear can be

  used in these cases.

Returns:

X – Square matrix of shape N x N. Solution to the Lyapunov equation

Return type:

TensorVariable

pytensor.tensor.slinalg.solve_triangular(a, b, *, trans=0, lower=False, unit_diagonal=False, check_finite=True, b_ndim=None)

Solve the equation a x = b for x, assuming a is a triangular matrix.

Parameters:

• a (TensorVariable) – Square input data

• b (TensorVariable) – Input data for the right hand side.

• lower (bool, optional) – Use only data contained in the lower triangle of a. Default is to use upper triangle.

• trans ({0, 1, 2, 'N', 'T', 'C'}, optional) – Type of system to solve: trans system 0 or ‘N’ a x = b 1 or ‘T’ a^T x = b 2 or ‘C’ a^H x = b

• unit_diagonal (bool, optional) – If True, diagonal elements of a are assumed to be 1 and will not be referenced.

• check_finite (bool, optional) – Whether to check that the input matrices contain only finite numbers. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.

• b_ndim (int) – Whether the core case of b is a vector (1) or matrix (2). This will influence how batched dimensions are interpreted.