pymor.algorithms.chol_qr

Module Contents

pymor.algorithms.chol_qr.shifted_chol_qr(A, product=None, maxiter=3, offset=0, orth_tol=None, check_finite=True, copy=True, product_norm=None)[source]

Orthonormalize a VectorArray using the shifted CholeskyQR algorithm.

This method computes a QR decomposition of a VectorArray via Cholesky factorizations of its Gramian matrix according to [FKN+20]. For ill-conditioned matrices, the Cholesky factorization will break down. In this case a diagonal shift will be applied to the Gramian.

  • shifted_chol_qr(A, maxiter=3, orth_tol=None) is equivalent to the shifted CholeskyQR3 algorithm (Algorithm 4.2 in [FKN+20]).

  • shifted_chol_qr(A, maxiter=np.inf, orth_tol=<some_number>) is equivalent to the shifted CholeskyQR algorithm (Algorithm 4.1 in [FKN+20]).

  • shifted_chol_qr(A, product=<some_product>, maxiter=3, orth_tol=None) is equivalent to the shifted CholeskyQR3 algorithm in an oblique inner product (Algorithm 5.1 in [FKN+20]).

Note

Note that the unshifted single iteration CholeskyQR algorithm is unstable. Setting maxiter=1 will perform the simple shifted CholeskyQR algorithm (Algorithm 2.1 in [FKN+20]) which is stable but leads to poor orthogonality in general, i.e.

\[\lVert Q^TQ - I \rVert_2 < 2,\]

(see Lemma 3.1 in [FKN+20]).

Parameters

A

The VectorArray which is to be orthonormalized.

product

The inner product Operator w.r.t. which to orthonormalize. If None, the Euclidean product is used.

maxiter

Maximum number of iterations. Defaults to 3.

offset

Assume that the first offset vectors are already orthonormal and apply the algorithm for the vectors starting at offset + 1.

orth_tol

If not None, check if the resulting VectorArray is really orthornormal and repeat the algorithm until the check passes or maxiter is reached.

check_finite

This argument is passed down to SciPy linalg functions. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.

copy

If True, create a copy of A instead of modifying A in-place.

product_norm

The spectral norm of the inner product. If None, it will be computed with eigs.

Returns

Q

The orthonormalized VectorArray.

R

The upper-triangular/trapezoidal matrix (if compute_R is True).