ldl_solve

→Tensor

ldl_solve(LD: Tensor, pivots: Tensor, B: Tensor)

source edit

Solve a symmetric linear system using an LDL factorization.

Given $A = L\,D\,L^\top$ produced by ldl_factor, solves

A\,X \,=\, B.

sytrf's $L$ is a product $P_1 L_1 P_2 L_2 \cdots$ of elimination steps, each preceded by its own row interchange. Moving every interchange to the left gives

A \,=\, P\,\tilde L\,D\,\tilde L^\top P^\top,

with $\tilde L$ unit-lower-triangular — $L$ with each column carrying the interchanges of the steps after it, as LAPACK's syconv arranges it — so the solve is a permutation, two triangular solves around the $1 \times 1$ / $2 \times 2$ blocks of $D$ (solved in closed form, as sytrs does), and the inverse permutation.

Parameters

LDTensor

Packed LDL factor from ldl_factor, shape (*, n, n). Only its lower triangle is read.

pivotsTensor

Integer pivot indices from ldl_factor, shape (*, n) (1-based, LAPACK convention). Both

1 \times 1

and

2 \times 2

Bunch-Kaufman blocks are supported.

BTensor

Right-hand side of shape (*, n, k) (or (n,)). Its batch dimensions broadcast with those of LD.

Returns

Tensor

Solution $X$ , shaped like B broadcast against LD's batch.

Raises

ValueError

If LD is not square, pivots does not match it, or a pivot is not one sytrf could have written.

TypeError

If pivots is not an integer tensor.

Notes

Supports indefinite symmetric $A$ (unlike Cholesky), so it is appropriate for KKT / saddle-point systems where Cholesky would fail. Cost per solve is $O(n^2 k)$ once the LDL factor is in hand: two triangular solves, and $O(nk)$ for $D$ . The pivots are replayed on the host in $O(n)$ per matrix; only when an interchange reaches earlier columns of $L$ is an $n \times n$ gather index built for $\tilde L$ ( $O(n^2)$ host work).

The interchanges are applied in the order sytrf performed them (LAPACK sytrs). Applying them as a single up-front row permutation of $B$ — as this function once did — is right only when no later step swaps rows an earlier column of $L$ has entries in, and silently wrong otherwise. The pivots are read on the host with tolist; the arithmetic stays on LD's device and is differentiable in LD and B.

Examples

>>> import lucid
>>> from lucid.linalg import ldl_factor, ldl_solve
>>> A = lucid.tensor([[4.0, 1.0], [1.0, 3.0]])
>>> LD, piv = ldl_factor(A)
>>> b = lucid.tensor([[5.0], [4.0]])
>>> ldl_solve(LD, piv, b)
tensor([[1.], [1.]])

>>> import lucid >>> from lucid.linalg import ldl_factor, ldl_solve >>> A = lucid.tensor([[4.0, 1.0], [1.0, 3.0]]) >>> LD, piv = ldl_factor(A) >>> b = lucid.tensor([[5.0], [4.0]]) >>> ldl_solve(LD, piv, b) tensor([[1.], [1.]])