`scipy.sparse.linalg._funm_multiply_krylov:funm_multiply_krylov`

source: build-install/usr/lib/python3.14/site-packages/scipy/sparse/linalg/_funm_multiply_krylov.py:149

Signature

def funm_multiply_krylov(f, A, b, *, assume_a = general, t = 1.0, atol = 0.0, rtol = 1e-06, restart_every_m = None, max_restarts = 20)

Summary

A restarted Krylov method for evaluating y = f(tA) b from ^[1] ^[2].

Parameters

f : callable: Callable object that computes the matrix function F = f(X).
A : {sparse array, ndarray, LinearOperator}: A real or complex N-by-N matrix. Alternatively, A can be a linear operator which can produce Ax using, e.g., scipy.sparse.linalg.LinearOperator.
b : ndarray: A vector to multiply the f(tA) with.
assume_a : str, optional: Indicate the structure of A. The algorithm will use this information to select the appropriated code path. The available options are 'hermitian'/'her' and 'general'/'gen'. If omitted, then it is assumed that A has a 'general' structure.
t : float, optional: The value to scale the matrix A with. The default is t = 1.0
atol, rtol : float, optional: Parameters for the convergence test. For convergence, norm(||y_k - y_k-1||) <= max(rtol*norm(b), atol) should be satisfied. The default is atol=0. and rtol=1e-6.
restart_every_m : int: If the iteration number reaches this value a restart is triggered. Larger values increase iteration cost but may be necessary for convergence. If omitted, min(20, n) is used.
max_restarts : int, optional: Maximum number of restart cycles. The algorithm will stop after max_restarts cycles even if the specified tolerance has not been achieved. The default is max_restarts=20

Returns

y : ndarray: The result of f(tA) b.

Notes

The convergence of the Krylov method heavily depends on the spectrum of A and the function f. With restarting, there are only formal proofs for functions of order 1 (e.g., exp, sin, cos) and Stieltjes functions ^[2] ^[3], while the general case remains an open problem.

Examples

import numpy as np
from scipy.sparse import csr_array
from scipy.sparse.linalg import funm_multiply_krylov
from scipy.linalg import expm, solve
A = csr_array([[3, 2, 0], [1, -1, 0], [0, 5, 1]], dtype=float)
b = np.array([2, 4, -1], dtype=float)
t = 0.1

✓

Compute ``y = exp(tA) b``.

y = funm_multiply_krylov(expm, A, b, t = t)

✓

✗

array([3.61649, 3.88422, 0.96073])

ref = expm(t * A.todense()) @ b
err = y - ref

✓

err

✗

array([-0.,  0., -0.])

Compute :math:`y = (A^3 - A) b`.

poly = lambda X : X @ X @ X - X
y = funm_multiply_krylov(poly, A, b)

✓

✗

array([132.,  24.,  70.])

ref = poly(A.todense()) @ b
err = y - ref

✓

err

✗

array([ 0.,  0., -0.])

Compute :math:`y = f(tA) b`, where :math:`f(X) = X^{-1}(e^{X} - I)`. This is known as the "phi function" from the exponential integrator literature.

phim_1 = lambda X : solve(X, expm(X) - np.eye(X.shape[0]))
y = funm_multiply_krylov(phim_1, A, b, t = t)

✓

✗

array([ 2.76984,  3.92769, -0.03111])

ref = phim_1(t * A.todense()) @ b
err = y - ref

✓

err

✗

array([-0.,  0., -0.])