pymor.reductors.h2

Reductors based on H2-norm.

Module Contents

class pymor.reductors.h2.GapIRKAReductor(fom, mu=None, solver_options=None)[source]

Bases: GenericIRKAReductor

Gap-IRKA reductor.

Parameters

fom

The full-order LTIModel to reduce.

mu

Parameter.

Methods

reduce

Reduce using gap-IRKA.

reduce(rom0_params, tol=0.0001, maxit=100, num_prev=1, conv_crit='sigma', projection='orth')[source]

Reduce using gap-IRKA.

See [BBG22] Algorithm 1.

Parameters

rom0_params

Can be:

  • order of the reduced model (a positive integer),

  • initial interpolation points (a 1D NumPy array),

  • dict with 'sigma', 'b', 'c' as keys mapping to initial interpolation points (a 1D NumPy array), right tangential directions (VectorArray from fom.input_space), and left tangential directions (VectorArray from fom.output_space), all of the same length (the order of the reduced model),

  • initial reduced-order model (LTIModel).

If the order of reduced model is given, initial interpolation data is generated randomly.

tol

Tolerance for the convergence criterion.

maxit

Maximum number of iterations.

num_prev

Number of previous iterations to compare the current iteration to. A larger number can avoid occasional cyclic behavior.

conv_crit

Convergence criterion:

  • 'sigma': relative change in interpolation points

  • 'htwogap': \(\mathcal{H}_2-gap\) distance of reduced-order models divided by \(\mathcal{L}_2\) norm of new reduced-order model

  • 'ltwo': relative \(\mathcal{L}_2\) distance of reduced-order models

projection

Projection method:

  • 'orth': projection matrix is orthogonalized with respect to the Euclidean inner product,

  • 'biorth': projection matrix is orthogonalized with respect to the E product.

Returns

rom

Reduced LTIModel model.

class pymor.reductors.h2.GenericIRKAReductor(fom, mu=None)[source]

Bases: pymor.core.base.BasicObject

Generic IRKA related reductor.

Parameters

fom

The full-order Model to reduce.

mu

Parameter values.

Methods

reconstruct

Reconstruct high-dimensional vector from reduced vector u.

reconstruct(u)[source]

Reconstruct high-dimensional vector from reduced vector u.

class pymor.reductors.h2.IRKAReductor(fom, mu=None)[source]

Bases: GenericIRKAReductor

Iterative Rational Krylov Algorithm reductor.

Parameters

fom

The full-order LTIModel to reduce.

mu

Parameter values.

Methods

reduce

Reduce using IRKA.

reduce(rom0_params, tol=0.0001, maxit=100, num_prev=1, force_sigma_in_rhp=False, projection='orth', conv_crit='sigma', compute_errors=False)[source]

Reduce using IRKA.

See [GAB08] (Algorithm 4.1) and [ABG10] (Algorithm 1).

Parameters

rom0_params

Can be:

  • order of the reduced model (a positive integer),

  • initial interpolation points (a 1D NumPy array),

  • dict with 'sigma', 'b', 'c' as keys mapping to initial interpolation points (a 1D NumPy array), right tangential directions (NumPy array of shape (len(sigma), fom.dim_input)), and left tangential directions (NumPy array of shape (len(sigma), fom.dim_input)),

  • initial reduced-order model (LTIModel).

If the order of reduced model is given, initial interpolation data is generated randomly.

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 biorthogonalized with respect to the E product

  • 'arnoldi': projection matrices are orthogonalized using the Arnoldi process (available only for SISO systems).

conv_crit

Convergence criterion:

  • 'sigma': relative change in interpolation points

  • 'h2': relative \(\mathcal{H}_2\) distance of reduced-order models

compute_errors

Should the relative \(\mathcal{H}_2\)-errors of intermediate reduced-order models be computed.

Warning

Computing \(\mathcal{H}_2\)-errors is expensive. Use this option only if necessary.

Returns

rom

Reduced LTIModel model.

class pymor.reductors.h2.OneSidedIRKAReductor(fom, version, mu=None)[source]

Bases: GenericIRKAReductor

One-Sided Iterative Rational Krylov Algorithm reductor.

Parameters

fom

The full-order LTIModel to reduce.

version

Version of the one-sided IRKA:

  • 'V': Galerkin projection using the input Krylov subspace,

  • 'W': Galerkin projection using the output Krylov subspace.

mu

Parameter values.

Methods

reduce

Reduce using one-sided IRKA.

reduce(rom0_params, tol=0.0001, maxit=100, num_prev=1, force_sigma_in_rhp=False, projection='orth', conv_crit='sigma', compute_errors=False)[source]

Reduce using one-sided IRKA.

Parameters

rom0_params

Can be:

  • order of the reduced model (a positive integer),

  • initial interpolation points (a 1D NumPy array),

  • dict with 'sigma', 'b', 'c' as keys mapping to initial interpolation points (a 1D NumPy array), right tangential directions (NumPy array of shape (len(sigma), fom.dim_input)), and left tangential directions (NumPy array of shape (len(sigma), fom.dim_input)),

  • initial reduced-order model (LTIModel).

If the order of reduced model is given, initial interpolation data is generated randomly.

tol

Tolerance for the largest change in interpolation points.

maxit

Maximum number of iterations.

num_prev

Number of previous iterations to compare the current iteration to. A larger number can avoid occasional cyclic behavior.

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 matrix is orthogonalized with respect to the Euclidean inner product,

  • 'Eorth': projection matrix is orthogonalized with respect to the E product.

conv_crit

Convergence criterion:

  • 'sigma': relative change in interpolation points,

  • 'h2': relative \(\mathcal{H}_2\) distance of reduced-order models.

compute_errors

Should the relative \(\mathcal{H}_2\)-errors of intermediate reduced-order models be computed.

Warning

Computing \(\mathcal{H}_2\)-errors is expensive. Use this option only if necessary.

Returns

rom

Reduced LTIModel model.

class pymor.reductors.h2.TFIRKAReductor(fom, mu=None)[source]

Bases: GenericIRKAReductor

Realization-independent IRKA reductor.

See [BG12].

Parameters

fom

TransferFunction or Model with a transfer_function attribute, with eval_tf and eval_dtf methods that should be defined at least over the open right half of the complex plane.

mu

Parameter values.

Methods

reconstruct

Reconstruct high-dimensional vector from reduced vector u.

reduce

Reduce using TF-IRKA.

reconstruct(u)[source]

Reconstruct high-dimensional vector from reduced vector u.

reduce(rom0_params, tol=0.0001, maxit=100, num_prev=1, force_sigma_in_rhp=False, conv_crit='sigma', compute_errors=False)[source]

Reduce using TF-IRKA.

Parameters

rom0_params

Can be:

  • order of the reduced model (a positive integer),

  • initial interpolation points (a 1D NumPy array),

  • dict with 'sigma', 'b', 'c' as keys mapping to initial interpolation points (a 1D NumPy array), right tangential directions (NumPy array of shape (len(sigma), fom.dim_input)), and left tangential directions (NumPy array of shape (len(sigma), fom.dim_input)),

  • initial reduced-order model (LTIModel).

If the order of reduced model is given, initial interpolation data is generated randomly.

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 poles in the left half-plane are reflected.

conv_crit

Convergence criterion:

  • 'sigma': relative change in interpolation points

  • 'h2': relative \(\mathcal{H}_2\) distance of reduced-order models

compute_errors

Should the relative \(\mathcal{H}_2\)-errors of intermediate reduced-order models be computed.

Warning

Computing \(\mathcal{H}_2\)-errors is expensive. Use this option only if necessary.

Returns

rom

Reduced LTIModel model.

class pymor.reductors.h2.TSIAReductor(fom, mu=None)[source]

Bases: GenericIRKAReductor

Two-Sided Iteration Algorithm reductor.

Parameters

fom

The full-order LTIModel to reduce.

mu

Parameter values.

Methods

reduce

Reduce using TSIA.

reduce(rom0_params, 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 [BKohlerS11].

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

rom0_params

Can be:

  • order of the reduced model (a positive integer),

  • initial interpolation points (a 1D NumPy array),

  • dict with 'sigma', 'b', 'c' as keys mapping to initial interpolation points (a 1D NumPy array), right tangential directions (NumPy array of shape (len(sigma), fom.dim_input)), and left tangential directions (NumPy array of shape (len(sigma), fom.dim_input)),

  • initial reduced-order model (LTIModel).

If the order of reduced model is given, initial interpolation data is generated randomly.

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 biorthogonalized with respect to the E product

conv_crit

Convergence criterion:

  • 'sigma': relative change in interpolation points

  • 'h2': relative \(\mathcal{H}_2\) distance of reduced-order models

compute_errors

Should the relative \(\mathcal{H}_2\)-errors of intermediate reduced-order models be computed.

Warning

Computing \(\mathcal{H}_2\)-errors is expensive. Use this option only if necessary.

Returns

rom

Reduced LTIModel.

class pymor.reductors.h2.VectorFittingReductor(s, Hs, weights=None, conjugate=True)[source]

Bases: pymor.core.base.BasicObject

Vector fitting reductor.

Only for single-input single-output (SISO) systems.

Parameters

s

Sampling points in the complex plane as a 1D NumPy array.

Hs

Transfer function values at the sampling points s as a 1D NumPy array. Alternatively, TransferFunction or Model with transfer_function attribute.

weights

Weights in the weighted least squares error as a 1D NumPy array. If not given, it is set to a vector of ones.

conjugate

Whether to include conjugated data.

Methods

reduce

Reduce using vector fitting.

reduce(r=None, lambdas=None, tol=0.0001, maxit=100)[source]

Reduce using vector fitting.

Based on [DrmavcGB15].

Parameters

r

Reduced order (if not given, it is len(lambdas)).

lambdas

Initial poles (if not given, it is set to -np.logspace(-1, 0, r)).

tol

Tolerance for the convergence criterion.

maxit

Maximum number of iterations.

Returns

rom

Reduced-order LTIModel.