InterpConstantOperator#

class InterpConstantOperator(training_parameters=None, entries=None, InterpolatorClass: type = None, fromblock=False)[source]#

Parametric constant operator \(\Ophat_{\ell}(\qhat,\u;\bfmu) = \chat(\bfmu) \in \RR^r\) where the parametric dependence is handled with elementwise interpolation.

\[\chat(\bfmu) = \textrm{interpolate}( (\bfmu_0,\chat^{(0)}),\ldots,(\bfmu_{s-1},\chat^{(s-1)}); \bfmu)\]

Here,

  • \(\bfmu_0,\ldots,\bfmu_{s-1}\in\RR^p\) are the (fixed) training parameter values, and

  • \(\chat^{(i)} = \chat(\bfmu_i) \in \RR^r\) is the operator vector evaluated at the training parameter values.

See opinf.operators.ConstantOperator.

Parameters:
training_parameterslist of s scalars or (p,) 1D ndarrays

Parameter values for which the operator vector is known or will be inferred from data. If not provided in the constructor, use set_training_parameters() later.

entrieslist of s ndarrays, or None

Operator vectors corresponding to the training_parameters. If not provided in the constructor, use set_entries() later.

InterpolatorClasstype or None

Class for the elementwise interpolation. Must obey the syntax

>>> interpolator = InterpolatorClass(data_points, data_values)
>>> interpolator_evaluation = interpolator(new_data_point)

This can be, e.g., a class from scipy.interpolate. If None (default), use scipy.interpolate.CubicSpline for one-dimensional parameters and scipy.interpolate.LinearNDInterpolator otherwise. If not provided in the constructor, use set_interpolator() later.

fromblockbool

If True, interpret entries as a horizontal concatenation of arrays; if False (default), interpret entries as a list of arrays.

Properties:
entries#

Operator matrices corresponding to the training parameters values, i.e., entries[i] is the operator matrix corresponding to the parameter value training_parameters[i].

interpolator#

Interpolator object for evaluating the operator at specified parameter values.

parameter_dimension#

Dimension \(p\) of the parameter vector \(\bfmu\) that the operator matrix depends on.

shape#

Shape of the operator matrix when evaluated at a parameter value.

state_dimension#

Dimension \(r\) of the state \(\qhat\) that the operator acts on.

training_parameters#

Parameter values where the operator matrix is known or will be inferred from data.

Methods:

apply

Apply the operator to the given state and input at the specified parameter value, \(\Ophat_\ell(\qhat,\u;\bfmu)\).

copy

Return a copy of the operator.

datablock

Return the data matrix block corresponding to the operator.

evaluate

Evaluate the operator at the given parameter value.

galerkin

Project this operator to a low-dimensional linear space.

jacobian

Construct the state Jacobian of the operator, \(\ddqhat\Ophat_\ell(\qhat,\u;\bfmu)\).

load

Load a parametric operator from an HDF5 file.

operator_dimension

Number of columns in the concatenated operator matrix.

save

Save the operator to an HDF5 file.

set_entries

Set the operator matrices at the training parameter values.

set_interpolator

Construct the interpolator for the operator matrix.

set_training_parameters

Set the training parameter values.

verify

Verify dimension attributes and evaluate().