pymor.algorithms.rand_la

Module Contents

class pymor.algorithms.rand_la.RandomizedRangeFinder(A, range_product=None, source_product=None, A_adj=None, power_iterations=0, block_size=None, failure_tolerance=1e-15, num_testvecs=20, lambda_min=None, iscomplex=False, qr_method='gram_schmidt', error_estimator='bs18')[source]

Bases: pymor.core.base.BasicObject

Adaptive randomized range approximation of A.

Algorithm 1 in [BS18] is used as the default error estimator.

Given the Operator A, the return value of this method is the VectorArray B with the property

\[\Vert A - P_{span(B)} A \Vert \leq tol\]

with a failure probability smaller than failure_tolerance, where the norm denotes the operator norm. The inner product of the range of A is given by range_product and the inner product of the source of A is given by source_product.

A faster alternative is the leave-one-out error estimator from [ET24] which can be activated by setting error_estimator='loo'. The leave-one-out error estimator estimates the Frobenius norm of the approximation error. Note that while it can be faster and more accurate in practice, this option does only support Euclidian inner products and there are no results on its failure probability so far.

Parameters

A

The Operator A.

source_product

Inner product Operator of the source of A.

range_product

Inner product Operator of the range of A.

power_iterations

Number of power iterations.

A_adj

Adjoint Operator to use for power iterations. If None the adjoint is computed using A, source_product and range_product. Set to A for a self for a known self-adjoint operator.

failure_tolerance

Maximum failure probability. Only needed for error_estimator=’bs18’.

num_testvecs

Number of test vectors. Only needed for error_estimator=’bs18’.

lambda_min

The smallest eigenvalue of source_product. Only needed for error_estimator=’bs18’. If None, the smallest eigenvalue is computed using scipy.

block_size

Number of basis vectors to add per iteration.

iscomplex

If True, the random vectors are chosen complex.

qr_method

QR method used for orthogonalization. Either ‘gram_schmidt’ (default) or ‘shifted_chol_qr’.

error_estimator

The error estimator used. Either ‘bs18’ (default) or ‘loo’.

Methods

find_range

Find the range of A.

find_range(basis_size=None, tol=None)[source]

Find the range of A.

Parameters

basis_size

Maximum dimension of range approximation.

tol

Error tolerance for the algorithm.

Returns

VectorArray which contains the basis, whose span approximates the range of A.

pymor.algorithms.rand_la.randomized_ghep(A, E=None, n=6, power_iterations=0, oversampling=20, single_pass=False, return_evecs=False)[source]

Approximates a few eigenvalues of a Hermitian linear Operator with randomized methods.

Approximates modes eigenvalues w with corresponding eigenvectors v which solve the eigenvalue problem

\[A v_i = w_i v_i\]

or the generalized eigenvalue problem

\[A v_i = w_i E v_i\]

if E is not None.

This method is an implementation of algorithm 6 and 7 in [SLK16].

Parameters

A

The Hermitian linear Operator for which the eigenvalues are to be computed.

E

The Hermitian Operator which defines the generalized eigenvalue problem.

n

The number of eigenvalues and eigenvectors which are to be computed. Defaults to 6.

power_iterations

The number of power iterations to increase the relative weight of the larger singular values. Ignored when single_pass is True.

oversampling

The number of samples that are drawn in addition to the desired basis size in the randomized range approximation process.

single_pass

If True, computes the GHEP where only one set of matvecs Ax is required, but at the expense of lower numerical accuracy. If False, the methods performs two sets of matvecs Ax (default).

return_evecs

If True, the eigenvectors are computed and returned. Defaults to False.

Returns

w

A 1D NumPy array which contains the computed eigenvalues.

V

A VectorArray which contains the computed eigenvectors.

pymor.algorithms.rand_la.randomized_svd(A, n, source_product=None, range_product=None, power_iterations=0, oversampling=20)[source]

Randomized SVD of an Operator based on [SHvBW21].

Viewing the Operator \(A\) as an \(m\) by \(n\) matrix, this methods computes and returns the randomized generalized singular value decomposition of A according [SHvBW21]:

\[A = U \Sigma V^{-1},\]

where the inner product on the range \(\mathbb{R}^m\) is given by

\[(x, y)_R = x^TRy\]

and the inner product on the source \(\mathbb{R}^n\) is given by

\[(x, y)_S = x^TSy.\]

Note that \(U\) is \(R\)-orthogonal, i.e. \(U^TRU=I\) and \(V\) is \(S\)-orthogonal, i.e. \(V^TSV=I\). In particular, V^{-1}=V^TS.

Parameters

A

The Operator for which the randomized SVD is to be computed.

n

The number of eigenvalues and eigenvectors which are to be computed.

source_product

Source product Operator \(S\) w.r.t. which the randomized SVD is computed.

range_product

Range product Operator \(R\) w.r.t. which the randomized SVD is computed.

power_iterations

The number of power iterations to increase the relative weight of the larger singular values.

oversampling

The number of samples that are drawn in addition to the desired basis size in the randomized range approximation process.

Returns

U

VectorArray containing the approximated left singular vectors.

s

One-dimensional NumPy array of the approximated singular values.

V

VectorArray containing the approximated right singular vectors.