Linear algebra module

Module for linear algebra related LazyLinOps.

lazylinop.wip.linalg.coshm(L, scale=1.0, nmax=8, backend='scipy', use_numba=False)

Constructs an hyperbolic cosine of linear operator L as a lazy linear operator C(L). It uses the equation 2 * cosh(z) = exp(scale * L) + exp(-scale * L). Of note, it is only an approximation (see nmax argument).

Args:

L: 2d array Linear operator. scale: float, optional Scale factor cosh(scale * L) (default is 1). nmax: int, optional Stop the serie expansion after nmax (default is 8). backend: str, optional It can be ‘scipy’ (default) to use scipy.linalg.coshm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of coshm(scale * L).

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.
ValueError: nmax must be >= 1.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import coshm
>>> scale = 0.01
>>> N = 10
>>> L = np.eye(N, N, k=0)
>>> E1 = coshm(L, scale=scale, nmax=32, backend='serie')
>>> E2 = sp.linalg.coshm(scale * L)
>>> E3 = coshm(eye(N, N, k=0), scale=scale, nmax=32, backend='serie')
>>> X = np.random.rand(N, 2 * N)
>>> np.allclose(E1 @ X, E2 @ X)
True
>>> np.allclose(E2 @ X, E3 @ X)
True

References:

See also scipy.linalg.coshm function.

lazylinop.wip.linalg.cosm(L, scale=1.0, nmax=8, backend='scipy')

Constructs a cosinus of linear operator L as a lazy linear operator C(L). It uses the equation expm(i * scale * L) = cos(scale * L) + i * sin(scale * L) where i^2 = -1 and returns real part. Of note, it is only an approximation (see nmax argument).

Args:

L: 2d array Linear operator. scale: float, optional Scale factor cosm(scale * L) (default is 1). nmax: int, optional Stop the serie expansion after nmax (default is 8). backend: str, optional It can be ‘scipy’ (default) to use scipy.linalg.cosm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of cosm(scale * L).

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.
ValueError: If L is a 2d array, backend must be ‘scipy’.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import cosm
>>> scale = 0.01
>>> coefficients = np.array([1.0, scale, 0.5 * scale ** 2])
>>> N = 10
>>> L = np.eye(N, N, k=0)
>>> E1 = cosm(L, scale=scale, nmax=4)
>>> E2 = sp.linalg.cosm(scale * L)
>>> np.allclose(E1.toarray(), E2)
>>> E3 = eye(N, N, k=0)
>>> X = np.random.rand(N, 2 * N)
>>> np.allclose(E2 @ X, E3 @ X)

References:

See also scipy.linalg.cosm function. See also expm().

lazylinop.wip.linalg.expm(L, scale=1.0, nmax=8, backend='scipy', use_numba=False)

Constructs exponentiation of linear operator L as a lazy linear operator E(L). Of note, it is only an approximation E(L) @ X ~= sum((scale * L)^i / factorial(i), i=0 to nmax) @ X.

Args:

L: 2d array: Linear operator.
scale: float, optional: Scale factor expm(scale * L) (default is 1).
nmax: int, optional: Stop the serie expansion after nmax (default is 8).
backend: str, optional: It can be ‘scipy’ (default) to use scipy.linalg.expm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of expm(scale * L).
use_numba: bool, optional: If True, use prange from Numba (default is False) to compute batch in parallel. It is useful only and only if the batch size and the length of a are big enough.

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.
ValueError: If L is a 2d array, backend must be ‘scipy’.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import expm
>>> scale = 0.01
>>> coefficients = np.array([1.0, scale, 0.5 * scale ** 2])
>>> N = 10
>>> L = np.eye(N, N, k=0)
>>> E1 = expm(L, scale=scale, nmax=4)
>>> E2 = sp.linalg.expm(scale * L)
>>> np.allclose(E1.toarray(), E2)
>>> E3 = eye(N, N, k=0)
>>> X = np.random.rand(N, 2 * N)
>>> np.allclose(E2 @ X, E3 @ X)

References:

See also scipy.linalg.expm function.

lazylinop.wip.linalg.logm(L, scale=1.0, nmax=8, backend='scipy', use_numba=False)

Constructs logarithm of linear operator L as a lazy linear operator Log(L). Of note, it is only an approximation Log(L) @ X ~= sum((-1)^(n + 1) * (L - Id)^n / n, n=1 to nmax) @ X.

Args:

L: 2d array Linear operator. scale: float, optional Scale factor logm(scale * L) (default is 1). nmax: int, optional Stop the serie expansion after nmax (default is 8). backend: str, optional It can be ‘scipy’ (default) to use scipy.linalg.logm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of logm(scale * L).

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.
ValueError: nmax must be >= 1.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import logm
>>> scale = 0.01
>>> N = 10
>>> E1 = logm(eye(N, N, k=0), scale=scale, nmax=4, backend='scipy')
>>> E2 = sp.linalg.logm(scale * np.eye(N, N, k=0))
>>> X = np.random.rand(N, 2 * N)
>>> np.allclose(E1 @ X, E2 @ X)

References:

See also scipy.linalg.logm function. See also logarithm of a matrix.

lazylinop.wip.linalg.sinhm(L, scale=1.0, nmax=8, backend='scipy', use_numba=False)

Constructs an hyperbolic cosine of linear operator L as a lazy linear operator S(L). It uses the equation 2 * sinh(z) = exp(scale * L) - exp(-scale * L). Of note, it is only an approximation (see nmax argument).

Args:

L: 2d array Linear operator. scale: float, optional Scale factor cosh(scale * L) (default is 1). nmax: int, optional Stop the serie expansion after nmax (default is 8). backend: str, optional It can be ‘scipy’ (default) to use scipy.linalg.sinhm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of sinhm(scale * L).

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import sinhm
>>> scale = 0.01
>>> N = 10
>>> L = np.eye(N, N, k=0)
>>> E1 = sinhm(L, scale=scale, nmax=32, backend='serie')
>>> E2 = sp.linalg.sinhm(scale * L)
>>> E3 = sinhm(eye(N, N, k=0), scale=scale, nmax=32, backend='serie')
>>> X = np.random.rand(N, 2 * N)
>>> np.allclose(E1 @ X, E2 @ X)
True
>>> np.allclose(E2 @ X, E3 @ X)
True

References:

See also scipy.linalg.sinhm function.

lazylinop.wip.linalg.sinm(L, scale=1.0, nmax=8, backend='scipy')

Constructs a cosinus of linear operator L as a lazy linear operator C(L). It uses the equation expm(i * scale * L) = cos(scale * L) + i * sin(scale * L) where i^2 = -1 and returns imaginary part. Of note, it is only an approximation (see nmax argument).

Args:

L: 2d array Linear operator. scale: float, optional Scale factor sinm(scale * L) (default is 1). nmax: int, optional Stop the serie expansion after nmax (default is 8). backend: str, optional It can be ‘scipy’ (default) to use scipy.linalg.sinm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of sinm(scale * L).

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.
ValueError: If L is a 2d array, backend must be ‘scipy’.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import sinm
>>> scale = 0.01
>>> coefficients = np.array([1.0, scale, 0.5 * scale ** 2])
>>> N = 10
>>> L = np.eye(N, N, k=0)
>>> E1 = sinm(L, scale=scale, nmax=4)
>>> E2 = sp.linalg.sinm(scale * L)
>>> np.allclose(E1.toarray(), E2)
>>> E3 = eye(N, N, k=0)
>>> X = np.random.rand(N, 2 * N)
>>> np.allclose(E2 @ X, E3 @ X)

References:

See also scipy.linalg.sinm function. See also expm().

lazylinop.wip.linalg.sqrtm(L, scale=1.0, nmax=8, backend='scipy')

Constructs square root of linear operator L as a lazy linear operator R(L). It uses the equation L^1/2 = sum((-1)^n * (1/2 n) * (Id - L)^n, n=0 to inf). Of note, it is only an approximation (see nmax argument).

Args:

L: 2d array Linear operator. scale: float, optional Scale factor cosh(scale * L) (default is 1). nmax: int, optional Stop the serie expansion after nmax (default is 8). backend: str, optional It can be ‘scipy’ (default) to use scipy.linalg.sqrtm function. nmax parameter is useless if backend is ‘scipy’. It can be ‘serie’ to use a serie expansion of sqrtm(scale * L).

Returns:

LazyLinOp

Raises:

ValueError: L @ x does not work because # of columns of L is not equal to the # of rows of x.
ValueError: backend value is either ‘scipy’ or ‘serie’.
ValueError: If L is a 2d array, backend must be ‘scipy’.

Examples:

>>> import numpy as np
>>> import scipy as sp
>>> from lazylinop import eye
>>> from lazylinop.wip.linear_algebra import sqrtm
>>> scale = 0.1
>>> N = 4
>>> L = np.eye(N, N, k=0)
>>> X = np.random.rand(N, 2 * N)
>>> E1 = sqrtm(L, scale=scale, nmax=256, backend='serie')
>>> E2 = sp.linalg.sqrtm(scale * L)
>>> np.allclose(E1 @ X, E2 @ X)
>>> L = eye(N, N, k=0)
>>> E3 = sqrtm(L, scale=scale, nmax=256, backend='serie')
>>> np.allclose(E2 @ X, E3 @ X)

References: