Complex and large numerical models require a substantial computational effort and large amounts of memory. Model order reduction methods reduce the size of the numerical models. The reduced models approximate the full system's response characteristics with less computational effort and memory requirements. Good reduced models are able to replace the original model.

However, the reduced model is not valid anymore when parameters of the original model such as material properties, the geometry or boundary conditions are changed. Applications, for which repeated model evaluations are required over a large range of parameter values are for example

  •     optimization
  •     uncertainty quantification
  •     control

  

In these applications, one is interested in having a reduced model that allows a rapid and accurate simulation over the range of the parameters in the so-called online phase. In order to obtain this reduced model, one is willing to forgo a large up-front cost, the so-called offline cost. 

A straightforward approach for these applications would be to generate a reduced order model (ROM) for each queried point in the parameter space. However, since generating a ROM is a costly task that requires assembling the full order model and performing computations with it, this is in general infeasible. Therefore, methods were sought to generate parametric ROMs, which gave rise to the field of parametric Model Order Reduction (pMOR).

Theory

We describe the dynamic system using a differential equation of second order

\begin{align} \mathbf{M}(\mathbf{p})\ddot{\mathbf{x}}\left(t, \mathbf{p}\right) + \mathbf{C}(\mathbf{p})\dot{\mathbf{x}}\left(t, \mathbf{p}\right)+\mathbf{K}(\mathbf{p})\mathbf{x}\left(t, \mathbf{p}\right) &= \mathbf{f}(t, \mathbf{p}) \\ y(t, \mathbf{p}) &= \mathbf{G}(\mathbf{p})\mathbf{x}\left(t, \mathbf{p}\right) \end{align}

with system mass, damping, and stiffness described by the time-invariant matrices \mathbf{M},\mathbf{C},\mathbf{K}\in\mathbb{R}^{n\times n}. The force acting on the system is denoted by \mathbf{f}\left(t\right)\in\mathbb{R}^{n}, the degrees of freedom \ddot{\mathbf{x}}, \dot{\mathbf{x}}, \mathbf{x}\in\mathbb{R}^{n} are acceleration, velocity, and displacement. The quantity of interest is obtained from the degrees of freedom via the output matrix \mathbf{G}(\mathbf{p}) \in \mathbb{R}^{n\times n}.

The system matrices, the force and the degrees of freedom depend on d parameters \mathbf{p}=[p_1, p_2, \dots, p_d]

Applying a Laplace transformation leads to the equation of motion in frequency domain

        

\begin{equation}         s^2\mathbf{M}(\mathbf{p})\mathbf{x}(s, \mathbf{p}) + s\mathbf{C}(\mathbf{p})\mathbf{x}(s, \mathbf{p}) + \mathbf{K}(\mathbf{p})\mathbf{x}(s, \mathbf{p}) = \mathbf{f}(s, \mathbf{p})         \end{equation}

with complex frequency s\in\mathbb{C}.

By inverting the dynamic stiffness matrix, we obtain the transfer function:

\begin{equation}         \mathbf{H}(s, \mathbf{p}) = \mathbf{G}(\mathbf{p}) \left( s^2\mathbf{M}(\mathbf{p}) + s\mathbf{C}(\mathbf{p}) + \mathbf{K}\right)^{-1}  \mathbf{f}(s, \mathbf{p})         \end{equation}

We assume that a coordinate transformation is used to reduce the full order model.The reduced quantities now also depend on the parameters \mathbf{p}=[p_1, p_2, \dots, p_d].

Parametric Model Order Reduction methods can be classified into global and local methods:

  • Global parametric Model Order Reduction
  • Local parametric Model Order Reduction
    • Interpolation of Reduced Basis
    • Interpolation of Reduced System Matrices
    • Interpolation of Transfer Function

Here, only global parametric Model Order Reduction and local parametric Model Order Reduction by interpolation of reduced system matrices are treated. For further information about the other methods refer to Benner et. al (2015) and references therein.

Sampling

For both, global and local approaches, a prior sampling is needed to generate a database of some reduced quantity. Two simple methods for sampling are Full Factorial Design and Random Sampling:

Full Factorial Design

For each parameter, a number of levels is defined (in this case 3 for p_1 and 3 for p_2) and then a regular grid is generated according to the number of levels. This method is optimal with regard to space filling, because the minimal distance between two samples is always the same. A disadvantage of the method is that the number of sampling points increases drastically in higher dimensions.

Random Sampling

An alternative is random sampling, for which the position of the samples is chosen randomly. Consequently, the number of samples can be chosen arbitrarily. However, due to the randomness a good space filling is not guaranteed. In fact, it can occur that there are large regions without any sample.


Please be aware of the notation used in the following: p_1 denotes the first parameter and \mathbf{p}_1 denotes the first sample point.

Global parametric Model Order Reduction

In Global pMOR, the reduced bases obtained in the sampling are concatenated to one global basis:

\begin{equation}     \mathbf{V} = [\mathbf{V}_1, \mathbf{V}_2, \dots, \mathbf{V}_K],     \end{equation}


with \mathbf{V}_1, \dots \mathbf{V}_K denoting the local bases corresponding to \mathbf{p}_1, \dots \mathbf{p}_K.

It is quite possible that the local bases have common components among each other, leading to a potentially rank-deficient global basis. To avoid this, the concatenation step is usually followed by an SVD to remove the rank-deficient components from \mathbf{V}.

Advantages of global pMOR are

  • the reduced order can be different for different sample points \mathbf{p}_k
  • easy to implement

Disadvantages of global pMOR are

  •   reduced size of the parametric reduced model grows with the number of samples
  •   the individual reduced bases \mathbf{V}_k for each sample point \mathbf{p}_k need to be generated on the same mesh

Local pMOR - Interpolation of Reduced Matrices

To avoid having to assemble the full order model and performing the projection in the online phase one can also  interpolate the reduced matrices.

However, for this approach a straightforward interpolation of the reduced matrices is not meaningful because they relate to different reduced bases. Let us illustrate this using the following example, where the reduced basis will be obtained by Modal Truncation:

l_x = 0.8, \, l_y = 1.0


l_x = 1.0, \, l_y = 0.95

It can be seen that a direct interpolation between the modes is not meaningful. Modes \mathbf{v}_2 and \mathbf{v}_3 for example change position and are partially mirrored. Mode \mathbf{v}_8  for l_x = 0.8, \, l_y = 1.0 is the same as mode \mathbf{v}_9 for l_x = 0.8, \, l_y = 1.0 but rotated. Similar relations can be found for the other modes as well.

To overcome this problem, the following approach was suggested in (Panzer et. al 2010

  1. Find a generalized coordinate system. For this purpose, find the most significant basis vectors by concatenating all bases and then performing a singular value decomposition:
\begin{equation}     [\mathbf{V}_1, \mathbf{V}_2, \dots \mathbf{V}_K] = \mathbf{U} \bm{\Sigma} \mathbf{Y}.     \end{equation}

    The most significant basis vectors are the first r columns in \mathbf{U} and denoted with \mathbf{R}:

    \begin{equation}     \mathbf{R} = \mathbf{U}(:,1:r).     \end{equation}

      2. Transform the individual reduced system matrices from their individual bases \mathbf{V}_k to the generalized coordinate system \mathbf{R}:

\begin{equation}  \tilde{\mathbf{K}}_r(\mathbf{p}_k) = \mathbf{T}_k^T \mathbf{K}_r(\mathbf{p}_k) \mathbf{T}_k, \quad \tilde{\mathbf{D}}_r(\mathbf{p}_k) = \mathbf{T}_k^T \mathbf{D}_r(\mathbf{p}_k) \mathbf{T}_k, \quad \tilde{\mathbf{M}}_r(\mathbf{p}_k) = \mathbf{T}_k^T \mathbf{M}_r(\mathbf{p}_k) \mathbf{T}_k, \quad \tilde{\mathbf{f}}_r(\mathbf{p}_k) = \mathbf{T}_k^T \mathbf{f}_r(\mathbf{p}_k), \quad \tilde{\mathbf{G}}_r(\mathbf{p}_k) = \mathbf{G}_r(\mathbf{p}_k) \mathbf{T}_k  \end{equation}

    with

\begin{equation}     \mathbf{T}_k = (\mathbf{R}^T\mathbf{V}_k)^{-1}.     \end{equation}

An example for this procedure is shown here.

Advantages of local parametric Model Order Reduction by Interpolation of Reduced System Matrices are

  • the size of the reduced system does not grow with the number of samples
  • efficient evaluation of the interpolated ROM in the online phase because the full order model is only needed in the offline phase

Disadvantages of local parametric Model Order Reduction by Interpolation of Reduced System Matrices are

  • in the original method, the same reduced size is required for each sampled parameter point \mathbf{p}_k. This problem was solved in (Geuss et. al 2014)
  • the individual reduced bases \mathbf{V}_k for each sample point \mathbf{p}_k need to be generated on the same mesh


Outlook: Surrogate Modeling for Frequency Response Functions

In the following slides we want to investigate the concept of surrogate modeling for frequency transfer functions. Surrogate models find widespread use in all kinds of applications that require the repeated evaluation of possibly expensive simulation models. Here we want to focus on a specific type of surrogate models that can be used to approximate the model response of structural dynamic models described through their frequency response functions (FRF).

As we know, we can interpret the FRF as a rational function in terms of the model parameters. This motivates us to approximate the FRFs using a rational approximation. The rational approximation \mathcal{R} (\mathbf{p}; \mathbf{a}, \mathbf{b}) is built through the ratio of two polynomials as follows

   

  \begin{equation*}         \mathcal{R} (\mathbf{p}; \mathbf{a}, \mathbf{b}) = \frac{\sum_{i=0}^{n_a - 1} a_i \psi_i (\mathbf{p})}{\sum_{i=0}^{n_b - 1} b_i \psi_i (\mathbf{p})}     \end{equation*}    

In there, \{a_i, i = 0, \ldots, n_a - 1 \} and \{b_i, i = 0, \ldots, n_b - 1 \} are the unknown coefficients in the expansions and the \psi_i denote a set of corresponding basis functions.

Our goal is now to determine the unknown coefficients \mathbf{a} and \mathbf{b}. Once the coefficients are known, the simple form of the model allows the fast repeated evaluation of the surrogate model.

In general, two broad categories of methods are available to determine the coefficients: intrusive and non-intrusive approaches. In (Schneider et. al 2020) a non-intrusive regression approach is presented, which we comprehensively present here. Assume now, that a sample of the input parameters \{\mathbf{p}_k, k = 1, \ldots, K \} and corresponding evaluations of a FRF at a specific frequency \{y_k = h (\mathbf{p}_k), k = 1, \ldots, K\} are available.

We define the following sample estimate of the modified mean-square error

\begin{equation}         \hat{\mathrm{err}} = \sum_{k=1}^{K} \left| h (\mathbf{p}_k) \cdot \sum_{i=0}^{n_b - 1} b_i \psi_i (\mathbf{p}_k) - \sum_{i=0}^{n_a - 1} a_i \psi_i (\mathbf{p}_k) \right|^2     \end{equation}

and find the optimal coefficients as the minimizers of the above error, i.e., 

    \begin{equation}         \{\widetilde{\mathbf{a}}, \widetilde{\mathbf{b}}\} = \underset{\mathbf{a}, \mathbf{b}}{\arg \min} \, \sum_{k=1}^{K} \left| h (\mathbf{p}_k) \cdot \sum_{i=0}^{n_b - 1} b_i \psi_i (\mathbf{p}_k) - \sum_{i=0}^{n_a - 1} a_i \psi_i (\mathbf{p}_k) \right|^2     \end{equation}

In order to solve the minimization problem, we set the derivatives of the objective functions with respect to the coefficients \mathbf{a} and \mathbf{b} to zero and solve for \mathbf{a} and \mathbf{b}.

The minimizer is the solution of the following homogeneous linear system of equations

\begin{equation} \begin{bmatrix} \bm{\Psi}_P^T \bm{\Psi}_P & - \bm{\Psi}_P^T \text{diag} \left( \mathbf{y} \right) \bm{\Psi}_Q \\ - \bm{\Psi}_Q^T \text{diag} \left( \mathbf{y}^\ast \right) \bm{\Psi}_P & \bm{\Psi}_Q^T \text{diag} \left( \mathbf{y} \circ \mathbf{y}^\ast \right) \bm{\Psi}_Q \end{bmatrix} \begin{bmatrix} \mathbf{a} \\ \mathbf{b} \end{bmatrix} = \begin{bmatrix} \mathbf{0} \\ \mathbf{0} \end{bmatrix} \, . \label{eq:pad_eq_sys_A} \end{equation}

Matrices \bm{\Psi}_P \in \mathbb{R}^{K\times n_p} and \bm{\Psi}_Q \in \mathbb{R}^{K\times n_q} have as (k,j)-element \psi_j (\mathbf{p}_k) and vector \mathbf{y} \in \mathbb{C}^K has as k-element the frequency transfer function evaluation h \left( \mathbf{p}_k \right). This homogeneous system of equations can be solved through the use of singular value decomposition. In the following we investigate the approach on a simple problem.

We choose to approximate the frequency transfer function 

\begin{equation} h (p) = \frac{\omega^2}{p - \omega^2 + \mathrm{i} \eta p} \end{equation}

with frequency \omega = 1 \mathrm{rad} \, \mathrm{Hz} and \eta = 0.01. We sample a number of points of p and solve for the coefficients in the rational approximation. The resulting approximation over p is depicted below.

References

Benner, P., Gugercin, S., & Willcox, K. (2015). A survey of projection-based model reduction methods for parametric dynamical systems. SIAM review, 57(4), 483-531.

Panzer, H., Mohring, J., Eid, R., & Lohmann, B. (2010). Parametric model order reduction by matrix interpolation.

Geuss, M., Panzer, H. K., Clifford, I. D., & Lohmann, B. (2014). Parametric model order reduction using pseudoinverses for the matrix interpolation of differently sized reduced models. IFAC Proceedings Volumes, 47(3), 9468-9473.

Schneider, F., Papaioannou, I., Ehre, M., & Straub, D. (2020). Polynomial chaos based rational approximation in linear structural dynamics with parameter uncertainties. Computers & Structures, 233, 106223.


  • Keine Stichwörter