`krypy.utils` - Krylov Subspace Utilities¶

The utils module provides helper functions for common tasks in the process of solving linear algebraic systems.

Collection of standard functions.

This method provides functions like inner products, norms, …

exception krypy.utils.ArgumentError¶

Bases: Exception

Raised when an argument is invalid.

Analogue to ValueError which is not used here in order to be able to distinguish between built-in errors and krypy errors.

exception krypy.utils.AssumptionError¶

Bases: Exception

Raised when an assumption is not satisfied.

Differs from ArgumentError in that all passed arguments are valid but computations reveal that assumptions are not satisfied and the result cannot be computed.

exception krypy.utils.ConvergenceError(msg, solver)¶

Bases: Exception

Raised when a method did not converge.

The ConvergenceError holds a message describing the error and the attribute solver through which the last approximation and other relevant information can be retrieved.

exception krypy.utils.LinearOperatorError¶

Bases: Exception

Raised when a LinearOperator cannot be applied.

exception krypy.utils.InnerProductError¶

Bases: Exception

Raised when the inner product is indefinite.

exception krypy.utils.RuntimeError¶

Bases: Exception

Raised for errors that do not fit in any other exception.

class krypy.utils.Arnoldi(A, v, maxiter=None, ortho='mgs', M=None, Mv=None, Mv_norm=None, ip_B=None)¶

Bases: object

Arnoldi algorithm.

Computes V and H such that \(AV_n=V_{n+1}\underline{H}_n\). If the Krylov subspace becomes A-invariant then V and H are truncated such that \(AV_n = V_n H_n\).

Parameters:

A – a linear operator that can be used with scipy’s aslinearoperator with shape==(N,N).
v – the initial vector with shape==(N,1).
maxiter – (optional) maximal number of iterations. Default: N.
ortho –
(optional) orthogonalization algorithm: may be one of
- 'mgs': modified Gram-Schmidt (default).
- 'dmgs': double Modified Gram-Schmidt.
- 'lanczos': Lanczos short recurrence.
- 'house': Householder.
M – (optional) a self-adjoint and positive definite preconditioner. If M is provided, then also a second basis \(P_n\) is constructed such that \(V_n=MP_n\). This is of importance in preconditioned methods. M has to be None if ortho=='house' (see B).
ip_B –
(optional) defines the inner product to use. See inner().

ip_B has to be None if ortho=='house'. It’s unclear to me (andrenarchy), how a variant of the Householder QR algorithm can be used with a non-Euclidean inner product. Compare http://math.stackexchange.com/questions/433644/is-householder-orthogonalization-qr-practicable-for-non-euclidean-inner-products

advance()¶: Carry out one iteration of Arnoldi.

get()¶

get_last()¶

class krypy.utils.BoundCG(evals, exclude_zeros=False)¶

Bases: object

CG residual norm bound.

Computes the \(\kappa\)-bound for the CG error \(A\)-norm when the eigenvalues of the operator are given, see [LieS13].

Parameters:

evals – an array of eigenvalues \(\lambda_1,\ldots,\lambda_N\in\mathbb{R}\). The eigenvalues will be sorted internally such that \(0=\lambda_1=\ldots=\lambda_{t-1}<\lambda_t\leq\ldots\lambda_N\) for \(t\in\mathbb{N}\).
steps – (optional) the number of steps \(k\) to compute the bound for. If steps is None (default), then \(k=N\) is used.

Returns:

array \([\eta_0,\ldots,\eta_k]\) with

\[\eta_n = 2 \left( \frac{\sqrt{\kappa_{\text{eff}}} - 1} {\sqrt{\kappa_{\text{eff}}} + 1} \right)^n \quad\text{for}\quad n\in\{0,\ldots,k\}\]

where \(\kappa_{\text{eff}}=\frac{\lambda_N}{\lambda_t}\).

Initialize with array/list of eigenvalues or Intervals object.

eval_step(step)¶: Evaluate bound for given step.

get_step(tol)¶: Return step at which bound falls below tolerance.

class krypy.utils.BoundMinres(evals)¶

Bases: object

MINRES residual norm bound.

Computes a bound for the MINRES residual norm when the eigenvalues of the operator are given, see [Gre97].

Parameters:

evals – an array of eigenvalues \(\lambda_1,\ldots,\lambda_N\in\mathbb{R}\). The eigenvalues will be sorted internally such that \(\lambda_1\leq\ldots\lambda_s<0=\lambda_{s+1}=\ldots=\lambda_{s+t-1}<\lambda_t\leq\ldots\lambda_N\) for \(s,t\in\mathbb{N}\) and \(s<t\).
steps – (optional) the number of steps \(k\) to compute the bound for. If steps is None (default), then \(k=N\) is used.

Returns:

array \([\eta_0,\ldots,\eta_k]\) with

\[\eta_n = 2 \left( \frac{ \sqrt{|\lambda_1\lambda_N|} - \sqrt{|\lambda_s\lambda_t|}} { \sqrt{|\lambda_1\lambda_N|} + \sqrt{|\lambda_s\lambda_t|}} \right)^{\left[\frac{n}{2}\right]} \quad\text{for}\quad n\in\{0,\ldots,k\}\]

if \(s>0\). If \(s=0\), i.e., if the eigenvalues are non-negative, then the result of bound_cg() is returned.

Initialize with array/list of eigenvalues or Intervals object.

static __new__(cls, evals)¶: Use BoundCG if all eigenvalues are non-negative.

eval_step(step)¶: Evaluate bound for given step.

get_step(tol)¶: Return step at which bound falls below tolerance.

exception krypy.utils.ConvergenceError(msg, solver)

Bases: Exception

Raised when a method did not converge.

The ConvergenceError holds a message describing the error and the attribute solver through which the last approximation and other relevant information can be retrieved.

class krypy.utils.Givens(x)¶

Bases: object

Compute Givens rotation for provided vector x.

Computes Givens rotation \(G=\begin{bmatrix}c&s\\-\overline{s}&c\end{bmatrix}\) such that \(Gx=\begin{bmatrix}r\\0\end{bmatrix}\).

apply(x)¶: Apply Givens rotation to vector x.

class krypy.utils.House(x)¶

Bases: object

Compute Householder transformation for given vector.

Initialize Householder transformation \(H\) such that \(Hx = \alpha \|x\|_2 e_1\) with \(|\alpha|=1\)

The algorithm is a combination of Algorithm 5.1.1 on page 236 and the treatment of the complex case in Section 5.1.13 on page 243 in Golub, Van Loan. Matrix computations. Fourth Edition. 2013.

apply(x)¶

Apply Householder transformation to vector x.

Applies the Householder transformation efficiently to the given vector.

matrix()¶

Build matrix representation of Householder transformation.

Builds the matrix representation \(H = I - \beta vv^*\).

Use with care! This routine may be helpful for testing purposes but should not be used in production codes for high dimensions since the resulting matrix is dense.

class krypy.utils.IdentityLinearOperator(shape)¶: Bases: krypy.utils.LinearOperator

class krypy.utils.LinearOperator(shape, dtype, dot=None, dot_adj=None)¶

Bases: object

Linear operator.

Is partly based on the LinearOperator from scipy (BSD License).

__add__(X)¶

__mul__(X)¶

__neg__()¶

__pow__(X)¶

__repr__()¶: Return repr(self).

__rmul__(X)¶

__sub__(X)¶

adj¶

dot(X)¶

dot_adj(X)¶

class krypy.utils.MatrixLinearOperator(A)¶

Bases: krypy.utils.LinearOperator

__repr__()¶: Return repr(self).

class krypy.utils.NormalizedRootsPolynomial(roots)¶

Bases: object

A polynomial with specified roots and p(0)=1.

Represents the polynomial

\[p(\lambda) = \prod_{i=1}^n \left(1-\frac{\lambda}{\theta_i}\right).\]

Parameters:	roots – array with roots \(\theta_1,\dots,\theta_n\) of the polynomial and `roots.shape==(n,)`.

__call__(points)¶

Evaluate polyonmial at given points.

Parameters:	points – a point \(x\) or array of points \(x_1,\dots,x_m\) with `points.shape==(m,)`.
Returns:	\(p(x)\) or array of shape `(m,)` with \(p(x_1),\dots,p(x_m)\).

minmax_candidates()¶

Get points where derivative is zero.

Useful for computing the extrema of the polynomial over an interval if the polynomial has real roots. In this case, the maximum is attained for one of the interval endpoints or a point from the result of this function that is contained in the interval.

class krypy.utils.Projection(X, Y=None, ip_B=None, orthogonalize=True, iterations=2)¶

Bases: object

Generic projection.

This class can represent any projection (orthogonal and oblique) on a N-dimensional Hilbert space. A projection is a linear operator \(P\) with \(P^2=P\). A projection is uniquely defined by its range \(\mathcal{V}:=\operatorname{range}(P)\) and its kernel \(\mathcal{W}:=\operatorname{ker}(P)\); this projection is called \(P_{\mathcal{V},\mathcal{W}}\).

Let X and Y be two full rank arrays with shape==(N,k) and let \(\mathcal{X}\oplus\mathcal{Y}^\perp=\mathbb{C}^N\) where \(\mathcal{X}:=\operatorname{colspan}(X)\) and \(\mathcal{Y}:=\operatorname{colspan}(Y)\). Then this class constructs the projection \(P_{\mathcal{X},\mathcal{Y}^\perp}\). The requirement \(\mathcal{X}\oplus\mathcal{Y}^\perp=\mathbb{C}^N\) is equivalent to \langle X,Y\rangle being nonsingular.

Parameters:

X – array with shape==(N,k) and \(\operatorname{rank}(X)=k\).
Y – (optional) None or array with shape==(N,k) and \(\operatorname{rank}(X)=k\). If Y is None then Y is set to X which means that the resulting projection is orthogonal.
ip_B – (optional) inner product, see inner(). None, a numpy.array or a LinearOperator is preferred due to the applicability of the proposed algorithms in [Ste11], see below.
orthogonalize – (optional) True orthogonalizes the bases provided in X and Y with respect to the inner product defined by ip_B. Defaults to True as the orthogonalization is suggested by the round-off error analysis in [Ste11].
iterations – (optional) number of applications of the projection. It was suggested in [Ste11] to use 2 iterations (default) in order to achieve high accuracy (“twice is enough” as in the orthogonal case).

This projection class makes use of the round-off error analysis of oblique projections in the work of Stewart [Ste11] and implements the algorithms that are considered as the most stable ones (e.g., the XQRY representation in [Ste11]).

apply(a, return_Ya=False)¶

Apply the projection to an array.

The computation is carried out without explicitly forming the matrix corresponding to the projection (which would be an array with shape==(N,N)).

See also _apply().

apply_adj(a)¶

apply_complement(a, return_Ya=False)¶

Apply the complementary projection to an array.

Parameters:	z – array with `shape==(N,m)`.
Returns:	\(P_{\mathcal{Y}^\perp,\mathcal{X}}z = z - P_{\mathcal{X},\mathcal{Y}^\perp} z\).

apply_complement_adj(a)¶

matrix()¶

Builds matrix representation of projection.

Builds the matrix representation \(P = X \langle Y,X\rangle^{-1} \langle Y, I_N\rangle\).

Use with care! This routine may be helpful for testing purposes but should not be used in production codes for high dimensions since the resulting matrix is dense.

operator()¶

Get a LinearOperator corresponding to apply().

Returns:	a LinearOperator that calls apply().

operator_complement()¶

Get a LinearOperator corresponding to apply_complement().

Returns:	a LinearOperator that calls apply_complement().

class krypy.utils.Timer¶

Bases: list

Measure execution time of multiple code blocks with with.

Example:

t = Timer()
with t:
    print('time me!')
print('don\'t time me!')
with t:
    print('time me, too!')
print(t)

Result:

time me!
don't time me!
time me, too!
[6.389617919921875e-05, 6.008148193359375e-05]

__enter__()¶

__exit__(a, b, c)¶

krypy.utils.angles(F, G, ip_B=None, compute_vectors=False)¶

Principal angles between two subspaces.

This algorithm is based on algorithm 6.2 in Knyazev, Argentati. Principal angles between subspaces in an A-based scalar product: algorithms and perturbation estimates. 2002. This algorithm can also handle small angles (in contrast to the naive cosine-based svd algorithm).

Parameters:

F – array with shape==(N,k).
G – array with shape==(N,l).
ip_B – (optional) angles are computed with respect to this inner product. See inner().
compute_vectors – (optional) if set to False then only the angles are returned (default). If set to True then also the principal vectors are returned.

Returns:

theta if compute_vectors==False
theta, U, V if compute_vectors==True

where

theta is the array with shape==(max(k,l),) containing the principal angles \(0\leq\theta_1\leq\ldots\leq\theta_{\max\{k,l\}}\leq \frac{\pi}{2}\).
U are the principal vectors from F with \(\langle U,U \rangle=I_k\).
V are the principal vectors from G with \(\langle V,V \rangle=I_l\).

The principal angles and vectors fulfill the relation \(\langle U,V \rangle = \begin{bmatrix} \cos(\Theta) & 0_{m,l-m} \\ 0_{k-m,m} & 0_{k-m,l-m} \end{bmatrix}\) where \(m=\min\{k,l\}\) and \(\cos(\Theta)=\operatorname{diag}(\cos(\theta_1),\ldots,\cos(\theta_m))\). Furthermore, \(\theta_{m+1}=\ldots=\theta_{\max\{k,l\}}=\frac{\pi}{2}\).

krypy.utils.arnoldi(*args, **kwargs)¶

krypy.utils.arnoldi_res(A, V, H, ip_B=None)¶

Measure Arnoldi residual.

Parameters:	A – a linear operator that can be used with scipy’s aslinearoperator with `shape==(N,N)`. V – Arnoldi basis matrix with `shape==(N,n)`. H – Hessenberg matrix: either \(\underline{H}_{n-1}\) with `shape==(n,n-1)` or \(H_n\) with `shape==(n,n)` (if the Arnoldi basis spans an A-invariant subspace). ip_B – (optional) the inner product to use, see `inner()`.
Returns:	either \(\\|AV_{n-1} - V_n \underline{H}_{n-1}\\|\) or \(\\|A V_n - V_n H_n\\|\) (in the invariant case).

krypy.utils.arnoldi_projected(H, P, k, ortho='mgs')¶

Compute (perturbed) Arnoldi relation for projected operator.

Assume that you have computed an Arnoldi relation

\[A V_n = V_{n+1} \underline{H}_n\]

where \(V_{n+1}\in\mathbb{C}^{N,n+1}\) has orthogonal columns (with respect to an inner product \(\langle\cdot,\cdot\rangle\)) and \(\underline{H}_n\in\mathbb{C}^{n+1,n}\) is an extended upper Hessenberg matrix.

For \(k<n\) you choose full rank matrices \(X\in\mathbb{C}^{n-1,k}\) and \(Y\in\mathbb{C}^{n,k}\) and define \(\tilde{X}:=A V_{n_1}X = V_n \underline{H}_{n-1} X\) and \(\tilde{Y}:=V_n Y\) such that \(\langle \tilde{Y}, \tilde{X} \rangle = Y^*\underline{H}_{n-1} X\) is invertible. Then the projections \(P\) and \(\tilde{P}\) characterized by

\(\tilde{P}x = x - \tilde{X} \langle \tilde{Y},\tilde{X} \rangle^{-1} \langle\tilde{Y},x\rangle\)
\(P = I - \underline{H}_{n-1}X (Y^*\underline{H}_{n-1}X)^{-1}Y^*\)

are well defined and \(\tilde{P}V_{n+1} = [V_n P, v_{n+1}]\) holds.

This method computes for \(i<n-k\) the Arnoldi relation

\[(\tilde{P}A + E_i) W_i = W_{i+1} \underline{G}_i\]

where \(W_{i+1}=V_n U_{i+1}\) has orthogonal columns with respect to \(\langle\cdot,\cdot\rangle\), \(\underline{G}_i\) is an extended upper Hessenberg matrix and \(E_i x = v_{n+1} F_i \langle W_i,x\rangle\) with \(F_i=[f_1,\ldots,f_i]\in\mathbb{C}^{1,i}\).

The perturbed Arnoldi relation can also be generated with the operator \(P_{V_n} \tilde{P} A\):

\[P_{V_n} \tilde{P} A W_i = W_{i+1} \underline{G}_i.\]

In a sense the perturbed Arnoldi relation is the best prediction for the behavior of the Krylov subspace \(K_i(\tilde{P}A,\tilde{P}v_1)\) that can be generated only with the data from \(K_{n+1}(A,v_1)\) and without carrying out further matrix-vector multiplications with A.

Parameters:

H – the extended upper Hessenberg matrix \(\underline{H}_n\) with shape==(n+1,n).
P – the projection \(P:\mathbb{C}^n\longrightarrow\mathbb{C}^n\) (has to be compatible with get_linearoperator()).
k – the dimension of the null space of P.

Returns:

U, G, F where

U is the coefficient matrix \(U_{i+1}\) with shape==(n,i+1),
G is the extended upper Hessenberg matrix \(\underline{G}_i\) with shape==(i+1,i),
F is the error matrix \(F_i\) with shape==(1,i).

krypy.utils.bound_perturbed_gmres(pseudo, p, epsilon, deltas)¶

Compute GMRES perturbation bound based on pseudospectrum

Computes the GMRES bound from [SifEM13].

krypy.utils.gap(lamda, sigma, mode='individual')¶

Compute spectral gap.

Useful for eigenvalue/eigenvector bounds. Computes the gap \(\delta\geq 0\) between two sets of real numbers lamda and sigma. The gap can be computed in several ways and may not exist, see the mode parameter.

Parameters:

lamda – a non-empty set \(\Lambda=\{\lambda_1,\ldots,\lambda_n\}\) given as a single real number or a list or numpy.array with real numbers.
sigma – a non-empty set \(\Sigma=\{\sigma_1,\ldots,\sigma_m\}\). See lamda.
mode –
(optional). Defines how the gap should be computed. May be one of
- 'individual' (default): \(\delta=\min_{\substack{i\in\{1,\ldots,n\}\\j\in\{1,\ldots,m\}}} |\lambda_i - \sigma_j|\). With this mode, the gap is always be defined.
- 'interval': determine the maximal \(\delta\) such that \(\Sigma\subset\mathbb{R}\setminus[\min_{\lambda\in\Lambda}\lambda-\delta,\max_{\lambda\in\Lambda}\lambda+\delta]\). If the gap does not exists, None is returned.

Returns:

\(\delta\) or None.

krypy.utils.get_linearoperator(shape, A, timer=None)¶: Enhances aslinearoperator if A is None.

krypy.utils.hegedus(A, b, x0, M=None, Ml=None, ip_B=None)¶

Rescale initial guess appropriately (Hegedüs trick).

The Hegedüs trick rescales the initial guess to \(\gamma_{\min} x_0\) such that

\[\|r_0\|_{M^{-1}} = \| M M_l (b - A \gamma_{\min} x_0) \|_{M^{-1}} = \min_{\gamma\in\mathbb{C}} \| M M_l (b - A \gamma x_0) \|_{M^{-1}} \leq \| M M_l b \|_{M^{-1}}.\]

This is achieved by \(\gamma_{\min} = \frac{\langle z, M M_l b \rangle_{M^{-1}}}{\|z\|_{M^{-1}}^2}\) for \(z=M M_l A x_0\) because then \(r_0=P_{z^\perp}b\). (Note that the right hand side of formula (5.8.16) in [LieS13] has to be complex conjugated.)

The parameters are the parameters you want to pass to gmres(), minres() or cg().

Returns:	the adapted initial guess with the above property.

krypy.utils.inner(X, Y, ip_B=None)¶

Euclidean and non-Euclidean inner product.

numpy.vdot only works for vectors and numpy.dot does not use the conjugate transpose.

Parameters:	X – numpy array with `shape==(N,m)` Y – numpy array with `shape==(N,n)` ip_B – (optional) May be one of the following `None`: Euclidean inner product. a self-adjoint and positive definite operator \(B\) (as `numpy.array` or `LinearOperator`). Then \(X^*B Y\) is returned. a callable which takes 2 arguments X and Y and returns \(\langle X,Y\rangle\).

Caution: a callable should only be used if necessary. The choice potentially has an impact on the round-off behavior, e.g. of projections.

Returns:	numpy array \(\langle X,Y\rangle\) with `shape==(m,n)`.

krypy.utils.ip_euclid(X, Y)¶

Euclidean inner product.

numpy.vdot only works for vectors and numpy.dot does not use the conjugate transpose.

Parameters:	X – numpy array with `shape==(N,m)` Y – numpy array with `shape==(N,n)`
Returns:	numpy array \(X^* Y\) with `shape==(m,n)`.

krypy.utils.norm(x, y=None, ip_B=None)¶

Compute norm (Euclidean and non-Euclidean).

Parameters:	x – a 2-dimensional `numpy.array`. y – a 2-dimensional `numpy.array`. ip_B – see `inner()`.

Compute \(\sqrt{\langle x,y\rangle}\) where the inner product is defined via ip_B.

krypy.utils.norm_MMlr(M, Ml, A, Mr, b, x0, yk, inner_product=<function ip_euclid>)¶

krypy.utils.norm_squared(x, Mx=None, inner_product=<function ip_euclid>)¶: Compute the norm^2 w.r.t. to a given scalar product.

krypy.utils.orthonormality(V, ip_B=None)¶

Measure orthonormality of given basis.

Parameters:	V – a matrix \(V=[v_1,\ldots,v_n]\) with `shape==(N,n)`. ip_B – (optional) the inner product to use, see `inner()`.
Returns:	\(\\| I_n - \langle V,V \rangle \\|_2\).

krypy.utils.qr(X, ip_B=None, reorthos=1)¶

QR factorization with customizable inner product.

Parameters:	X – array with `shape==(N,k)` ip_B – (optional) inner product, see `inner()`. reorthos – (optional) numer of reorthogonalizations. Defaults to 1 (i.e. 2 runs of modified Gram-Schmidt) which should be enough in most cases (TODO: add reference).
Returns:	Q, R where \(X=QR\) with \(\langle Q,Q \rangle=I_k\) and R upper triangular.

krypy.utils.ritz(H, V=None, hermitian=False, type='ritz')¶

Compute several kinds of Ritz pairs from an Arnoldi/Lanczos relation.

This function computes Ritz, harmonic Ritz or improved harmonic Ritz values and vectors with respect to the Krylov subspace \(K_n(A,v)\) from the extended Hessenberg matrix \(\underline{H}_n\) generated with n iterations the Arnoldi algorithm applied to A and v.

Parameters:

H – Hessenberg matrix from Arnoldi/Lanczos algorithm.
V –
(optional) Arnoldi/Lanczos vectors, \(V\in\mathbb{C}^{N,n+1}\). If provided, the Ritz vectors are also returned. The Arnoldi vectors have to form an orthonormal basis with respect to an inner product.

Caution: if you are using the Lanzcos or Gram-Schmidt Arnoldi algorithm without reorthogonalization, then the orthonormality of the basis is usually lost. For accurate results it is advisable to use the Householder Arnoldi (ortho='house') or modified Gram-Schmidt with reorthogonalization (ortho='dmgs').
hermitian – (optional) if set to True the matrix \(H_n\) must be Hermitian. A Hermitian matrix \(H_n\) allows for faster and often more accurate computation of Ritz pairs.
type –
(optional) type of Ritz pairs, may be one of 'ritz', 'harmonic' or 'harmonic_like'. Two choices of Ritz pairs fit in the following description:

Given two n-dimensional subspaces \(X,Y\subseteq \mathbb{C}^N\), find a basis \(z_1,\ldots,z_n\) of \(X\) and \(\theta_1,\ldots,\theta_n\in\mathbb{C}\) such that \(A z_i - \theta_i z_i \perp Y\) for all \(i\in\{1,\ldots,n\}\).

In this setting the choices are
- 'ritz': regular Ritz pairs, i.e. \(X=Y=K_n(A,v)\).
- 'harmonic': harmonic Ritz pairs, i.e.
  
  \(X=K_n(A,v)\) and \(Y=AK_n(A,v)\).
- 'harmonic_improved': the returned vectors U (and V, if requested) are the same as with type='harmonic'. The theta array contains the improved Ritz values \(\theta_i = u_i^* H_n u_i\), cf. section 2 in Morgan, Zeng. Harmonic Projection Methods for Large Non-symmetric Eigenvalue Problems. 1998. It can be shown that the residual norm of improved Ritz pairs is always less than or equal to the residual norm of the harmonic Ritz pairs. However, the improved Ritz pairs do not fit into the framework above since the orthogonality condition is lost.

Returns:

If V is not None then theta, U, resnorm, Z is returned.
If V is None then theta, U, resnorm is returned.

Where

theta are the Ritz values \([\theta_1,\ldots,\theta_n]\).
U are the coefficients of the Ritz vectors in the Arnoldi basis, i.e. \(z_i=Vu_i\) where \(u_i\) is the i-th column of U.
resnorm is a residual norm vector.
Z are the actual Ritz vectors, i.e. Z=dot(V,U).

krypy.utils.shape_vec(x)¶: Take a (n,) ndarray and return it as (n,1) ndarray.

krypy.utils.shape_vecs(*args)¶

Reshape all ndarrays with shape==(n,) to shape==(n,1).

Recognizes ndarrays and ignores all others.

krypy.utils.strakos(n, l_min=0.1, l_max=100, rho=0.9)¶

Return the Strakoš matrix.

See [Str92].

krypy.utils - Krylov Subspace Utilities¶

`krypy.utils` - Krylov Subspace Utilities¶