Linear Algebra API#

Linear algebra routines operate on SpaceCore spaces and operators. Preconditions are mathematical: square, Hermitian, positive-definite, and residual statements refer to the domain and codomain inner products, not only coordinate arrays.

Linear solves#

`spacecore.linalg.cg`	Solve \(A x = b\) by conjugate gradients.
`spacecore.linalg.CGResult`	Store the result returned by `cg()`.

cg solves A x = b for Hermitian positive-definite A : X -> X with residuals measured in X.norm.
CGResult stores x, convergence flag, iteration count, residual norm, and status details.

Least squares#

`spacecore.linalg.lsqr`	Solve \(\min_x \\|A x - b\\|\) by LSQR.
`spacecore.linalg.LSQRResult`	Store the result returned by `lsqr()`.

lsqr solves least-squares problems for A : X -> Y using X.inner and Y.inner.
LSQRResult stores the solution, convergence data, and residual diagnostics.

Eigenvalue and spectral methods#

`spacecore.linalg.lanczos_smallest`	Approximate the smallest eigenpair of a Hermitian operator.
`spacecore.linalg.power_iteration`	Estimate the dominant eigenpair of a self-adjoint action.
`spacecore.linalg.LanczosResult`	Store the result returned by `lanczos_smallest()`.
`spacecore.linalg.PowerIterationResult`	Store the result returned by `power_iteration()`.

lanczos_smallest approximates the smallest eigenpair of a Hermitian A : X -> X.
power_iteration estimates the dominant eigenpair of a self-adjoint action or quadratic-form Hessian.

Matrix functions#

`spacecore.linalg.expm_multiply`	Compute \(\exp(t A) v\) by Krylov projection.
`spacecore.linalg.ExpmMultiplyResult`	Store the result returned by `expm_multiply()`.

expm_multiply approximates exp(t A) v for square Hermitian A : X -> X using Krylov projection.

Autodoc#

spacecore.linalg.cg(A, b, *, x0=None, tol=1e-06, atol=0.0, maxiter=None, check_every=64)[source]#

Solve \(A x = b\) by conjugate gradients.

Require A to be square in the SpaceCore sense (A.domain == A.codomain), Hermitian, and positive-definite with respect to A.domain.inner. The implementation uses only LinOp.apply() and the domain-space inner product; it never materializes a dense matrix.

Parameters:

A (LinOp) – Linear operator that must be Hermitian positive-definite with respect to A.domain.inner. A.domain must equal A.codomain, including the underlying space type and inner-product geometry. An operator that is provably non-self-adjoint in this geometry (A.is_hermitian() is False) is rejected at entry with a ValueError. Operators whose Hermiticity is unknown (A.is_hermitian() is None, e.g. matrix-free operators) are accepted unchecked; positive-definiteness is likewise not validated, so indefinite or otherwise unsuitable operators can still diverge or produce NaN outputs without an explicit error.
b (array-like) – Right-hand side in A.codomain.
x0 (array-like or None, optional) – Initial guess in A.domain. Default is the zero vector.
tol (float, optional) – Relative tolerance on the linear-system residual. result.converged is True when the residual norm is below atol + tol * norm(b). Default is 1e-6.
atol (float, optional) – Absolute residual tolerance. Default is 0.0.
maxiter (int or None, optional) – Maximum number of iterations. Default is prod(A.domain.shape).
check_every (int, optional) – Refresh convergence diagnostics every this many iterations and always on the final iteration. Default is DEFAULT_CONVERGENCE_CHECK_INTERVAL.

Returns:

Named tuple with fields:

x: approximate solution in A.domain
converged: whether the requested tolerance was met
num_iters: number of iterations executed
residual_norm: final residual norm

Return type:

CGResult

Raises:

TypeError – If A is not a LinOp.
ValueError – If A is not square or if iteration parameters are invalid.

See also

lsqr: Solve least-squares systems for rectangular operators.
lanczos_smallest: Approximate the smallest eigenpair of a Hermitian operator.

Notes

The residual norm is compared with \(\text{atol} + \text{tol} \| b \|\) only every check_every iterations, and always on the final iteration. This keeps convergence checks out of the hot loop while remaining compatible with JAX JIT control flow. maxiter and check_every should be treated as static JAX arguments.

Iteration also stops when no numerically useful CG update remains: either the squared residual is at machine-precision scale or the curvature inner(p, A p) is nonpositive/tiny relative to the residual scale. The residual is refreshed before this early exit, so converged still reflects the returned iterate.

For complex operators, residual norms and step sizes are computed from the real part of A.domain.inner(x, y). SpaceCore’s complex inner-product convention conjugates the first argument; custom Space subclasses must follow that convention for CG to converge correctly.

Inner products and norms use A.domain.inner and A.domain.norm. The method is correct on non-Euclidean geometries when the space supplies Riesz maps and A is Hermitian positive-definite in that geometry.

References

Hestenes, M. R. and Stiefel, E., “Methods of Conjugate Gradients for Solving Linear Systems,” J. Res. Natl. Bur. Stand., 49 (1952), 409-436.

Examples

Solve a small positive-definite system.

>>> import numpy as np
>>> import spacecore as sc
>>> ctx = sc.Context(sc.NumpyOps(), dtype=np.float64)
>>> X = sc.DenseCoordinateSpace((3,), ctx)
>>> M = ctx.asarray([[4.0, 1.0, 0.0], [1.0, 3.0, 1.0], [0.0, 1.0, 2.0]])
>>> A = sc.DenseLinOp(M, X, X, ctx)
>>> b = ctx.asarray([1.0, 2.0, 3.0])
>>> result = sc.cg(A, b, tol=1e-10)
>>> bool(result.converged)
True
>>> np.allclose(A.apply(result.x), b)
True

class spacecore.linalg.CGResult(x, converged, num_iters, residual_norm)[source]#

Bases: NamedTuple

Store the result returned by cg().

Parameters:

x (array-like) – Approximate solution in A.domain.
converged (bool-like) – Whether the final residual norm satisfied the requested tolerance.
num_iters (int-like) – Number of conjugate-gradient iterations executed.
residual_norm (scalar) – Norm of the final residual in A.codomain.

x: Any#: Alias for field number 0

converged: Any#: Alias for field number 1

num_iters: Any#: Alias for field number 2

residual_norm: Any#: Alias for field number 3

count(value, /)#: Return number of occurrences of value.

index(value, start=0, stop=sys.maxsize, /)#

Return first index of value.

Raises ValueError if the value is not present.

spacecore.linalg.lsqr(A, b, *, x0=None, tol=1e-06, atol=0.0, maxiter=None, check_every=64, residual_mode='exact')[source]#

Solve \(\min_x \|A x - b\|\) by LSQR.

Allow A to map between distinct domain and codomain spaces. The method uses LinOp.apply() for forward products and A.H.apply for adjoint products, so the normal equations are represented implicitly and no dense matrix is formed.

Parameters:

A (LinOp) – Linear operator with possibly distinct domain and codomain. For square A (A.domain == A.codomain), cg() is usually preferred when A is also Hermitian positive-definite.
b (array-like) – Right-hand side in A.codomain.
x0 (array-like or None, optional) – Initial guess in A.domain. Default is the zero vector.
tol (float, optional) – Relative tolerance for the normal-equation residual norm(A.H @ (A @ x - b)). result.converged is True when that residual is below atol + tol * norm(b). Default is 1e-6.
atol (float, optional) – Absolute tolerance for the normal-equation residual. Default is 0.0.
maxiter (int or None, optional) – Maximum number of iterations. Default is prod(A.domain.shape).
check_every (int, optional) – Refresh residual diagnostics every this many iterations and always on the final iteration. Default is DEFAULT_CONVERGENCE_CHECK_INTERVAL.
residual_mode ({"exact", "recurrence"}, optional) – Residual diagnostic mode. "exact" preserves the historical behavior: every diagnostic refresh recomputes A @ x - b and A.H @ (A @ x - b). This costs one additional forward application and one additional adjoint application on each check iteration, so small check_every values, especially check_every=1, can substantially increase runtime for expensive operators. Use larger values such as check_every=10 or check_every=20 when exact diagnostics are needed for matrix-free, PDE, neural-network, GPU, or JAX workloads. "recurrence" uses LSQR scalar recurrences for both returned residual diagnostics and avoids those extra applications.

Returns:

Named tuple with fields:

x: approximate least-squares solution in A.domain
converged: whether the requested tolerance was met
num_iters: number of iterations executed
residual_norm: final residual norm or recurrence estimate
normal_residual_norm: final normal-equation residual norm or recurrence estimate

Return type:

LSQRResult

Raises:

TypeError – If A is not a LinOp.
ValueError – If iteration parameters are invalid or residual_mode is unknown.

See also

cg: Solve square Hermitian positive-definite systems.
power_iteration: Estimate a dominant eigenpair.

Notes

Convergence is tested using \(\|A^*(A x - b)\| < \text{atol} + \text{tol}\|b\|\). In residual_mode="exact", exact residual diagnostics are refreshed only every check_every iterations, and always on the final iteration. Each refresh performs one additional forward product and one additional adjoint product beyond the LSQR recurrence itself.

In residual_mode="recurrence", residual_norm is the standard LSQR estimate abs(phi_bar) and normal_residual_norm is the LSQR scalar estimate alpha * abs(tau) with tau = s * phi. These estimates avoid extra operator applications during checks, including the final check. This function is JIT-compatible on the JAX backend when maxiter, check_every, and residual_mode are static arguments.

The normal-equation residual can be much smaller than the solution error for ill-conditioned A. For ill-conditioned problems, use a tighter tol or check the residual and solution quality directly.

Works on real and complex operators. For complex operators, A.H uses the conjugate adjoint.

Inner products and norms use A.domain.inner / A.domain.norm for domain-space quantities and A.codomain.norm for least-squares residuals. The method is therefore correct on non-Euclidean geometries when the spaces provide Riesz maps and A.rapply is the true metric adjoint.

References

Paige, C. C. and Saunders, M. A., “LSQR: An Algorithm for Sparse Linear Equations and Sparse Least Squares,” ACM Trans. Math. Soft., 8 (1982), 43-71.

Examples

Solve a small overdetermined least-squares problem.

>>> import numpy as np
>>> import spacecore as sc
>>> ctx = sc.Context(sc.NumpyOps(), dtype=np.float64)
>>> X = sc.DenseCoordinateSpace((2,), ctx)
>>> Y = sc.DenseCoordinateSpace((3,), ctx)
>>> M = ctx.asarray([[1.0, 0.0], [0.0, 1.0], [1.0, 1.0]])
>>> A = sc.DenseLinOp(M, X, Y, ctx)
>>> b = ctx.asarray([1.0, 2.0, 3.0])
>>> result = sc.lsqr(A, b, tol=1e-10)
>>> np.allclose(result.x, [1.0, 2.0])
True

class spacecore.linalg.LSQRResult(x, converged, num_iters, residual_norm, normal_residual_norm)[source]#

Bases: NamedTuple

Store the result returned by lsqr().

Parameters:

x (array-like) – Approximate least-squares solution in A.domain.
converged (bool-like) – Whether the normal-equation residual satisfied the requested tolerance.
num_iters (int-like) – Number of LSQR iterations executed.
residual_norm (scalar) – Norm of A x - b in A.codomain in exact mode, or the LSQR recurrence estimate in recurrence mode.
normal_residual_norm (scalar) – Norm of A.H @ (A x - b) in A.domain in exact mode, or the LSQR recurrence estimate in recurrence mode.

x: Any#: Alias for field number 0

converged: Any#: Alias for field number 1

num_iters: Any#: Alias for field number 2

residual_norm: Any#: Alias for field number 3

normal_residual_norm: Any#: Alias for field number 4

count(value, /)#: Return number of occurrences of value.

index(value, start=0, stop=sys.maxsize, /)#

Return first index of value.

Raises ValueError if the value is not present.

spacecore.linalg.lanczos_smallest(A, initial_vector, *, max_iter=100, tol=1e-06, check_every=64)[source]#

Approximate the smallest eigenpair of a Hermitian operator.

The operator is supplied as a square LinOp in the SpaceCore sense (A.domain == A.codomain), and initial_vector is an element of A.domain. The implementation keeps fixed-size coordinate arrays for JAX compatibility, safely handles zero initial vectors, and refines the returned eigenvalue with the Rayleigh quotient of the reconstructed Ritz vector in the original space.

Mathematically, Lanczos builds an orthonormal Krylov basis V for span{v, A v, A^2 v, ...} and a tridiagonal projection \(T_k = V^* A V\). The returned vector is the Ritz vector reconstructed in the original coordinates, and the returned scalar is the Rayleigh quotient \(\langle x, A x \rangle_X / \langle x, x \rangle_X\).

Parameters:

A (LinOp) – Linear operator that must be Hermitian/self-adjoint with respect to A.domain.inner. A.domain must equal A.codomain, including the underlying space type and inner-product geometry. Operators with structurally unknown Hermiticity (A.is_hermitian() returns None) are accepted on trust; the caller is responsible for ensuring Hermiticity. Non-Hermitian inputs produce undefined results.
initial_vector (array-like) – Starting vector in A.domain. If it is numerically zero, the algorithm falls back to a deterministic coordinate vector.
max_iter (int, optional) – Maximum Krylov dimension. Must be a Python int rather than a traced JAX scalar; under jax.jit it is treated as a static argument and changing it triggers retracing. Default is 100.
tol (float, optional) – Tolerance used for two purposes. Iteration stops at a check point when the off-diagonal Lanczos coefficient falls below tol; the returned converged flag is True when the Ritz residual estimate is below tol. Default is 1e-6.
check_every (int, optional) – Refresh the breakdown-based stopping decision every this many iterations and always on the final iteration. Default is DEFAULT_CONVERGENCE_CHECK_INTERVAL.

Returns:

Named tuple with fields:

eigenvalue: smallest Ritz eigenvalue estimate
eigenvector: associated Ritz vector in A.domain
residual_norm: standard Ritz residual estimate
krylov_dim: actual Krylov dimension reached
converged: whether residual_norm < tol

Return type:

LanczosResult

Raises:

TypeError – If A is not a LinOp.
ValueError – If A is not square, is known to be non-Hermitian, or if max_iter is invalid.

See also

power_iteration: Estimate the dominant eigenpair.
expm_multiply: Apply a matrix exponential using the Lanczos basis.

Notes

The residual estimate is computed from the tridiagonal recurrence as \(\beta_m |y_{m-1}|\). Callers that need the true residual can evaluate A.apply(eigenvector) - eigenvalue * eigenvector once more in the original space.

The “smallest Ritz value” is the smallest eigenvalue of the projected tridiagonal matrix, not necessarily a good approximation of the smallest eigenvalue of A. Convergence to the actual smallest eigenvalue requires the bottom of the spectrum to be separated and the initial vector to have nonzero projection onto the corresponding eigenspace. For clustered low eigenvalues, increase max_iter or use multiple initial vectors.

Hermiticity is enforced only when it can be structurally verified: known non-Hermitian operators raise ValueError. Operators with unknown structure, such as many matrix-free operators and operators on custom spaces, are trusted.

This function is JIT-compatible on the JAX backend when max_iter and check_every are static arguments. For plain VectorSpace domains, Euclidean reorthogonalization is vectorized; custom spaces use Space.inner() to preserve the declared geometry.

Inner products and norms use A.domain.inner and A.domain.norm. The method is correct on non-Euclidean geometries when the space supplies Riesz maps and A is self-adjoint in that geometry.

References

Lanczos, C., “An Iteration Method for the Solution of the Eigenvalue Problem of Linear Differential and Integral Operators,” J. Res. Natl. Bur. Stand., 45 (1950), 255-282.

Examples

Approximate the smallest eigenpair of a diagonal operator.

>>> import numpy as np
>>> import spacecore as sc
>>> ctx = sc.Context(sc.NumpyOps(), dtype=np.float64)
>>> X = sc.DenseCoordinateSpace((3,), ctx)
>>> A = sc.DiagonalLinOp(ctx.asarray([1.0, 2.0, 4.0]), X, ctx)
>>> result = sc.lanczos_smallest(A, ctx.asarray([1.0, 1.0, 1.0]), max_iter=3)
>>> np.allclose(result.eigenvalue, 1.0)
True

class spacecore.linalg.LanczosResult(eigenvalue, eigenvector, residual_norm, krylov_dim, converged)[source]#

Bases: NamedTuple

Store the result returned by lanczos_smallest().

Parameters:

eigenvalue (scalar) – Ritz approximation to the smallest eigenvalue.
eigenvector (array-like) – Ritz vector in A.domain.
residual_norm (scalar) – Standard Ritz residual estimate.
krylov_dim (int-like) – Krylov dimension reached before breakdown or max_iter.
converged (bool-like) – Whether residual_norm < tol.

eigenvalue: Any#: Alias for field number 0

eigenvector: Any#: Alias for field number 1

residual_norm: Any#: Alias for field number 2

krylov_dim: Any#: Alias for field number 3

converged: Any#: Alias for field number 4

count(value, /)#: Return number of occurrences of value.

index(value, start=0, stop=sys.maxsize, /)#

Return first index of value.

Raises ValueError if the value is not present.

spacecore.linalg.power_iteration(A, *, x0=None, tol=1e-06, maxiter=None, check_every=64)[source]#

Estimate the dominant eigenpair of a self-adjoint action.

Accept a square LinOp or a QuadraticForm exposing hess_apply. Public dispatch converts either input into a fixed self-adjoint action before entering the numerical loop. Power iteration still requires hess_apply for quadratic forms because the vector update is x_next = normalize(H x). Optional two-sided scalar diagnostics such as hess_quad(x, Hx=None) can improve Rayleigh quotient evaluation, but they cannot replace the Hessian-vector action. “Dominant” means largest eigenvalue in absolute value, not necessarily the largest positive eigenvalue.

Parameters:

A (LinOp or QuadraticForm) – Square operator or quadratic form whose dominant eigenpair, largest in absolute value, is sought. Linear-operator inputs must satisfy A.domain == A.codomain; this includes the underlying space type and inner-product geometry. Quadratic-form inputs must provide hess_apply. If they also expose hess_quad(x, Hx=None), it is used for the Rayleigh quotient and must be compatible with being called as hess_quad(x, Hx=Hx). The Hx argument is the cached Hessian-vector product already computed by power iteration. If the quadratic form exposes hess_residual_norm(x, Hx, eigenvalue), it is used for the residual diagnostic. Otherwise generic space inner-product and norm diagnostics are used. For spectral-norm estimates of a rectangular operator, pass A.H @ A.
x0 (array-like or None, optional) – Initial vector in the action domain. Default is a normalized all-ones vector in the domain geometry.
tol (float, optional) – Residual-norm tolerance. result.converged is True when norm(A @ x - lambda * x) < tol. Default is 1e-6.
maxiter (int or None, optional) – Maximum number of iterations. Default is prod(A.domain.shape).
check_every (int, optional) – Accepted for backward compatibility. Residual diagnostics are now refreshed every iteration because the loop already carries A @ x; this argument is ignored and may be removed in a future release.

Returns:

Named tuple with fields:

eigenvalue: Rayleigh-quotient eigenvalue estimate
eigenvector: normalized eigenvector estimate
converged: whether residual_norm < tol
num_iters: number of iterations executed
residual_norm: norm of A x - eigenvalue * x

Return type:

PowerIterationResult

Raises:

TypeError – If A is neither a LinOp nor a QuadraticForm.
ValueError – If a linear-operator input is not square, is known to be non-Hermitian, or if iteration parameters are invalid.

See also

lanczos_smallest: Approximate the smallest eigenpair of a Hermitian operator.
cg: Solve Hermitian positive-definite systems.

Notes

The residual-based stopping criterion uses \(\|A x - \lambda x\|\) and is refreshed every iteration. The check_every argument is accepted for backward compatibility but is no longer used. This function is JIT-compatible on the JAX backend when maxiter is static.

Inner products and norms use domain.inner and domain.norm through the normalized self-adjoint action. The method is correct on non-Euclidean geometries when the space supplies Riesz maps and the action is self-adjoint in that geometry.

For operators with eigenvalues of mixed sign, the dominant eigenvalue is the one with largest absolute value, which may be negative. Convergence requires that this eigenvalue be separated from the rest in absolute value. If the dominant modulus is degenerate, for example both lambda and -lambda have maximum modulus, the iteration may oscillate between subspaces.

For most dense vector-space problems, the generic Rayleigh quotient real(inner(x, Hx)) is already cheap. Specialized quadratic-form scalar diagnostics are mainly useful when a subclass can evaluate <x, Hx> or the residual norm more accurately or with less overhead than reconstructing the scalar from generic space operations. A specialized hess_quad must accept the cached Hessian-vector product as hess_quad(x, Hx=Hx); implementations may ignore Hx or use it to avoid recomputation.

Examples

Estimate the largest eigenvalue of a diagonal operator.

>>> import numpy as np
>>> import spacecore as sc
>>> ctx = sc.Context(sc.NumpyOps(), dtype=np.float64)
>>> X = sc.DenseCoordinateSpace((3,), ctx)
>>> A = sc.DiagonalLinOp(ctx.asarray([1.0, 3.0, 2.0]), X, ctx)
>>> result = sc.power_iteration(A, maxiter=20, tol=1e-10)
>>> np.allclose(result.eigenvalue, 3.0)
True

class spacecore.linalg.PowerIterationResult(eigenvalue, eigenvector, converged, num_iters, residual_norm)[source]#

Bases: NamedTuple

Store the result returned by power_iteration().

Parameters:

eigenvalue (scalar) – Rayleigh-quotient estimate of the dominant eigenvalue.
eigenvector (array-like) – Normalized eigenvector estimate in the operator domain.
converged (bool-like) – Whether the residual norm satisfied tol.
num_iters (int-like) – Number of power iterations executed.
residual_norm (scalar) – Norm of A x - eigenvalue * x.

eigenvalue: Any#: Alias for field number 0

eigenvector: Any#: Alias for field number 1

converged: Any#: Alias for field number 2

num_iters: Any#: Alias for field number 3

residual_norm: Any#: Alias for field number 4

count(value, /)#: Return number of occurrences of value.

index(value, start=0, stop=sys.maxsize, /)#

Return first index of value.

Raises ValueError if the value is not present.

spacecore.linalg.expm_multiply(A, v, t=1.0, *, max_iter=30, tol=1e-10)[source]#

Compute \(\exp(t A) v\) by Krylov projection.

Require A to be square in the SpaceCore sense (A.domain == A.codomain) and Hermitian with respect to A.domain.inner. The method builds a Lanczos basis and applies the exponential of the small tridiagonal projection, avoiding dense materialization of A.

Parameters:

A (LinOp) – Linear operator that must be Hermitian/self-adjoint with respect to A.domain.inner. A.domain must equal A.codomain, including the underlying space type and inner-product geometry. Operators with structurally unknown Hermiticity (A.is_hermitian() returns None) are accepted on trust; the caller is responsible for ensuring Hermiticity. Non-Hermitian inputs produce undefined results.
v (array-like) – Initial vector in A.domain.
t (float or complex, optional) – Scalar multiplier on A. Complex values require a complex-valued ctx.dtype such as complex64 or complex128. Using a complex t with a real-valued context produces backend-dependent results. Default is 1.0.
max_iter (int, optional) – Maximum Krylov dimension. Values around 20-50 are usually sufficient when \(|t|\|A\|\) is moderate. Must be a Python int rather than a traced JAX scalar; under jax.jit it is treated as a static argument and changing it triggers retracing. Default is 30.
tol (float, optional) – Tolerance used both for Lanczos breakdown and for the convergence flag: result.converged is True when the projected exponential residual estimate is below tol. Default is 1e-10.

Returns:

Result vector in A.domain, the Krylov dimension used, the standard estimate abs(beta[m] * phi[m - 1]), and a convergence flag.

Return type:

ExpmMultiplyResult

Raises:

TypeError – If A is not a LinOp.
ValueError – If A is not square, is known to be non-Hermitian, or if max_iter is invalid.

See also

lanczos_smallest: Build the related Hermitian Krylov projection.
power_iteration: Estimate a dominant eigenpair.

Notes

The projected exponential is computed as \(\exp(t T) e_0\) using an eigendecomposition of the small real symmetric tridiagonal matrix T. This is JIT-compatible on the JAX backend when max_iter is static.

Hermiticity is enforced only when it can be structurally verified: known non-Hermitian operators raise ValueError. Operators with unknown structure, such as many matrix-free operators and operators on custom spaces, are trusted.

The returned residual estimate is \(|\beta_m \phi_{m-1}|\), where phi is the projected exponential vector. Callers that need the true residual can perform one additional operator application.

Examples

Apply the exponential of a diagonal operator.

>>> import numpy as np
>>> import spacecore as sc
>>> ctx = sc.Context(sc.NumpyOps(), dtype=np.float64)
>>> X = sc.DenseCoordinateSpace((2,), ctx)
>>> A = sc.DiagonalLinOp(ctx.asarray([0.0, 1.0]), X, ctx)
>>> v = ctx.asarray([2.0, 3.0])
>>> result = sc.expm_multiply(A, v, t=0.5, max_iter=5)
>>> np.allclose(result.result, [2.0, 3.0 * np.exp(0.5)], atol=1e-10)
True

class spacecore.linalg.ExpmMultiplyResult(result, krylov_dim, residual_estimate, converged)[source]#

Bases: NamedTuple

Store the result returned by expm_multiply().

Parameters:

result (array-like) – Vector in the domain of the input operator approximating exp(t * A) @ v.
krylov_dim (int-like) – Actual Krylov dimension reached before breakdown or max_iter.
residual_estimate (scalar) – Projected exponential residual estimate abs(beta[m] * phi[m - 1]).
converged (bool-like) – Boolean indicating whether residual_estimate < tol.

result: Any#: Alias for field number 0

krylov_dim: Any#: Alias for field number 1

residual_estimate: Any#: Alias for field number 2

converged: Any#: Alias for field number 3

count(value, /)#: Return number of occurrences of value.

index(value, start=0, stop=sys.maxsize, /)#

Return first index of value.

Raises ValueError if the value is not present.