InterpROM

class sharpy.rom.utils.librom_interp.InterpROM(SS, VV=None, WWT=None, Vref=None, WTref=None, method_proj=None)[source]

State-space 1D interpolation class.

This class allows interpolating from a list of state-space models, SS.

State-space models are required to have the same number of inputs and outputs and need to have the same number of states.

For state-space interpolation, state-space models also need to be defined over the same set of generalised coordinates. If this is not the case, the projection matrices W and V used to produce the ROMs, ie

\[\mathbf{A}_{proj} = \mathbf{W}^\top \mathbf{A V}\]

where A is the full-states matrix, also need to be provided. This will allow projecting the state-space models onto a common set of generalised coordinates before interpoling.

For development purposes, the method currently creates a hard copy of the projected matrices into the self.AA, self.BB, self.CC lists

Inputs:

  • SS: list of state-space models (instances of libss.StateSpace class)

  • VV: list of V matrices used to produce SS. If None, it is assumed that ROMs are defined over the same basis

  • WWT: list of W^T matrices used to derive the ROMs.

  • Vref, WTref: reference subspaces for projection. Some methods neglect this input (e.g. panzer)

  • method_proj: method for projection of state-space models over common coordinates. Available options are:

    • leastsq: find left/right projectors using least squares approx. Suitable for all basis.

    • strongMAC: strong Modal Assurance Criterion [4] enforcement for general basis. See Ref. [3], Eq. (7)

    • strongMAC_BT: strong Modal Assurance Criterion [4] enforcement for basis obtained by Balanced Truncation. Equivalent to strongMAC

    • maraniello_BT: this is equivalent to strongMAC and strongMAC_BT but avoids inversions. However, performance are the same as other strongMAC approaches - it works only when basis map the same subspaces

    • weakMAC_right_orth: weak MAC enforcement [1,3] for state-space models with right orthonoraml basis, i.e. V.T V = I. This is like Ref. [1], but implemented only on one side.

    • weakMAC: implementation of weak MAC enforcement for a general system. The method orthonormalises the right basis (V) and then solves the orthogonal Procrustes problem.

    • for orthonormal basis (V.T V = I): !!! These methods are not tested !!!

      • panzer: produces a new reference point based on svd [2]

      • amsallem: project over Vref,WTref [1]

References:

[1] D. Amsallem and C. Farhat, An online method for interpolating linear parametric reduced-order models, SIAM J. Sci. Comput., 33 (2011), pp. 2169–2198.

[2] Panzer, J. Mohring, R. Eid, and B. Lohmann, Parametric model order reduction by matrix interpolation, at–Automatisierungstechnik, 58 (2010), pp. 475–484.

[3] Mahony, R., Sepulchre, R. & Absil, P. -a., 2004. Riemannian Geometry of Grassmann Manifolds with a View on Algorithmic Computation. Acta Applicandae Mathematicae, 80(2), pp.199–220.

[4] Geuss, M., Panzer, H. & Lohmann, B., 2013. On parametric model order reduction by matrix interpolation. 2013 European Control Conference (ECC), pp.3433–3438.

project()[source]

Project the state-space models onto the generalised coordinates of state-space model IImap