pymor.reductors package¶
Submodules¶
basic module¶
-
class
pymor.reductors.basic.
GenericPGReductor
(d, W, V, bases_are_biorthonormal, vector_ranged_operators=('initial_data', ), product=None)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Generic Petrov-Galerkin reductor.
Replaces each
Operator
of the givenDiscretization
with the projection onto the span of the given projection matrices.Parameters
- d
- The
Discretization
which is to be reduced. - W
VectorArray
containing the left projection matrix.- V
VectorArray
containing the right projection matrix.- bases_are_biorthonormal
- Indicate whether or not V and W are biorthonormal w.r.t.
product
. - vector_ranged_operators
- List of keys in
d.operators
for which the correspondingOperator
should be biorthogonally projected (i.e. operators which map to vectors in contrast to bilinear forms which map to functionals). - product
- Inner product for the projection of the
Operators
given byvector_ranged_operators
.
Methods
GenericPGReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
()[source]¶ Perform the Petrov-Galerkin projection.
Returns
The reduced
Discretization
.
-
class
pymor.reductors.basic.
GenericRBReductor
(d, RB=None, basis_is_orthonormal=None, vector_ranged_operators=('initial_data', ), product=None)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Generic reduced basis reductor.
Replaces each
Operator
of the givenDiscretization
with the Galerkin projection onto the span of the given reduced basis.Parameters
- d
- The
Discretization
which is to be reduced. - RB
VectorArray
containing the reduced basis on which to project.- basis_is_orthonormal
- If
RB
is specified, indicate whether or not the basis is orthonormal w.r.t.product
. - vector_ranged_operators
- List of keys in
d.operators
for which the correspondingOperator
should be orthogonally projected (i.e. operators which map to vectors in contrast to bilinear forms which map to functionals). - product
- Inner product for the orthonormalization of
RB
and the projection of theOperators
given byvector_ranged_operators
.
Methods
GenericRBReductor
extend_basis
,reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
extend_basis
(U, method='gram_schmidt', pod_modes=1, pod_orthonormalize=True, orthonormal=None, copy_U=True)[source]¶ Extend basis by new vectors.
Parameters
- U
VectorArray
containing the new basis vectors.- method
Basis extension method to use. The following methods are available:
trivial: Vectors in U
are appended to the basis. Duplicate vectors in the sense ofalmost_equal
are removed.gram_schmidt: New basis vectors are orthonormalized w.r.t. to the old basis using the gram_schmidt
algorithm.pod: Append the first POD modes of the defects of the projections of the vectors in U onto the existing basis (e.g. for use in POD-Greedy algorithm). Warning
In case of the
'gram_schmidt'
and'pod'
extension methods, the existing reduced basis is assumed to be orthonormal w.r.t. the given inner product.- pod_modes
- In case
method == 'pod'
, the number of POD modes that shall be appended to the basis. - pod_orthonormalize
- If
True
andmethod == 'pod'
, re-orthonormalize the new basis vectors obtained by the POD in order to improve numerical accuracy. - orthonormal
- If
method == 'trivial'
, set this toTrue
to indicate that the basis will remain orthonormal after extending. - copy_U
- If
copy_U
isFalse
, the new basis vectors might be removed fromU
.
Raises
- ExtensionError
- Raised when the selected extension method does not yield a basis of increased dimension.
-
reduce
(dim=None)[source]¶ Perform the reduced basis projection.
Parameters
- dim
- If specified, the desired reduced state dimension. Must not be larger than the current reduced basis dimension.
Returns
The reduced
Discretization
.
bt module¶
-
class
pymor.reductors.bt.
BRBTReductor
(d, gamma, solver_options=None)[source]¶ Bases:
pymor.reductors.bt.GenericBTReductor
Bounded Real (BR) Balanced Truncation reductor.
See [A05] (Section 7.5.3) and [OJ88].
Parameters
- d
- The system which is to be reduced.
- gamma
- Upper bound for the
-norm.
- solver_options
- The solver options to use to solve the positive Riccati equations.
Methods
BRBTReductor
error_bounds
,gramians
GenericBTReductor
reconstruct
,reduce
,sv_U_V
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.bt.
BTReductor
(d)[source]¶ Bases:
pymor.reductors.bt.GenericBTReductor
Standard (Lyapunov) Balanced Truncation reductor.
See Section 7.3 in [A05].
Parameters
- d
- The system which is to be reduced.
Methods
BTReductor
error_bounds
,gramians
GenericBTReductor
reconstruct
,reduce
,sv_U_V
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.bt.
GenericBTReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Generic Balanced Truncation reductor.
Parameters
- d
- The system which is to be reduced.
Methods
GenericBTReductor
error_bounds
,gramians
,reconstruct
,reduce
,sv_U_V
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r=None, tol=None, projection='bfsr')[source]¶ Generic Balanced Truncation.
Parameters
- r
- Order of the reduced model if
tol
isNone
, maximum order iftol
is specified. - tol
- Tolerance for the error bound if
r
isNone
. - projection
Projection method used:
'sr'
: square root method'bfsr'
: balancing-free square root method (default, since it avoids scaling by singular values and orthogonalizes the projection matrices, which might make it more accurate than the square root method)'biorth'
: like the balancing-free square root method, except it biorthogonalizes the projection matrices (usinggram_schmidt_biorth
)
Returns
- rd
- Reduced system.
-
class
pymor.reductors.bt.
LQGBTReductor
(d, solver_options=None)[source]¶ Bases:
pymor.reductors.bt.GenericBTReductor
Linear Quadratic Gaussian (LQG) Balanced Truncation reductor.
See Section 3 in [MG91].
Parameters
- d
- The system which is to be reduced.
- solver_options
- The solver options to use to solve the Riccati equations.
Methods
LQGBTReductor
error_bounds
,gramians
GenericBTReductor
reconstruct
,reduce
,sv_U_V
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
coercive module¶
-
class
pymor.reductors.coercive.
CoerciveRBEstimator
(residual, residual_range_dims, coercivity_estimator)[source]¶ Bases:
pymor.core.interfaces.ImmutableInterface
Instantiated by
CoerciveRBReductor
.Not to be used directly.
Methods
CoerciveRBEstimator
estimate
,restricted_to_subbasis
ImmutableInterface
generate_sid
,with_
,__setattr__
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
-
class
pymor.reductors.coercive.
CoerciveRBReductor
(d, RB=None, basis_is_orthonormal=None, vector_ranged_operators=('initial_data', ), product=None, coercivity_estimator=None)[source]¶ Bases:
pymor.reductors.basic.GenericRBReductor
Reduced Basis reductor for
StationaryDiscretizations
with coercive linear operator.The only addition to
GenericRBReductor
is an error estimator which evaluates the dual norm of the residual with respect to a given inner product. For the reduction of the residual we useResidualReductor
for improved numerical stability [BEOR14].Parameters
- d
- The
Discretization
which is to be reduced. - RB
VectorArray
containing the reduced basis on which to project.- basis_is_orthonormal
- If
RB
is specified, indicate whether or not the basis is orthonormal w.r.t.product
. - vector_ranged_operators
- List of keys in
d.operators
for which the correspondingOperator
should be orthogonally projected (i.e. operators which map to vectors in contrast to bilinear forms which map to functionals). - product
- Inner product for the orthonormalization of
RB
, the projection of theOperators
given byvector_ranged_operators
and for the computation of Riesz representatives of the residual. IfNone
, the Euclidean product is used. - coercivity_estimator
None
or aParameterFunctional
returning a lower bound for the coercivity constant of the given problem. Note that the computed error estimate is only guaranteed to be an upper bound for the error when an appropriate coercivity estimate is specified.
Methods
GenericRBReductor
extend_basis
,reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.coercive.
SimpleCoerciveRBEstimator
(estimator_matrix, coercivity_estimator)[source]¶ Bases:
pymor.core.interfaces.ImmutableInterface
Instantiated by
SimpleCoerciveRBReductor
.Not to be used directly.
Methods
SimpleCoerciveRBEstimator
estimate
,restricted_to_subbasis
ImmutableInterface
generate_sid
,with_
,__setattr__
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
-
class
pymor.reductors.coercive.
SimpleCoerciveRBReductor
(d, RB=None, basis_is_orthonormal=None, vector_ranged_operators=('initial_data', ), product=None, coercivity_estimator=None)[source]¶ Bases:
pymor.reductors.basic.GenericRBReductor
Reductor for linear
StationaryDiscretizations
with affinely decomposed operator and rhs.Note
The reductor
CoerciveRBReductor
can be used for arbitrary coerciveStationaryDiscretizations
and offers an improved error estimator with better numerical stability.The only addition is to
GenericRBReductor
is an error estimator, which evaluates the norm of the residual with respect to a given inner product.Parameters
- d
- The
Discretization
which is to be reduced. - RB
VectorArray
containing the reduced basis on which to project.- basis_is_orthonormal
- If
RB
is specified, indicate whether or not the basis is orthonormal w.r.t.product
. - vector_ranged_operators
- List of keys in
d.operators
for which the correspondingOperator
should be orthogonally projected (i.e. operators which map to vectors in contrast to bilinear forms which map to functionals). - product
- Inner product for the orthonormalization of
RB
, the projection of theOperators
given byvector_ranged_operators
and for the computation of Riesz representatives of the residual. IfNone
, the Euclidean product is used. - coercivity_estimator
None
or aParameterFunctional
returning a lower bound for the coercivity constant of the given problem. Note that the computed error estimate is only guaranteed to be an upper bound for the error when an appropriate coercivity estimate is specified.
Methods
GenericRBReductor
extend_basis
,reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
h2 module¶
-
class
pymor.reductors.h2.
IRKAReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Iterative Rational Krylov Algorithm reductor.
Parameters
- d
LTISystem
.
Methods
IRKAReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r, sigma=None, b=None, c=None, rd0=None, tol=0.0001, maxit=100, num_prev=1, force_sigma_in_rhp=False, projection='orth', use_arnoldi=False, conv_crit='sigma', compute_errors=False)[source]¶ Reduce using IRKA.
See [GAB08] (Algorithm 4.1) and [ABG10] (Algorithm 1).
Parameters
- r
- Order of the reduced order model.
- sigma
Initial interpolation points (closed under conjugation).
If
None
, interpolation points are log-spaced between 0.1 and 10. Ifsigma
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be a one-dimensional array-like of lengthr
.sigma
andrd0
cannot both be notNone
.- b
Initial right tangential directions.
If
None
, if is chosen as all ones. Ifb
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be aVectorArray
of lengthr
fromd.B.source
.b
andrd0
cannot both be notNone
.- c
Initial left tangential directions.
If
None
, if is chosen as all ones. Ifc
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be aVectorArray
of lengthr
fromd.C.range
.c
andrd0
cannot both be notNone
.- rd0
Initial reduced order model.
If
None
, thensigma
,b
, andc
are used. Otherwise, it needs to be anLTISystem
of orderr
and it is used to constructsigma
,b
, andc
.- tol
- Tolerance for the convergence criterion.
- maxit
- Maximum number of iterations.
- num_prev
- Number of previous iterations to compare the current iteration to. Larger number can avoid occasional cyclic behavior of IRKA.
- force_sigma_in_rhp
- If
False
, new interpolation are reflections of the current reduced order model’s poles. Otherwise, only poles in the left half-plane are reflected. - projection
Projection method:
'orth'
: projection matrices are orthogonalized with respect to the Euclidean inner product'biorth'
: projection matrices are biorthogolized with respect to the E product
- use_arnoldi
- Should the Arnoldi process be used for rational interpolation. Available only for SISO systems. Otherwise, it is ignored.
- conv_crit
Convergence criterion:
'sigma'
: relative change in interpolation points'h2'
: relativedistance of reduced-order models
- compute_errors
Should the relative
-errors of intermediate reduced order models be computed.
Warning
Computing
-errors is expensive. Use this option only if necessary.
Returns
- rd
- Reduced
LTISystem
model.
-
class
pymor.reductors.h2.
TF_IRKAReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Realization-independent IRKA reductor.
See [BG12].
Parameters
- d
- Discretization with
eval_tf
andeval_dtf
methods.
Methods
TF_IRKAReductor
reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r, sigma=None, b=None, c=None, rd0=None, tol=0.0001, maxit=100, num_prev=1, force_sigma_in_rhp=False, conv_crit='sigma')[source]¶ Reduce using TF-IRKA.
Parameters
- r
- Order of the reduced order model.
- sigma
Initial interpolation points (closed under conjugation).
If
None
, interpolation points are log-spaced between 0.1 and 10. Ifsigma
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be a one-dimensional array-like of lengthr
.sigma
andrd0
cannot both be notNone
.- b
Initial right tangential directions.
If
None
, if is chosen as all ones. Ifb
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be aNumPy array
of shape(m, r)
.b
andrd0
cannot both be notNone
.- c
Initial left tangential directions.
If
None
, if is chosen as all ones. Ifc
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be aNumPy array
of shape(p, r)
.c
andrd0
cannot both be notNone
.- rd0
Initial reduced order model.
If
None
, thensigma
,b
, andc
are used. Otherwise, it needs to be anLTISystem
of orderr
and it is used to constructsigma
,b
, andc
.- tol
- Tolerance for the convergence criterion.
- maxit
- Maximum number of iterations.
- num_prev
- Number of previous iterations to compare the current iteration to. Larger number can avoid occasional cyclic behavior of TF-IRKA.
- force_sigma_in_rhp
- If
False
, new interpolation are reflections of the current reduced order model’s poles. Otherwise, only the poles in the left half-plane are reflected. - conv_crit
Convergence criterion:
'sigma'
: relative change in interpolation points'h2'
: relativedistance of reduced-order models
Returns
- rd
- Reduced
LTISystem
model.
-
class
pymor.reductors.h2.
TSIAReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Two-Sided Iteration Algorithm reductor.
Parameters
- d
LTISystem
.
Methods
TSIAReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(rd0, tol=0.0001, maxit=100, num_prev=1, projection='orth', conv_crit='sigma', compute_errors=False)[source]¶ Reduce using TSIA.
See [XZ11] (Algorithm 1) and [BKS11].
In exact arithmetic, TSIA is equivalent to IRKA (under some assumptions on the poles of the reduced model). The main difference in implementation is that TSIA computes the Schur decomposition of the reduced matrices, while IRKA computes the eigenvalue decomposition. Therefore, TSIA might behave better for non-normal reduced matrices.
Parameters
- rd0
- Initial reduced order model.
- tol
- Tolerance for the convergence criterion.
- maxit
- Maximum number of iterations.
- num_prev
- Number of previous iterations to compare the current iteration to. Larger number can avoid occasional cyclic behavior of TSIA.
- projection
Projection method:
'orth'
: projection matrices are orthogonalized with respect to the Euclidean inner product'biorth'
: projection matrices are biorthogolized with respect to the E product
- conv_crit
Convergence criterion:
'sigma'
: relative change in interpolation points'h2'
: relativedistance of reduced-order models
- compute_errors
Should the relative
-errors of intermediate reduced order models be computed.
Warning
Computing
-errors is expensive. Use this option only if necessary.
Returns
- rd
- Reduced
LTISystem
.
-
pymor.reductors.h2.
_convergence_criterion
(data, conv_crit)[source]¶ Compute the convergence criterion for given data.
interpolation module¶
-
class
pymor.reductors.interpolation.
DelayBHIReductor
(d)[source]¶ Bases:
pymor.reductors.interpolation.GenericBHIReductor
Bitangential Hermite interpolation for delay systems.
Parameters
Methods
GenericBHIReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.interpolation.
GenericBHIReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Generic bitangential Hermite interpolation reductor.
This is a generic reductor for reducing any linear
InputStateOutputSystem
with the transfer function which can be written in the generalized coprime factorizationas in [BG09]. The interpolation here is limited to only up to the first derivative. Hence, interpolation points are assumed to be pairwise distinct.
Parameters
- d
- Discretization.
Methods
GenericBHIReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(sigma, b, c, projection='orth')[source]¶ Bitangential Hermite interpolation.
Parameters
- sigma
- Interpolation points (closed under conjugation), list of
length
r
. - b
- Right tangential directions,
VectorArray
of lengthr
fromself.d.input_space
. - c
- Left tangential directions,
VectorArray
of lengthr
fromself.d.output_space
. - projection
Projection method:
'orth'
: projection matrices are orthogonalized with respect to the Euclidean inner product'biorth'
: projection matrices are biorthogolized with respect to the E product
Returns
- rd
- Reduced discretization.
-
class
pymor.reductors.interpolation.
LTI_BHIReductor
(d)[source]¶ Bases:
pymor.reductors.interpolation.GenericBHIReductor
Bitangential Hermite interpolation for
LTISystems
.Parameters
- d
LTISystem
.
Methods
LTI_BHIReductor
reduce
,reduce_arnoldi
GenericBHIReductor
reconstruct
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(sigma, b, c, projection='orth', use_arnoldi=False)[source]¶ Bitangential Hermite interpolation.
Parameters
- sigma
- Interpolation points (closed under conjugation), list of
length
r
. - b
- Right tangential directions,
VectorArray
of lengthr
fromself.d.input_space
. - c
- Left tangential directions,
VectorArray
of lengthr
fromself.d.output_space
. - projection
Projection method:
'orth'
: projection matrices are orthogonalized with respect to the Euclidean inner product'biorth'
: projection matrices are biorthogolized with respect to the E product
- use_arnoldi
- Should the Arnoldi process be used for rational interpolation. Available only for SISO systems. Otherwise, it is ignored.
Returns
- rd
- Reduced discretization.
-
reduce_arnoldi
(sigma, b, c)[source]¶ Bitangential Hermite interpolation for SISO
LTISystems
.Parameters
- sigma
- Interpolation points (closed under conjugation), list of
length
r
. - b
- Right tangential directions,
VectorArray
of lengthr
fromself.d.B.source
. - c
- Left tangential directions,
VectorArray
of lengthr
fromself.d.C.range
.
Returns
- rd
- Reduced
LTISystem
model.
-
class
pymor.reductors.interpolation.
SO_BHIReductor
(d)[source]¶ Bases:
pymor.reductors.interpolation.GenericBHIReductor
Bitangential Hermite interpolation for second-order systems.
Parameters
Methods
GenericBHIReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.interpolation.
TFInterpReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Loewner bitangential Hermite interpolation reductor.
See [BG12].
Parameters
- d
- Discretization with
eval_tf
andeval_dtf
methods.
Methods
TFInterpReductor
reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(sigma, b, c)[source]¶ Realization-independent tangential Hermite interpolation.
Parameters
- sigma
- Interpolation points (closed under conjugation), list of
length
r
. - b
- Right tangential directions,
NumPy array
of shape(d.m, r)
. - c
- Left tangential directions,
NumPy array
of shape(d.p, r)
.
Returns
- lti
LTISystem
interpolating the transfer function ofd
.
parabolic module¶
-
class
pymor.reductors.parabolic.
ParabolicRBEstimator
(residual, residual_range_dims, initial_residual, initial_residual_range_dims, coercivity_estimator)[source]¶ Bases:
pymor.core.interfaces.ImmutableInterface
Instantiated by
ParabolicRBReductor
.Not to be used directly.
Methods
ParabolicRBEstimator
estimate
,restricted_to_subbasis
ImmutableInterface
generate_sid
,with_
,__setattr__
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
-
class
pymor.reductors.parabolic.
ParabolicRBReductor
(d, RB=None, basis_is_orthonormal=None, product=None, coercivity_estimator=None)[source]¶ Bases:
pymor.reductors.basic.GenericRBReductor
Reduced Basis Reductor for parabolic equations.
This reductor uses
GenericRBReductor
for the actual RB-projection. The only addition is the assembly of an error estimator which bounds the discrete l2-in time / energy-in space error similar to [GP05], [HO08] as follows:Here,
denotes the norm induced by the problem’s mass matrix (e.g. the L^2-norm) and
is an arbitrary energy norm w.r.t. which the space operator
is coercive, and
is a lower bound for its coercivity constant. Finally,
denotes the implicit Euler timestepping residual for the (fixed) time step size
,
where
denotes the mass operator and
the source term. The dual norm of the residual is computed using the numerically stable projection from [BEOR14].
Parameters
- d
- The
InstationaryDiscretization
which is to be reduced. - RB
VectorArray
containing the reduced basis on which to project.- basis_is_orthonormal
- Indicate whether or not the basis is orthonormal w.r.t.
product
. - product
- The energy inner product
Operator
w.r.t. which the reduction error is estimated andRB
is orthonormalized. - coercivity_estimator
None
or aParameterFunctional
returning a lower boundfor the coercivity constant of
d.operator
w.r.t.product
.
Methods
GenericRBReductor
extend_basis
,reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
residual module¶
-
class
pymor.reductors.residual.
ImplicitEulerResidualOperator
(operator, mass, rhs, dt, name=None)[source]¶ Bases:
pymor.operators.basic.OperatorBase
Instantiated by
ImplicitEulerResidualReductor
.Methods
Attributes
-
apply
(U, U_old, mu=None)[source]¶ Apply the operator to a
VectorArray
.Parameters
- U
VectorArray
of vectors to which the operator is applied.- mu
- The
Parameter
for which to evaluate the operator.
Returns
VectorArray
of the operator evaluations.
-
-
class
pymor.reductors.residual.
ImplicitEulerResidualReductor
(RB, operator, mass, dt, rhs=None, product=None)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Reduced basis residual reductor with mass operator for implicit Euler timestepping.
Given an operator, mass and a functional, the concatenation of residual operator with the Riesz isomorphism is given by:
riesz_residual.apply(U, U_old, mu) == product.apply_inverse(operator.apply(U, mu) + 1/dt*mass.apply(U, mu) - 1/dt*mass.apply(U_old, mu) - rhs.as_vector(mu))
This reductor determines a low-dimensional subspace of the image of a reduced basis space under
riesz_residual
usingestimate_image_hierarchical
, computes an orthonormal basisresidual_range
of this range space and then returns the Petrov-Galerkin projectionprojected_riesz_residual == riesz_residual.projected(range_basis=residual_range, source_basis=RB)
of the
riesz_residual
operator. Given reduced basis coefficient vectorsu
andu_old
, the dual norm of the residual can then be computed asprojected_riesz_residual.apply(u, u_old, mu).l2_norm()
Moreover, a
reconstruct
method is provided such thatresidual_reductor.reconstruct(projected_riesz_residual.apply(u, u_old, mu)) == riesz_residual.apply(RB.lincomb(u), RB.lincomb(u_old), mu)
Parameters
- operator
- See definition of
riesz_residual
. - mass
- The mass operator. See definition of
riesz_residual
. - dt
- The time step size. See definition of
riesz_residual
. - rhs
- See definition of
riesz_residual
. IfNone
, zero right-hand side is assumed. - RB
VectorArray
containing a basis of the reduced space onto which to project.- product
- Inner product
Operator
w.r.t. which to compute the Riesz representatives.
Methods
ImplicitEulerResidualReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.residual.
NonProjectedImplicitEulerResidualOperator
(operator, mass, rhs, dt, product)[source]¶ Bases:
pymor.reductors.residual.ImplicitEulerResidualOperator
Instantiated by
ImplicitEulerResidualReductor
.Not to be used directly.
Methods
Attributes
-
apply
(U, U_old, mu=None)[source]¶ Apply the operator to a
VectorArray
.Parameters
- U
VectorArray
of vectors to which the operator is applied.- mu
- The
Parameter
for which to evaluate the operator.
Returns
VectorArray
of the operator evaluations.
-
-
class
pymor.reductors.residual.
NonProjectedResidualOperator
(operator, rhs, riesz_representatives, product)[source]¶ Bases:
pymor.reductors.residual.ResidualOperator
Instantiated by
ResidualReductor
.Not to be used directly.
Methods
Attributes
-
apply
(U, mu=None)[source]¶ Apply the operator to a
VectorArray
.Parameters
- U
VectorArray
of vectors to which the operator is applied.- mu
- The
Parameter
for which to evaluate the operator.
Returns
VectorArray
of the operator evaluations.
-
-
class
pymor.reductors.residual.
ResidualOperator
(operator, rhs, name=None)[source]¶ Bases:
pymor.operators.basic.OperatorBase
Instantiated by
ResidualReductor
.Methods
Attributes
-
apply
(U, mu=None)[source]¶ Apply the operator to a
VectorArray
.Parameters
- U
VectorArray
of vectors to which the operator is applied.- mu
- The
Parameter
for which to evaluate the operator.
Returns
VectorArray
of the operator evaluations.
-
-
class
pymor.reductors.residual.
ResidualReductor
(RB, operator, rhs=None, product=None, riesz_representatives=False)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Generic reduced basis residual reductor.
Given an operator and a right-hand side, the residual is given by:
residual.apply(U, mu) == operator.apply(U, mu) - rhs.as_range_array(mu)
When operator maps to functionals instead of vectors, we are interested in the Riesz representative of the residual:
residual.apply(U, mu) == product.apply_inverse(operator.apply(U, mu) - rhs.as_range_array(mu))
Given a basis
RB
of a subspace of the source space ofoperator
, this reductor usesestimate_image_hierarchical
to determine a low-dimensional subspace containing the image of the subspace underresidual
(resp.riesz_residual
), computes an orthonormal basisresidual_range
for this range space and then returns the Petrov-Galerkin projectionprojected_residual == project(residual, range_basis=residual_range, source_basis=RB)
of the residual operator. Given a reduced basis coefficient vector
u
, w.r.t.RB
, the (dual) norm of the residual can then be computed asprojected_residual.apply(u, mu).l2_norm()
Moreover, a
reconstruct
method is provided such thatresidual_reductor.reconstruct(projected_residual.apply(u, mu)) == residual.apply(RB.lincomb(u), mu)
Parameters
- RB
VectorArray
containing a basis of the reduced space onto which to project.- operator
- See definition of
residual
. - rhs
- See definition of
residual
. IfNone
, zero right-hand side is assumed. - product
- Inner product
Operator
w.r.t. which to orthonormalize and w.r.t. which to compute the Riesz representatives in caseoperator
maps to functionals. - riesz_representatives
- If
True
compute the Riesz representative of the residual.
Methods
ResidualReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
sobt module¶
-
class
pymor.reductors.sobt.
GenericSOBTpvReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Generic Second-Order Balanced Truncation position/velocity reductor.
See [RS08].
Parameters
- d
- The system which is to be reduced.
Methods
GenericSOBTpvReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r, projection='bfsr')[source]¶ Reduce using GenericSOBTpv.
Parameters
- r
- Order of the reduced model.
- projection
Projection method used:
'sr'
: square root method'bfsr'
: balancing-free square root method (default, since it avoids scaling by singular values and orthogonalizes the projection matrices, which might make it more accurate than the square root method)'biorth'
: like the balancing-free square root method, except it biorthogonalizes the projection matrices
Returns
- rd
- Reduced system.
-
class
pymor.reductors.sobt.
SOBTReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Second-Order Balanced Truncation reductor.
See [CLVV06].
Parameters
- d
- The system which is to be reduced.
Methods
SOBTReductor
reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r, projection='bfsr')[source]¶ Reduce using SOBT.
Parameters
- r
- Order of the reduced model.
- projection
Projection method used:
'sr'
: square root method'bfsr'
: balancing-free square root method (default, since it avoids scaling by singular values and orthogonalizes the projection matrices, which might make it more accurate than the square root method)'biorth'
: like the balancing-free square root method, except it biorthogonalizes the projection matrices
Returns
- rd
- Reduced system.
-
class
pymor.reductors.sobt.
SOBTfvReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
Free-velocity Second-Order Balanced Truncation reductor.
See [MS96].
Parameters
- d
- The system which is to be reduced.
Methods
SOBTfvReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r, projection='bfsr')[source]¶ Reduce using SOBTfv.
Parameters
- r
- Order of the reduced model.
- projection
Projection method used:
'sr'
: square root method'bfsr'
: balancing-free square root method (default, since it avoids scaling by singular values and orthogonalizes the projection matrices, which might make it more accurate than the square root method)'biorth'
: like the balancing-free square root method, except it biorthogonalizes the projection matrices
Returns
- rd
- Reduced system.
-
class
pymor.reductors.sobt.
SOBTpReductor
(d)[source]¶ Bases:
pymor.reductors.sobt.GenericSOBTpvReductor
Second-Order Balanced Truncation position reductor.
See [RS08].
Parameters
- d
- The system which is to be reduced.
Methods
GenericSOBTpvReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.sobt.
SOBTpvReductor
(d)[source]¶ Bases:
pymor.reductors.sobt.GenericSOBTpvReductor
Second-Order Balanced Truncation position-velocity reductor.
See [RS08].
Parameters
- d
- The system which is to be reduced.
Methods
GenericSOBTpvReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.sobt.
SOBTvReductor
(d)[source]¶ Bases:
pymor.reductors.sobt.GenericSOBTpvReductor
Second-Order Balanced Truncation velocity reductor.
See [RS08].
Parameters
- d
- The system which is to be reduced.
Methods
GenericSOBTpvReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
class
pymor.reductors.sobt.
SOBTvpReductor
(d)[source]¶ Bases:
pymor.reductors.sobt.GenericSOBTpvReductor
Second-Order Balanced Truncation velocity-position reductor.
See [RS08].
Parameters
- d
- The system which is to be reduced.
Methods
GenericSOBTpvReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
sor_irka module¶
-
class
pymor.reductors.sor_irka.
SOR_IRKAReductor
(d)[source]¶ Bases:
pymor.core.interfaces.BasicInterface
SOR-IRKA reductor.
Parameters
- d
- SecondOrderSystem.
Methods
SOR_IRKAReductor
reconstruct
,reduce
BasicInterface
disable_logging
,enable_logging
,has_interface_name
,implementor_names
,implementors
Attributes
BasicInterface
logger
,logging_disabled
,name
,uid
-
reduce
(r, sigma=None, b=None, c=None, rd0=None, tol=0.0001, maxit=100, num_prev=1, force_sigma_in_rhp=False, projection='orth', use_arnoldi=False, conv_crit='sigma', compute_errors=False, irka_options=None)[source]¶ Reduce using SOR-IRKA.
It uses IRKA as the intermediate reductor, to reduce from 2r to r poles. See Section 5.3.2 in [W12].
Parameters
- r
- Order of the reduced order model.
- sigma
Initial interpolation points (closed under conjugation).
If
None
, interpolation points are log-spaced between 0.1 and 10. Ifsigma
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be a one-dimensional array-like of lengthr
.sigma
andrd0
cannot both be notNone
.- b
Initial right tangential directions.
If
None
, if is chosen as all ones. Ifb
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be aVectorArray
of lengthr
fromd.B.source
.b
andrd0
cannot both be notNone
.- c
Initial left tangential directions.
If
None
, if is chosen as all ones. Ifc
is anint
, it is used as a seed to generate it randomly. Otherwise, it needs to be aVectorArray
of lengthr
fromd.Cp.range
.c
andrd0
cannot both be notNone
.- rd0
Initial reduced order model.
If
None
, thensigma
,b
, andc
are used. Otherwise, it needs to be anLTISystem
of orderr
and it is used to constructsigma
,b
, andc
.- tol
- Tolerance for the convergence criterion.
- maxit
- Maximum number of iterations.
- num_prev
- Number of previous iterations to compare the current iteration to. Larger number can avoid occasional cyclic behavior of IRKA.
- force_sigma_in_rhp
- If
False
, new interpolation are reflections of the current reduced order model’s poles. Otherwise, only the poles in the left half-plane are reflected. - projection
Projection method:
'orth'
: projection matrices are orthogonalized with respect to the Euclidean inner product'biorth'
: projection matrices are biorthogolized with respect to the E product
- conv_crit
Convergence criterion:
'sigma'
: relative change in interpolation points'h2'
: relativedistance of reduced-order models
- compute_errors
Should the relative
-errors of intermediate reduced order models be computed.
Warning
Computing
-errors is expensive. Use this option only if necessary.
- irka_options
- Dict of options for IRKAReductor.reduce.
Returns
- rd
- Reduced
LTISystem
model.