fips.linalg.KroneckerOperator.solve

KroneckerOperator.solve(B)

Solve self @ X = B through the factors.

Each factor is factored once (Cholesky, or a symmetric LDL^T with a logged warning when roundoff has tipped it off positive-definite) and the factorization is kept. When the rows are the whole product, in any order, the inverse is \(D^{-1} P (C_1^{-1} \otimes \dots \otimes C_K^{-1}) P^T D^{-1}\).

When the rows are a strict subset of the product, with d entries of the product missing, the inverse of the kept rows’ correlation is the Schur complement \(G_{kk} - G_{kd} G_{dd}^{-1} G_{dk}\) of the full product’s inverse \(G\), at the cost of d products with the factors’ inverses and one d x d solve. That is the inverse in exact arithmetic, but in floating point it loses accuracy when the factors are ill-conditioned (the complement subtracts large terms of \(G\)). So the solve checks its result: when the backward error of any column, \(\|R y - c\| / (\|R\| \|y\| + \|c\|)\) for the rows’ correlation \(R\), is larger than n * eps (the bound a dense Cholesky solve meets), it logs a warning and solves the dense matrix instead. When the complement would hold more than the dense matrix (N d > n^2), the dense matrix is solved from the start.

Parameters:

B (array-like, SciPy sparse matrix or Operator) – The right-hand side: a vector, or a matrix with as many rows as the operator.

Returns:

X, dense, with the shape of B.

Return type:

ndarray

Raises:

numpy.linalg.LinAlgError – If a standard deviation is zero: the operator is singular.