pymor.operators.numpy

Operators based on NumPy arrays.

This module provides the following NumPy-based Operators:

• NumpyMatrixOperator wraps a 2D NumPy array as an Operator.

• NumpyMatrixBasedOperator should be used as base class for all Operators which assemble into a NumpyMatrixOperator.

• NumpyGenericOperator wraps an arbitrary Python function between NumPy arrays as an Operator.

• NumpyCirculantOperator matrix-free operator of a circulant matrix.

• NumpyToeplitzOperator matrix-free operator of a Toeplitz matrix.

• NumpyHankelOperator matrix-free operator of a Hankel matrix.

Module Contents

class pymor.operators.numpy.NumpyCirculantOperator(c, solver=None, name=None)

Bases: pymor.operators.interface.Operator, pymor.core.cache.CacheableObject

Matrix-free representation of a (block) circulant matrix by a NumPy array.

A (block) circulant matrix is a special kind of (block) Toeplitz matrix which is (block) square and completely determined by its first (matrix-valued) column via

\[\begin{split}C = \begin{bmatrix} c_1 & c_n & c_{n-1} & \cdots & \cdots & c_2 \\ c_2 & c_1 & c_n & & & \vdots \\ c_3 & c_2 & \ddots & & & \vdots \\ \vdots & & & \ddots & c_n & c_{n-1} \\ \vdots & & & c_2 & c_1 & c_n \\ c_n & \cdots & \cdots & c_3 & c_2 & c_1 \end{bmatrix} \in \mathbb{C}^{n*p \times n*m},\end{split}\]

where the so-called circulant vector \(c \in \mathbb{C}^{\times n\times p\times m}\) denotes the first (matrix-valued) column of the matrix. The matrix \(C\) as seen above is not explicitly constructed, only c is stored. Efficient matrix-vector multiplications are realized with DFT in the class’ apply method. See [GVL13] Chapter 4.8.2. for details.

Parameters:

• c – The NumPy array of shape (n) or (n, p, m) that defines the circulant vector.

• name – Name of the operator.

Methods

apply	Apply the operator to a VectorArray.
apply_adjoint	Apply the adjoint operator.

cache_region = 'memory'

linear = True

apply(U, mu=None)

Apply the operator to a VectorArray.

Parameters:

• U – VectorArray of vectors to which the operator is applied.

• mu – The parameter values for which to evaluate the operator.

Returns:

|VectorArray| of the operator evaluations.

apply_adjoint(V, mu=None)

Apply the adjoint operator.

For any given linear Operator op, parameter values mu and VectorArrays U, V in the source resp. range we have:

op.apply_adjoint(V, mu).dot(U) == V.inner(op.apply(U, mu))

Thus, when op is represented by a matrix M, apply_adjoint is given by left-multiplication of (the complex conjugate of) M with V.

Parameters:

• V – VectorArray of vectors to which the adjoint operator is applied.

• mu – The parameter values for which to apply the adjoint operator.

Returns:

|VectorArray| of the adjoint operator evaluations.

class pymor.operators.numpy.NumpyGenericOperator(mapping, adjoint_mapping=None, dim_source=1, dim_range=1, linear=False, parameters={}, solver=None, name=None)

Bases: pymor.operators.interface.Operator

Wraps an arbitrary Python function between NumPy arrays as an Operator.

Parameters:

• mapping –

  The function to wrap. If parameters is None, the function is of the form mapping(U) and is expected to be vectorized. In particular:

  mapping(U).shape == (dim_range,) + U.shape[1:].
  

  If parameters is not None, the function has to have the signature mapping(U, mu).

• adjoint_mapping –

  The adjoint function to wrap. If parameters is None, the function is of the form adjoint_mapping(U) and is expected to be vectorized. In particular:

  adjoint_mapping(U).shape == (dim_source,) + U.shape[1:].
  

  If parameters is not None, the function has to have the signature adjoint_mapping(U, mu).

• dim_source – Dimension of the operator’s source.

• dim_range – Dimension of the operator’s range.

• linear – Set to True if the provided mapping and adjoint_mapping are linear.

• parameters – The Parameters the operator depends on.

• solver – The Solver for the operator.

• name – Name of the operator.

Methods

apply	Apply the operator to a VectorArray.
apply_adjoint	Apply the adjoint operator.

apply(U, mu=None)

Apply the operator to a VectorArray.

Parameters:

• U – VectorArray of vectors to which the operator is applied.

• mu – The parameter values for which to evaluate the operator.

Returns:

|VectorArray| of the operator evaluations.

apply_adjoint(V, mu=None)

Apply the adjoint operator.

For any given linear Operator op, parameter values mu and VectorArrays U, V in the source resp. range we have:

op.apply_adjoint(V, mu).dot(U) == V.inner(op.apply(U, mu))

Thus, when op is represented by a matrix M, apply_adjoint is given by left-multiplication of (the complex conjugate of) M with V.

Parameters:

• V – VectorArray of vectors to which the adjoint operator is applied.

• mu – The parameter values for which to apply the adjoint operator.

Returns:

|VectorArray| of the adjoint operator evaluations.

class pymor.operators.numpy.NumpyHankelOperator(c, r=None, solver=None, name=None)

Bases: pymor.operators.interface.Operator

Matrix-free representation of a finite dimensional Hankel operator by a NumPy array.

A finite dimensional Hankel operator can be represented by a (block) matrix with constant anti-diagonals (anti-diagonal blocks), i.e.:

\[\begin{split}H = \begin{bmatrix} c_1 & c_2 & c_3 & \cdots & \cdots & r_1 \\ c_2 & c_3 & & & & \vdots \\ c_3 & & & & & \vdots \\ \vdots & & & & & r_{n-2} \\ \vdots & & & & r_{n-2} & r_{n-1} \\ c_m & \cdots & \cdots & r_{n-2} & r_{n-1} & r_n \end{bmatrix} \in \mathbb{C}^{n*p \times k*m},\end{split}\]

where \(c\in\mathbb{C}^{s\times p\times m}\) and \(r\in\mathbb{C}^{k\times p\times m}\) denote the first (matrix-valued) column and last (matrix-valued) row of the (block) Hankel matrix, respectively. The matrix \(H\) as seen above is not explicitly constructed, only the arrays c and r are stored. Efficient matrix-vector multiplications are realized with DFT in the class’ apply method (see [MSKC21] Algorithm 3.1. for details).

Parameters:

• c – The NumPy array of shape (n) or (n, p, m) that defines the first column of the (block) Hankel matrix.

• r – The NumPy array of shape (k,) or (k, p, m) that defines the last row of the (block) Hankel matrix. If supplied, its first entry r[0] has to be equal to c[-1]. Defaults to None. If r is None, the behavior of scipy.linalg.hankel is mimicked which sets r to zero (except for the first entry).

• name – Name of the operator.

Methods

apply	Apply the operator to a VectorArray.
apply_adjoint	Apply the adjoint operator.

linear = True

apply(U, mu=None)

Apply the operator to a VectorArray.

Parameters:

• U – VectorArray of vectors to which the operator is applied.

• mu – The parameter values for which to evaluate the operator.

Returns:

|VectorArray| of the operator evaluations.

apply_adjoint(V, mu=None)

Apply the adjoint operator.

For any given linear Operator op, parameter values mu and VectorArrays U, V in the source resp. range we have:

op.apply_adjoint(V, mu).dot(U) == V.inner(op.apply(U, mu))

Thus, when op is represented by a matrix M, apply_adjoint is given by left-multiplication of (the complex conjugate of) M with V.

Parameters:

• V – VectorArray of vectors to which the adjoint operator is applied.

• mu – The parameter values for which to apply the adjoint operator.

Returns:

|VectorArray| of the adjoint operator evaluations.

class pymor.operators.numpy.NumpyMatrixBasedOperator

Bases: pymor.operators.interface.Operator

Base class for operators which assemble into a NumpyMatrixOperator.

sparse

True if the operator assembles into a sparse matrix, False if the operator assembles into a dense matrix, None if unknown.

Methods

apply	Apply the operator to a VectorArray.
apply_adjoint	Apply the adjoint operator.
as_range_array	Return a VectorArray representation of the operator in its range space.
as_source_array	Return a VectorArray representation of the operator in its source space.
assemble	Assemble the operator for given parameter values.
export_matrix	Save the matrix of the operator to a file.

linear = True

apply(U, mu=None)

Apply the operator to a VectorArray.

Parameters:

• U – VectorArray of vectors to which the operator is applied.

• mu – The parameter values for which to evaluate the operator.

Returns:

|VectorArray| of the operator evaluations.

apply_adjoint(V, mu=None)

Apply the adjoint operator.

For any given linear Operator op, parameter values mu and VectorArrays U, V in the source resp. range we have:

op.apply_adjoint(V, mu).dot(U) == V.inner(op.apply(U, mu))

Thus, when op is represented by a matrix M, apply_adjoint is given by left-multiplication of (the complex conjugate of) M with V.

Parameters:

• V – VectorArray of vectors to which the adjoint operator is applied.

• mu – The parameter values for which to apply the adjoint operator.

Returns:

|VectorArray| of the adjoint operator evaluations.

as_range_array(mu=None)

Return a VectorArray representation of the operator in its range space.

In the case of a linear operator with NumpyVectorSpace as source, this method returns for given parameter values mu a VectorArray V in the operator’s range, such that

V.lincomb(U.to_numpy()) == self.apply(U, mu)

for all VectorArrays U.

Parameters:

mu – The parameter values for which to return the VectorArray representation.

Returns:

V – The VectorArray defined above.

as_source_array(mu=None)

Return a VectorArray representation of the operator in its source space.

In the case of a linear operator with NumpyVectorSpace as range, this method returns for given parameter values mu a VectorArray V in the operator’s source, such that

self.range.make_array(V.inner(U)) == self.apply(U, mu)

for all VectorArrays U.

Parameters:

mu – The parameter values for which to return the VectorArray representation.

Returns:

V – The VectorArray defined above.

assemble(mu=None)

Assemble the operator for given parameter values.

The result of the method strongly depends on the given operator. For instance, a matrix-based operator will assemble its matrix, a LincombOperator will try to form the linear combination of its operators, whereas an arbitrary operator might simply return a FixedParameterOperator. The only assured property of the assembled operator is that it no longer depends on a Parameter.

If the operator has a Solver, the assembled Operator will be equipped with the same Solver.

Parameters:

mu – The parameter values for which to assemble the operator.

Returns:

Parameter-independent, assembled |Operator|.

export_matrix(filename, matrix_name=None, output_format='matlab', mu=None)

Save the matrix of the operator to a file.

Parameters:

• filename – Name of output file.

• matrix_name – The name, the output matrix is given. (Comment field is used in case of Matrix Market output_format.) If None, the Operator’s name is used.

• output_format – Output file format. Either matlab or matrixmarket.

• mu – The parameter values to assemble the to be exported matrix for.

class pymor.operators.numpy.NumpyMatrixOperator(matrix, solver=None, name=None, default_dense_solver=ScipyLUSolveSolver, default_sparse_solver=ScipySpSolveSolver)

Bases: NumpyMatrixBasedOperator

Wraps a 2D NumPy array or SciPy spmatrix as an Operator.

Note

In the case of a NumPy array, the apply_inverse method by default uses check_finite=True and check_cond=True. Setting them to False (e.g., via defaults) can significantly speed up the computation, especially for smaller matrices.

Parameters:

• matrix – The NumPy array or SciPy spmatrix which is to be wrapped.

• solver – The Solver for the operator.

• name – Name of the operator.

Methods

apply	Apply the operator to a VectorArray.
apply_adjoint	Apply the adjoint operator.
as_range_array	Return a VectorArray representation of the operator in its range space.
as_source_array	Return a VectorArray representation of the operator in its source space.
assemble	Assemble the operator for given parameter values.
from_file

apply(U, mu=None)

Apply the operator to a VectorArray.

Parameters:

• U – VectorArray of vectors to which the operator is applied.

• mu – The parameter values for which to evaluate the operator.

Returns:

|VectorArray| of the operator evaluations.

apply_adjoint(V, mu=None)

Apply the adjoint operator.

For any given linear Operator op, parameter values mu and VectorArrays U, V in the source resp. range we have:

op.apply_adjoint(V, mu).dot(U) == V.inner(op.apply(U, mu))

Thus, when op is represented by a matrix M, apply_adjoint is given by left-multiplication of (the complex conjugate of) M with V.

Parameters:

• V – VectorArray of vectors to which the adjoint operator is applied.

• mu – The parameter values for which to apply the adjoint operator.

Returns:

|VectorArray| of the adjoint operator evaluations.

as_range_array(mu=None)

Return a VectorArray representation of the operator in its range space.

In the case of a linear operator with NumpyVectorSpace as source, this method returns for given parameter values mu a VectorArray V in the operator’s range, such that

V.lincomb(U.to_numpy()) == self.apply(U, mu)

for all VectorArrays U.

Parameters:

mu – The parameter values for which to return the VectorArray representation.

Returns:

V – The VectorArray defined above.

as_source_array(mu=None)

Return a VectorArray representation of the operator in its source space.

In the case of a linear operator with NumpyVectorSpace as range, this method returns for given parameter values mu a VectorArray V in the operator’s source, such that

self.range.make_array(V.inner(U)) == self.apply(U, mu)

for all VectorArrays U.

Parameters:

mu – The parameter values for which to return the VectorArray representation.

Returns:

V – The VectorArray defined above.

assemble(mu=None)

Assemble the operator for given parameter values.

The result of the method strongly depends on the given operator. For instance, a matrix-based operator will assemble its matrix, a LincombOperator will try to form the linear combination of its operators, whereas an arbitrary operator might simply return a FixedParameterOperator. The only assured property of the assembled operator is that it no longer depends on a Parameter.

If the operator has a Solver, the assembled Operator will be equipped with the same Solver.

Parameters:

mu – The parameter values for which to assemble the operator.

Returns:

Parameter-independent, assembled |Operator|.

classmethod from_file(path, key=None, solver=None, name=None)

class pymor.operators.numpy.NumpyToeplitzOperator(c, r=None, solver=None, name=None)

Bases: pymor.operators.interface.Operator

Matrix-free representation of a finite dimensional Toeplitz matrix by a NumPy array.

A finite dimensional Toeplitz operator can be represented by a (block) matrix with constant diagonals (diagonal blocks), i.e.:

\[\begin{split}T = \begin{bmatrix} c_1 & r_2 & r_3 & \cdots & \cdots & r_n \\ c_2 & c_1 & r_2 & & & \vdots \\ c_3 & c_2 & \ddots & & & \vdots \\ \vdots & & & \ddots & r_2 & r_3 \\ \vdots & & & c_2 & c_1 & r_2 \\ c_m & \cdots & \cdots & c_3 & c_2 & c_1 \end{bmatrix} \in \mathbb{C}^{n*p \times k*m},\end{split}\]

where \(c\in\mathbb{C}^{n\times p\times m}\) and \(r\in\mathbb{C}^{k\times p\times m}\) denote the first (matrix-valued) column and first (matrix-valued) row of the (block) Toeplitz matrix, respectively. The matrix \(T\) as seen above is not explicitly constructed, only the arrays c and r are stored. The operator’s apply method takes advantage of the fact that any (block) Toeplitz matrix can be embedded in a larger (block) circulant matrix to leverage efficient matrix-vector multiplications with DFT.

Parameters:

• c – The NumPy array of shape either (n) or (n, p, m) that defines the first column of the (block) Toeplitz matrix.

• r – The NumPy array of shape (k,) or (k, p, m) that defines the first row of the Toeplitz matrix. If supplied, its first entry r[0] has to be equal to c[0]. Defaults to None. If r is None, the behavior of scipy.linalg.toeplitz is mimicked which sets r = c.conj() (except for the first entry).

• name – Name of the operator.

Methods

apply	Apply the operator to a VectorArray.
apply_adjoint	Apply the adjoint operator.

linear = True

apply(U, mu=None)

Apply the operator to a VectorArray.

Parameters:

• U – VectorArray of vectors to which the operator is applied.

• mu – The parameter values for which to evaluate the operator.

Returns:

|VectorArray| of the operator evaluations.

apply_adjoint(V, mu=None)

Apply the adjoint operator.

For any given linear Operator op, parameter values mu and VectorArrays U, V in the source resp. range we have:

op.apply_adjoint(V, mu).dot(U) == V.inner(op.apply(U, mu))

Thus, when op is represented by a matrix M, apply_adjoint is given by left-multiplication of (the complex conjugate of) M with V.

Parameters:

• V – VectorArray of vectors to which the adjoint operator is applied.

• mu – The parameter values for which to apply the adjoint operator.

Returns:

|VectorArray| of the adjoint operator evaluations.