# Development of a new vibroacoustic superelement for heavy-fluid cavities

G Buron¹,², F Thouverez¹, L Jézéquel¹, A Beley² and F Thévenon²

¹ Laboratoire de Tribologie et Dynamique des Systèmes, École Centrale de Lyon, Écully
69134, France

² Ansys France, Villeurbanne 69100, France

E-mail: guillaume.buron@ec-lyon.fr, guillaume.buron@ansys.com

**Abstract.** We propose a Reduced Order Model (ROM) of heavy-fluid cavities with the objective to build a superelement in a finite element substructuring context. We consider the non-symmetric displacement-pressure finite element formulation. While there are numerous fluid-structure ROM in the literature, few are generalized to superelements. Typical vibroacoustic substructuring methods involve modal synthesis using uncoupled rigid-wall fluid modes and *in vacuo* structure modes. Vibroacoustic domains are composed of a structural domain and a fluid domain coupled on a fluid-structure boundary. Each domain can be modeled by symmetric monophysics formulations. Only the fluid-structure boundary requires a non-symmetric multiphysics formulation. Being defined on a submanifold, this subproblem is typically significantly smaller than the other subproblems. We propose to condense the uncoupled subdomains onto the fluid-structure boundary and use a Petrov-Galerkin procedure to obtain a reduced order representation of this loaded boundary. The captured boundary dynamics is then propagated to each uncoupled subdomain. Galerkin procedures are also applied to the uncoupled subproblems to increase the reduction bases. At the end, a Petrov-Galerkin procedure is applied to the full order model to build the superelement. The proposed superelement generation strategy is used to study an industrial water tank in a seismic analysis context. Results show a significant improvement over typical uncoupled superelements at the same reduction basis size. While the uncoupled superelements cannot describe the dynamics of the model beyond the weakly-coupled sloshing regime, the proposed superelement accurately captures the strongly-coupled dynamics in the studied frequency bandwidth.

## 1. Introduction

Industrial applications nowadays often involve the study of structural dynamics in a multiphysics environment. From the automotive industry to civil engineering, vibroacoustics is a classical example of such coupling. Vibroacoustics is commonly modeled in commercial software by strongly coupled finite element formulations, implemented in monolithic solvers. The use of finite element industrial models with millions of degrees of freedom is becoming more and more common. With the digital twin paradigm, the trend is toward increasingly larger problem sizes and more complex models. Furthermore, fluid cavities typically involve large acoustic domains of numerous unknowns. In this context, ROM are required not only to speed up computations, but also to enable a substructuring approach where complex models are partitioned into multiple components for design and efficiency reasons. To that extent, superelements are defined with full order interfaces for loading and component assembly while their internal domains are reduced

<figure><img src="image_1.png" /></figure>

Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
by technics such as component mode synthesis that are not a priori specialized on a prescribed load-case.

Despite several vibroacoustic ROM proposed in the literature, the generalization to vibroacoustic superelements is currently an active subject of research. Typical vibroacoustic ROM are based on uncoupled monophysics fluid and structural modes. Commonly, fluid added mass effects on the structure are considered [8], but other extensions are made like static coupling between the structure and fluid modes [7] or structural added mass effects on the fluid [10]. ROM and superelements can be found in [9] for several vibroacoustic formulations. More complex approaches based on the enrichment of uncoupled bases have also been proposed [3][1][2].

In this study, we consider the vibroacoustic subproblem of heavy-fluid cavities with sloshing effects. We propose a vibroacoustic superelement centered around multiphysics boundary modes of the Fluid-Structure Interaction (FSI) subdomain. This superelement is used to study an industrial water tank in a seismic analysis context. Results are compared to typical uncoupled superelements, with the full order model as reference.

## Notations

In this study, we write scalars $x$, vectors $\mathbf{x}$ and matrices $\mathbf{X}$ with the indices $\gamma_k(\cdot)^\beta_\alpha$. $\alpha$ refers to a specific geometric domain. $\beta$ refers to a specific physics. $\gamma$ refers to the solution of a direct or transposed non-symmetric problem (respectively right $\gamma = \mathbf{r}$ or left $\gamma = \mathbf{l}$). $k$ refers to the index of an element in a sequence. When two indices are juxtaposed without parenthesis, it refers to the input and output of the transformation. For example, $\mathbf{A}^{12}$ is a transformation coupling the physics 1 and 2, while $\mathbf{A}_{12}$ couples the geometric domains 1 and 2. When parenthesis are used, we defines the specific meaning. We do not write the zero vector $\mathbf{0}$ and zero matrix $\mathbf{O}$ in block vectors and matrices to avoid a cluttered notation. The interior of a set is denoted by int.

## 2. Full order model

Several formulations of a vibroacoustic problem are commonly used. In this study, we consider a displacement and pressure (u, p) two-field formulation. Figure 1 describes a vibroacoustic domain of dimensionality $d$ represented by the open set $\Omega^{(up)}$. Let $\Omega^u \subset \mathbb{R}^d$ and $\Omega^p \subset \mathbb{R}^d$ be open sets that represents respectively the structural and fluid domains. The fluid-structure interaction occurs at the $(d-1)$-manifold boundary $\Gamma_b := \partial\Omega^u \cap \partial\Omega^p$. For the usual 3D case $(d=3)$, $\Omega^{(up)}$ is a volume while the FSI boundary $\Gamma_b$ is a surface. The substructuring master domain $\Omega_m^u \subset \Omega^u$ is defined for external loading or structural assembly with other components. Compressibility of the fluid is considered. A fluid free submanifold $\Gamma_{free}^p$ is present and sloshing effects are considered.

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Vibroacoustic domain $\Omega^{(up)}$.</figcaption></figure>

### 2.1. Governing equations

We consider the non-symmetric finite element (u, p) formulation described in [9] and implemented in the Ansys Mechanical APDL solver [4]. The problem is assumed to be conservative. The Full
Order Model (FOM) is the discrete problem

$$ \begin{bmatrix} \mathbf{M}^{uu} \\ -\rho_0 \mathbf{K}^{\mathrm{up}\top} \end{bmatrix} \begin{Bmatrix} \ddot{\mathbf{u}} \\ \ddot{\mathbf{p}} \end{Bmatrix} + \begin{bmatrix} \mathbf{K}^{uu} & \mathbf{K}^{\mathrm{up}} \\ \mathbf{K}^{\mathrm{pp}} & \mathbf{K}^{\mathrm{pp}} \end{bmatrix} \begin{Bmatrix} \mathbf{u} \\ \mathbf{p} \end{Bmatrix} = \begin{Bmatrix} \mathbf{f}^u \\ \mathbf{f}^p \end{Bmatrix} \quad (1) $$

$$ \Leftrightarrow: \mathbf{M} \begin{Bmatrix} \ddot{\mathbf{u}} \\ \ddot{\mathbf{p}} \end{Bmatrix} + \mathbf{K} \begin{Bmatrix} \mathbf{u} \\ \mathbf{p} \end{Bmatrix} = \mathbf{f}, \quad (2) $$

of size $n := n^u + n^p$, with $n^u$ displacement unknowns $\mathbf{u}$ and $n^p$ pressure unknowns $\mathbf{p}$. The $\mathbf{K}$ and $\mathbf{M}$ matrices are sparse thanks to the finite element formulation. Sloshing effects due to gravity are accounted for in the fluid matrix $\mathbf{M}^{\mathrm{pp}}$. $\rho_0$ is the mean fluid density.

## 2.2. Left/right eigenspaces

The direct $(\mathbf{K}, \mathbf{M})$ (3) and transposed $(\mathbf{K}^{\top}, \mathbf{M}^{\top})$ (4) eigenproblems define the modes of the FOM (1), with $\mathbf{\Lambda}$ its real spectrum, and $^{\ell}\mathbf{\Phi}$ and $^{r}\mathbf{\Phi}$ its real left and right eigenvectors. From the sequence of increasingly larger eigenvalues, we define $\mathbf{\Lambda}$ as the spectral matrix diag($\omega^2, \dots, \omega^2$). From the corresponding sequence of column eigenvectors, we define $^{\ell}\mathbf{\Phi}$ as the modal matrix $[\mathbf{\iota}^1 \dots \mathbf{\iota}^n]$.

$$ \mathbf{K}^{\mathrm{r}} \mathbf{\Phi} = \mathbf{M}^{\mathrm{r}} \mathbf{\Phi} \mathbf{\Lambda} \Leftrightarrow: \begin{bmatrix} \mathbf{K}^{\mathrm{uu}} & \mathbf{K}^{\mathrm{up}} \\ \mathbf{K}^{\mathrm{pp}} & \mathbf{K}^{\mathrm{pp}} \end{bmatrix} \begin{bmatrix} \mathbf{r}^{\mathrm{\Phi}^u} \\ \mathbf{r}^{\mathrm{\Phi}^p} \end{bmatrix} = \begin{bmatrix} \mathbf{M}^{\mathrm{uu}} \\ -\rho_0 \mathbf{K}^{\mathrm{up}\top} \end{bmatrix} \begin{bmatrix} \mathbf{r}^{\mathrm{\Phi}^u} \\ \mathbf{r}^{\mathrm{\Phi}^p} \end{bmatrix} \mathbf{\Lambda}, \quad (3) $$

$$ \mathbf{K}^{\mathrm{T}} \mathbf{\ell} \mathbf{\Phi} = \mathbf{M}^{\mathrm{T}} \mathbf{\ell} \mathbf{\Phi} \mathbf{\Lambda} \Leftrightarrow: \begin{bmatrix} \mathbf{K}^{\mathrm{uu}} & \mathbf{K}^{\mathrm{pp}} \\ \mathbf{K}^{\mathrm{up}\top} & \mathbf{K}^{\mathrm{pp}} \end{bmatrix} \begin{bmatrix} \mathbf{\ell}^{\mathrm{\Phi}^u} \\ \mathbf{\ell}^{\mathrm{\Phi}^p} \end{bmatrix} = \begin{bmatrix} \mathbf{M}^{\mathrm{uu}} & -\rho_0 \mathbf{K}^{\mathrm{up}} \\ \mathbf{M}^{\mathrm{pp}} & \mathbf{M}^{\mathrm{pp}} \end{bmatrix} \begin{bmatrix} \mathbf{\ell}^{\mathrm{\Phi}^u} \\ \mathbf{\ell}^{\mathrm{\Phi}^p} \end{bmatrix} \mathbf{\Lambda}. \quad (4) $$

The eigenvectors may be normalized to enforce either biorthonormalization $^{\ell}\mathbf{\Phi}^{\mathrm{T}} \mathbf{r}\mathbf{\Phi} = \mathbf{I}$ or $\mathbf{M}$-biorthogonality $^{\ell}\mathbf{\Phi}^{\mathrm{T}} \mathbf{M}^{\mathrm{r}}\mathbf{\Phi} = \mathbf{I}$.

The specific structure of problems (4-3) causes the spectrum to be real and leads to a relation between the left and right eigenvectors. For this vibroacoustic formulation, it can be proven that a left eigenvector $^{\ell}\mathbf{\phi}$ is deduced from its corresponding right eigenvector $^{r}\mathbf{\phi}$ by

$$ \begin{Bmatrix} \mathbf{\ell}^{\mathrm{\Phi}^u} \\ \mathbf{\ell}^{\mathrm{\Phi}^p} \end{Bmatrix} = \begin{bmatrix} \mathbf{\iota}^{\omega^2 \mathbf{I}} & \mathbf{\iota}^{\mathbf{r}^{\mathrm{\Phi}^u}} \\ \mathbf{\iota}^{\mathbf{r}^{\mathrm{\Phi}^p}} & \mathbf{\iota}^{\mathbf{r}^{\mathrm{\Phi}^p}} \end{bmatrix} \begin{Bmatrix} \mathbf{r}^{\mathrm{\Phi}^u} \\ \mathbf{r}^{\mathrm{\Phi}^p} \end{Bmatrix}. \quad (5) $$

While this structure enables symmetrization approaches as discussed in [1][2], we consider the non-symmetric problem directly.

## 2.3. Singularities

An in vacuo vibroacoustic component may have a singular FOM (1) due to rigid body motions. However, in a substructuring context, boundary conditions or assembly on the substructuring subdomain $\Gamma_m^u$ remove the rigid body motions in most cases.

In the presence of a fluid free surface $\Gamma_{\mathrm{free}}^p$, a fluid nullspace mode $(^{\ell}\mathbf{\phi}, ^{r}\mathbf{\phi})$ will appear, leading to a singular problem (1). It corresponds to an uniform pressure field inside the acoustic domain and the corresponding static deformation of the structure. It can be determined as the solution of

$$ \begin{bmatrix} \mathbf{K}^{\mathrm{uu}} & \mathbf{K}^{\mathrm{up}} \\ \mathbf{K}^{\mathrm{pp}} & \mathbf{K}^{\mathrm{pp}} \end{bmatrix} \begin{bmatrix} \mathbf{r}^{\mathrm{\Phi}^u} \\ \mathbf{1} \end{bmatrix} = \mathbf{0}, \quad (6) $$

with the trivial left eigenvector

$$ ^{\ell}\mathbf{\phi} = \begin{bmatrix} \mathbf{0} \\ \mathbf{1} \end{bmatrix} \quad (7) $$
due to no coupling $p \to u$ at order 0 of eigenproblem (4). With this procedure, a normalization step will be necessary to enforce biorthogonality.

In this study, we regularize ill-conditioned eigenproblems with a small initial spectral shift. When we study an ill-conditioned linear problem or when we need an explicit pseudo-inverse of a singular matrix, we use the bordering regularization technic [5]. These methods preserve matrix sparsity.

## 2.4. Discrete domain subdivision

When building a superelement, the chosen master domain $\Omega_m^u$ stays in physical space. The rest of the domain $\Omega^{(up)} \setminus \Omega_m^u$, called condensed domain, will be reduced to generalized unknowns. The mesh used for the discrete problem (1) approximates the space corresponding to $\Omega^{(up)}$. We partition the problem (1) by defining nodal subsets of this mesh.

### 2.4.1. Condensed-master subdivision
Let $C_m^u$ be the set of master nodes corresponding to $\Omega_m^u$, $C_c^u$ the nodes corresponding to $\Omega^u \setminus \Omega_m^u$ and $C_c^p$ the ones corresponding to $\Omega^p$. From (1), this leads to the partitioned problem

$$ \begin{bmatrix} \mathbf{M}_{cc}^{uu} & \mathbf{M}_{cm}^{uu} \\ -\rho_0 \mathbf{K}_{cc}^{up\top} & \mathbf{M}_{cc}^{pp} \\ \mathbf{M}_{cm}^{uu\top} & \mathbf{M}_{mm}^{uu} \end{bmatrix} \begin{Bmatrix} \ddot{u}_c \\ \ddot{p}_c \\ \ddot{u}_m \end{Bmatrix} + \begin{bmatrix} \mathbf{K}_{cc}^{uu} & \mathbf{K}_{cc}^{up} & \mathbf{K}_{cm}^{uu} \\ \mathbf{K}_{cc}^{pp} & \mathbf{K}_{cc}^{up} & \mathbf{K}_{cm}^{uu} \\ \mathbf{K}_{cm}^{uu\top} & \mathbf{K}_{mc}^{up} & \mathbf{K}_{mm}^{uu} \end{bmatrix} \begin{Bmatrix} \mathbf{u}_c \\ \mathbf{p}_c \\ \mathbf{u}_m \end{Bmatrix} = \begin{Bmatrix} \mathbf{f}_c^u \\ \mathbf{f}_c^p \\ \mathbf{f}_m^u \end{Bmatrix}. \quad (8) $$

### 2.4.2. Internal-boundary-master subdivision
We further partition the condensed domain into multiphysics and monophysics subdomains. The FSI boundary nodes corresponding to $\Gamma_b \setminus \partial\Omega_m^u$ are $\Sigma_b := C_c^u \cap C_c^p$. Then, internal nodes of each monophysics subdomain belong to $C_i^u := C_c^u \setminus \Sigma_b$ and $C_i^p := C_c^p \setminus \Sigma_b$. From (1), this leads to the partitioned problem

$$ \begin{bmatrix} \mathbf{M}_{ii}^{uu} & \mathbf{M}_{ib}^{uu} & \mathbf{M}_{im}^{uu} \\ \mathbf{M}_{ib}^{uu\top} & \mathbf{M}_{bb}^{pp} & \mathbf{M}_{bb}^{pp} \\ \mathbf{M}_{im}^{uu\top} & \mathbf{M}_{bm}^{pp\top} & \mathbf{M}_{mm}^{uu} \end{bmatrix} \begin{Bmatrix} \ddot{u}_i \\ \ddot{p}_i \\ \ddot{u}_b \\ \ddot{p}_b \\ \ddot{u}_m \end{Bmatrix} + \begin{bmatrix} \mathbf{K}_{ii}^{uu} & \mathbf{K}_{ib}^{pp} & \mathbf{K}_{ib}^{uu} & \mathbf{K}_{im}^{uu} \\ \mathbf{K}_{ib}^{uu\top} & \mathbf{K}_{ii}^{pp} & \mathbf{K}_{bb}^{uu} & \mathbf{K}_{bb}^{up} \\ \mathbf{K}_{im}^{uu\top} & \mathbf{K}_{mi}^{pp\top} & \mathbf{K}_{bb}^{pp} & \mathbf{K}_{mb}^{up} \\ \mathbf{K}_{bm}^{uu\top} & \mathbf{K}_{mi}^{up} & \mathbf{K}_{bm}^{uu\top} & \mathbf{K}_{mb}^{up} \end{bmatrix} \begin{Bmatrix} \mathbf{u}_i \\ \mathbf{p}_i \\ \mathbf{u}_b \\ \mathbf{p}_b \\ \mathbf{u}_m \end{Bmatrix} = \begin{Bmatrix} \mathbf{f}_i^u \\ \mathbf{f}_i^p \\ \mathbf{f}_b^u \\ \mathbf{f}_b^p \\ \mathbf{f}_m^u \end{Bmatrix}. \quad (9) $$

Note that the internal subdomains $C_i^u$ and $C_i^p$ are not coupled, thus their corresponding subproblems in (9) are monophysics and symmetric. Only the FSI boundary $\Sigma_b$ subproblem is multiphysics and non-symmetric. Being defined on a $(d-1)$-manifold, this subproblem is typically significantly smaller than the condensed domain $C_c^u \cup C_c^p$ subproblem from (8) that is defined on a $d$-manifold.

## 3. Reduced order model

In this study, we consider the fixed-interface approach to build a superelement. Thus, a reduced order representation of the condensed domain $C_c^u \cup C_c^p$ is needed with the $\mathbf{u}_m = \mathbf{0}$ boundary condition. Similarly to the Hurty/Craig-Bampton method [6] used in structural dynamics, the master domain dynamics is propagated to the condensed domain through constraint modes.

### 3.1. Uncoupled superelement

The most common vibroacoustic ROM method is to use monophysics modes from the uncoupled structural and fluid subproblems for the reduction. This classical strategy, leading to an uncoupled superelement, is compared to our approach. We use two technics to reduce the condensed domain $C_c^u \cup C_c^p$: one using in vacuo structure modes and another one using modes of the structure loaded with the fluid. In both approaches, rigid-wall fluid modes are used.
### 3.1.1. Uncoupled reduction
From partitioned FOM problem (8), applying the fixed-interface $\mathbf{u}_m := \mathbf{0}$ and rigid-wall $\mathbf{u}_c := \mathbf{0}$ conditions leads to the uncoupled acoustic eigenproblem

$$ \mathbf{K}_{cc}^{pp} \mathbf{\Phi}_{c}^{p} = \mathbf{M}_{cc}^{pp} \mathbf{\Phi}_{c}^{p} \mathbf{\Lambda}_{c}^{p}. \quad (10) $$

Similarly to the FOM eigenproblem (3), eigenproblem (10) is ill-conditioned due to a fluid nullspace mode. It is regularized by an initial spectral shift $\mathcal{K}_{cc}^{pp} := \mathbf{K}_{cc}^{pp} - s_0 \mathbf{M}_{cc}^{pp}$ where $s_0$ is chosen smaller than the lowest non-zero eigenvalue of (10).

From partitioned FOM problem (8), applying the fixed-interface $\mathbf{u}_m := \mathbf{0}$ and in vacuo structure $\mathbf{p}_c := \mathbf{0}$ conditions leads to the uncoupled structural eigenproblem

$$ \mathbf{K}_{cc}^{uu} \mathbf{\Phi}_{c}^{u} = \mathbf{M}_{cc}^{uu} \mathbf{\Phi}_{c}^{u} \mathbf{\Lambda}_{c}^{u}. \quad (11) $$

### 3.1.2. Added mass effects
Assuming static response of the fluid, fluid added mass effects can be accounted for several equivalent ways. Similar to [8], the uncoupled fluid subproblem is regularized by the condition $\mathbf{p} = \mathbf{0}$ on the corresponding fluid free surface $C_{cfree}^p := C_c^p \setminus (\text{int} C_c^p \cup \Sigma_b)$. The added mass matrix from the interior fluid domain $C_{c_i}^p := C_c^p \setminus C_{cfree}^p$ is $\mathbf{M}_{cc}^{(p \to u)} = -\rho_0 \mathbf{K}_{cc_i}^{\text{up}} \mathbf{K}_{c_i c_i}^{\text{pp}-1} \mathbf{K}_{cc_i}^{\text{up}\top}$, leading to the loaded uncoupled structural eigenproblem

$$ \mathbf{K}_{cc}^{uu} \mathbf{\Phi}_{c}^{u} = (\mathbf{M}_{cc}^{uu} + \mathbf{M}_{cc}^{(p \to u)}) \mathbf{\Phi}_{c}^{u} \mathbf{\Lambda}_{c}^{u}. \quad (12) $$

### 3.1.3. Monophysics superelement generalization
Considering the monophysics approach to the reduction of the condensed domain, we also take a monophysics approach to the constraint modes between the condensed and master domains, similar to [3]. For that purpose, only the structural constraint modes (13) need to be computed.

$$ \mathbf{\Psi}_{cm}^{u} = -\mathbf{K}_{cc}^{uu-1} \mathbf{K}_{cm}^{uu}. \quad (13) $$

The bases $\mathbf{\Phi}_{c}^{u}$, $\mathbf{\Phi}_{c}^{u}$ and $\mathbf{\Phi}_{c}^{p}$ are truncated up to an arbitrary eigenfrequency threshold we call cutoff frequency, or based on a target ROM size. From the truncated bases $\tilde{\mathbf{\Phi}}_{c}^{u}$, $\tilde{\mathbf{\Phi}}_{c}^{u}$, $\tilde{\mathbf{\Phi}}_{c}^{p}$ and constraint modes $\mathbf{\Psi}_{cm}^{u}$, we build the reduction bases, leading to the approximation

$$ \left\{ \begin{matrix} \mathbf{u}_c \\ \mathbf{p}_c \\ \mathbf{u}_m \end{matrix} \right\} \approx \left[ \begin{matrix} \tilde{\mathbf{\Phi}}_{c}^{u} & \mathbf{\Psi}_{cm}^{u} \\ \mathbf{\Phi}_{c}^{p} & \mathbf{I} \end{matrix} \right] \left\{ \begin{matrix} \mathbf{q}_{c}^{u} \\ \mathbf{q}_{c}^{p} \\ \mathbf{u}_{m} \end{matrix} \right\} := \mathbf{\Theta} \left\{ \begin{matrix} \mathbf{q}_{c}^{u} \\ \mathbf{q}_{c}^{p} \\ \mathbf{u}_{m} \end{matrix} \right\}, \quad (14) $$

and the similar loaded structure approximation with $\tilde{\mathbf{\Phi}}_{c}^{u}$ defining the $\mathbf{\Theta}$ reduction basis.

The ROM matrices can be computed by the Galerkin method with the reduction basis $\mathbf{\Theta}$ or $\mathbf{\Theta}$. We call respectively uncoupled ROM and loaded uncoupled ROM these problems of size $\tilde{n}_c^u + \tilde{n}_c^p + n_m^u$ and $\tilde{n}_c^u + \tilde{n}_c^p + n_m^u$. However, we choose to propagate the singularity of the FOM (1) to the ROM in order to avoid making any assumption in this study on the choice of regularization technic. For this purpose, we increase the reduction bases with the non-symmetric fluid nullspace mode (7-6), leading, without loss of generality, to a Petrov-Galerkin method.

## 3.2. Coupled superelement

While the monophysics modes from $\mathbf{\Theta}$ (14) or $\mathbf{\Theta}$ will be coupled during the Galerkin projection, we make the assumption that accounting for some multiphysics information in the reduction basis itself should improve the ROM performance for strongly-coupled models. Considering the partitioning (9), we propose to compute multiphysics modes of the FSI boundary $\Sigma_b$ while keeping a monophysics approach to the uncoupled internal structural $C_i^u$ and fluid $C_i^p$ subdomains.
### 3.2.1. Domain subdivision-based reduction
From partitioned FOM problem (9), applying the fixed-interface $\mathbf{u}_m := \mathbf{0}$ condition leads to the coupled subproblem

$$ \begin{bmatrix} \mathbf{M}_{ii}^{uu} & & \mathbf{M}_{ib}^{uu} & & \mathbf{M}_{ib}^{pp} \\ & \mathbf{M}_{ii}^{pp} & & \mathbf{M}_{ib}^{pp} & \\ & & \mathbf{M}_{bb}^{uu} & & \mathbf{M}_{bb}^{pp} \\ & & & -\rho_0 \mathbf{K}_{bb}^{up\top} & \mathbf{M}_{bb}^{pp} \\ \mathbf{M}_{ib}^{uu\top} & & & & \mathbf{M}_{ib}^{pp\top} \\ & & & & -\rho_0 \mathbf{K}_{bb}^{up\top} & \mathbf{M}_{bb}^{pp} \end{bmatrix} \begin{Bmatrix} \ddot{u}_i \\ \ddot{p}_i \\ \ddot{u}_b \\ \ddot{p}_b \end{Bmatrix} + \begin{bmatrix} \mathbf{K}_{ii}^{uu} & & \mathbf{K}_{ib}^{uu} & & \mathbf{K}_{ib}^{pp} \\ & \mathbf{K}_{ii}^{pp} & & \mathbf{K}_{ib}^{pp} & \\ & & \mathbf{K}_{bb}^{uu} & & \mathbf{K}_{bb}^{up} \\ & & & \mathbf{K}_{bb}^{pp} & \\ \mathbf{K}_{ib}^{uu\top} & & & & \mathbf{K}_{ib}^{pp\top} \\ & & & & -\rho_0 \mathbf{K}_{bb}^{up\top} & \mathbf{M}_{bb}^{pp} \end{bmatrix} \begin{Bmatrix} \mathbf{u}_i \\ \mathbf{p}_i \\ \mathbf{u}_b \\ \mathbf{p}_b \end{Bmatrix} = \begin{Bmatrix} \mathbf{f}_i^u \\ \mathbf{f}_i^p \\ \mathbf{f}_b^u \\ \mathbf{f}_b^p \end{Bmatrix}. \quad (15) $$

Assuming static response of the internal structural $C_i^u$ and fluid $C_i^p$ subdomains. The internal unknowns are condensed by

$$ \left\{ \begin{matrix} \mathbf{u}_i \\ \mathbf{p}_i \\ \mathbf{u}_b \\ \mathbf{p}_b \end{matrix} \right\} \approx \begin{bmatrix} -\mathbf{K}_{ii}^{uu-1} \mathbf{K}_{ib}^{uu} & & & \\ & \mathbf{I} & & \\ & & -\mathbf{K}_{ii}^{pp-1} \mathbf{K}_{ib}^{pp} & \\ & & & \mathbf{I} \end{bmatrix} \left\{ \begin{matrix} \mathbf{u}_b \\ \mathbf{p}_b \end{matrix} \right\} =: \begin{bmatrix} \mathbf{\Psi}_{ib}^u & & \\ & \mathbf{\Psi}_{ib}^p & \\ & & \mathbf{I} \end{bmatrix} \left\{ \begin{matrix} \mathbf{u}_b \\ \mathbf{p}_b \end{matrix} \right\}. \quad (16) $$

This is equivalent to studying the FSI boundary with the added mass of the internal structural $C_i^u$ and fluid $C_i^p$ subdomains. After applying by Galerkin projection the reduction (16) to problem (15), we deduce the non-symmetric loaded FSI boundary eigenproblem

$$ \underline{\mathbf{K}}_{bb} \begin{bmatrix} \mathbf{r} \underline{\mathbf{\Phi}}_b^u \\ \mathbf{r} \underline{\mathbf{\Phi}}_b^p \end{bmatrix} = \underline{\mathbf{M}}_{bb} \begin{bmatrix} \mathbf{r} \underline{\mathbf{\Phi}}_b^u \\ \mathbf{r} \underline{\mathbf{\Phi}}_b^p \end{bmatrix} \underline{\mathbf{\Lambda}}_b, \quad (17) $$

with the left modes $\mathbf{r}\underline{\mathbf{\Phi}}_b$ deduced with a relation similar to (5).

Exploiting the fact that $C_i^u \cap C_i^p = \emptyset$, the FSI boundary is simultaneously loaded by the internal structural and fluid subdomains through the monophysics constraint modes $\mathbf{\Psi}_{ib}^u$ and $\mathbf{\Psi}_{ib}^p$. Because the FSI boundary $\Sigma_b$ is a submanifold only comprised of the wetted face of the structure and the corresponding fluid layer, the eigenproblem (17) is typically much smaller than the eigenproblem corresponding to problem (15). This is an opportunity to capture efficiently in $\underline{\mathbf{\Lambda}}_b$ and $(\mathbf{r}\underline{\mathbf{\Phi}}_b, \mathbf{r}\underline{\mathbf{\Phi}}_b)$ some of the multiphysics dynamics of problem (15).

Internal subdomain modes are also computed. From partitioned problem (15), applying both the fixed-interface $\mathbf{u}_m := \mathbf{0}$ and the fixed FSI boundary $\mathbf{u}_b := \mathbf{0}$ and $\mathbf{p}_b := \mathbf{0}$ conditions leads to the monophysics internal eigenproblems

$$ \mathbf{K}_{ii}^{uu} \mathbf{\Phi}_i^u = \mathbf{M}_{ii}^{uu} \mathbf{\Phi}_i^u \underline{\mathbf{\Lambda}}_i^u, \quad (18) $$

$$ \mathbf{K}_{ii}^{pp} \mathbf{\Phi}_i^p = \mathbf{M}_{ii}^{pp} \mathbf{\Phi}_i^p \underline{\mathbf{\Lambda}}_i^p. \quad (19) $$

### 3.2.2. Multiphysics superelement generalization
Considering our multiphysics approach to the reduction of the condensed domain, we also take a multiphysics approach to the constraint modes between the condensed and master domains, similar to [1][2]. To that extent, we compute the multiphysics left and right constraint modes

$$ \begin{bmatrix} \ell \mathbf{\Psi}_{im}^{uu} \\ \ell \mathbf{\Psi}_{im}^{pu} \\ \ell \mathbf{\Psi}_{bm}^{uu} \\ \ell \mathbf{\Psi}_{bm}^{pu} \end{bmatrix} = - \begin{bmatrix} \mathbf{K}_{ii}^{uu} & & \mathbf{K}_{ib}^{uu} & & \mathbf{K}_{ib}^{pp} \\ & \mathbf{K}_{ib}^{pp} & & \mathbf{K}_{ib}^{pp} & \\ & & \mathbf{K}_{bb}^{uu} & & \mathbf{K}_{bb}^{up} \\ & & & \mathbf{K}_{bb}^{pp} & \\ \mathbf{K}_{ib}^{uu\top} & & & & \mathbf{K}_{ib}^{pp\top} \\ & & & & -\rho_0 \mathbf{K}_{bb}^{up\top} & \mathbf{M}_{bb}^{pp} \end{bmatrix}^{-\top} \begin{bmatrix} \mathbf{K}_{im}^{uu} \\ \mathbf{K}_{mi}^{up\top} \\ \mathbf{K}_{bm}^{uu} \\ \mathbf{K}_{mb}^{up\top} \end{bmatrix}, \quad (20) $$

$$ \begin{bmatrix} \mathbf{r} \mathbf{\Psi}_{im}^{uu} \\ \mathbf{r} \mathbf{\Psi}_{im}^{pu} \\ \mathbf{r} \mathbf{\Psi}_{bm}^{uu} \\ \mathbf{r} \mathbf{\Psi}_{bm}^{pu} \end{bmatrix} = - \begin{bmatrix} \mathbf{K}_{ii}^{uu} & & \mathbf{K}_{ib}^{uu} & & \mathbf{K}_{ib}^{pp} \\ & \mathbf{K}_{ii}^{pp} & & \mathbf{K}_{ib}^{pp} & \\ & & \mathbf{K}_{bb}^{uu} & & \mathbf{K}_{bb}^{up} \\ & & & \mathbf{K}_{bb}^{pp} & \\ \mathbf{K}_{ib}^{uu\top} & & & & \mathbf{K}_{ib}^{pp\top} \\ & & & & -\rho_0 \mathbf{K}_{bb}^{up\top} & \mathbf{M}_{bb}^{pp} \end{bmatrix}^{-1} \begin{bmatrix} \mathbf{K}_{im}^{uu} \\ \mathbf{K}_{bm}^{uu} \end{bmatrix}. \quad (21) $$
The bases $\tilde{\Phi}_i^u, \tilde{\Phi}_i^p$ and $(^\ell\tilde{\Phi}_b, {}^r\tilde{\Phi}_b)$ are truncated either by applying a chosen frequency cutoff or based on a chosen ROM size. From the truncated bases $\tilde{\Phi}_i^u, \tilde{\Phi}_i^p$ and $(^\ell\tilde{\Phi}_b, {}^r\tilde{\Phi}_b)$, and constraint modes (16), (20) and (21), we build the reduction bases, leading to the approximation

$$ \left\{ \begin{matrix} \mathbf{u}_i \\ \mathbf{p}_i \\ \mathbf{u}_b \\ \mathbf{p}_b \\ \mathbf{u}_m \end{matrix} \right\} \approx \left[ \begin{matrix} \tilde{\Phi}_i^u & \Psi_{ib}^u {}^r\tilde{\Phi}_b^u & {}^r\Psi_{im}^{uu} \\ \tilde{\Phi}_i^p & \Psi_{ib}^p {}^r\tilde{\Phi}_b^p & {}^r\Psi_{im}^{pu} \\ {}^r\tilde{\Phi}_b^u & {}^r\Psi_{bm}^{uu} & {}^r\Psi_{bm}^{pu} \\ {}^r\tilde{\Phi}_b^p & {}^r\Psi_{bm}^{pu} & \mathbf{I} \\ \end{matrix} \right] \left\{ \begin{matrix} \mathbf{q}_i^u \\ \mathbf{q}_i^p \\ \underline{\mathbf{q}}_b \\ \mathbf{u}_m \end{matrix} \right\} := {}^r\mathfrak{K} \left\{ \begin{matrix} \mathbf{q}_i^u \\ \mathbf{q}_i^p \\ \underline{\mathbf{q}}_b \\ \mathbf{u}_m \end{matrix} \right\}, \quad (22) $$

with the corresponding left reduction basis

$$ {}^{\ell}\mathfrak{K} := \left[ \begin{matrix} \tilde{\Phi}_i^u & \Psi_{ib}^{\ell} \tilde{\Phi}_b^u & {}^{\ell}\Psi_{im}^{uu} \\ \tilde{\Phi}_i^p & \Psi_{ib}^p \tilde{\Phi}_b^p & {}^{\ell}\Psi_{im}^{pu} \\ {}^{\ell}\tilde{\Phi}_b^u & {}^{\ell}\Psi_{bm}^{uu} & {}^{\ell}\Psi_{bm}^{pu} \\ {}^{\ell}\tilde{\Phi}_b^p & {}^{\ell}\Psi_{bm}^{pu} & \mathbf{I} \end{matrix} \right]. \quad (23) $$

By the Petrov-Galerkin method with reduction bases ($^{\ell}\mathfrak{K}$, $^{r}\mathfrak{K}$), we compute the ROM of size $\tilde{n}_i^u + \tilde{n}_i^p + \tilde{n}_b + n_m^u$ that we call *coupled* ROM. Without loss of generality, we also increase the reduction bases with the non-symmetric fluid nullspace mode (7-6) to propagate the FOM singularity to the ROM.

The superelement we propose could be viewed as a nested generalization of fixed-interface substructuring to vibroacoustic problems. At the vibroacoustic component scale, this is a Petrov-Galerkin generalization of the fixed-interface method. At the condensed domain scale, with the FSI boundary considered a pseudo master domain, it is two Galerkin structural and fluid fixed-interface methods with a further multiphysics interface reduction.

# 4. Numerical testing

## 4.1. Experimental procedure

We use as reference results the FOM truncated modes $\tilde{\Lambda}$ and $(^\ell\tilde{\Phi}, {}^r\tilde{\Phi})$ from eigenproblems (4) and (3). We perform modal analyses after reduction on the uncoupled ROM, the loaded uncoupled ROM and the proposed coupled ROM. Below a given cutoff frequency, we consider a ROM accurate if both its spectrum does not deviate significantly from the FOM spectrum, and its eigenvectors expanded in the FOM space span a similar subspace as the FOM eigenvectors. We choose not to compute explicitly the matrix of the frequency response functions because it becomes too computationally intensive even for relatively small models. It is also challenging to interpret the results for numerous degrees of freedom and having a very high modal density.

The typical indicator of eigenvector similarity is the Modal Assurance Criterion (MAC). Two eigenspaces can be compared by quantifying with the MAC the similarity of each of their corresponding eigenvectors. However, in the general case of a model with high modal density and with some degenerated modes, implementing this procedure is not trivial. Because some modes may be missing or spurious in one of the eigenspace, the mode numbers are not equivalent and cannot be used directly for identifying the eigenvectors to compare. Although modes can be matched between the two result sources based on the nearest eigenvalues, additional considerations must be taken for degenerated modes. For each repeated eigenvalue, using the MAC requires studying the entire modal subspace. Degeneracy is almost always due to geometric
symmetries. When the mesh structure breaks these symmetries, the computed eigenvalues of a degenerated mode can significantly drift apart. At high modal densities, it can be impractical to identify in a robust way all the repeated eigenvalues and the corresponding modal subspaces. Vibroacoustic cavities often have geometric symmetries. Therefore, as an indicator of eigenspace similarity that does not require mode matching, we propose a generalization of the subspace similarity indicator proposed in [11]. Instead of always considering two bases of the same size, we allow for one basis to be larger. We interpret the indicator as the capability of a basis to span the subspace of the reference basis. In details, this indicator quantifies the loss of orthonormality of the reference basis after application of the projector corresponding to the basis to characterize. The studied bases need to be orthonormal, thus, in the general case of a non-symmetric and non-biorthonormal transformation ($^{ℓ}$T, $^{r}$T), we compute its orthogonal projector $\mathcal{P}$($^{ℓ}$T, $^{r}$T) as

$$ \mathcal{P}(^{\ell}\mathbf{T}, ^{r}\mathbf{T}) = ^{r}\mathbf{T} (^{\ell}\mathbf{T}^{\top} ^{r}\mathbf{T})^{-1} ^{\ell}\mathbf{T}^{\top}. \qquad (24) $$

As a first approach, we investigate if the reduced space of each superelement is a good approximation of the considered FOM eigenspace. While it does not guarantee an accurate dynamics after reduction, the reduction bases should at least accurately span the FOM eigenspace.

When comparing a ROM modal subspace to the FOM modal subspace, the truncation strategy must be chosen with care. Because mode numbers are not in general equivalent, truncating the two eigenspaces based on size may lead to omitting some degenerated eigenvectors at the end of a basis. This artificially lowers the similarity indicator. A frequency truncation circumvents that limitation. However, if there is a mode near the cutoff frequency with an overestimated ROM frequency, it will be omitted. Thus, we arbitrary choose to relax the cutoff frequencies of our data points by +5% to only penalize bases with significantly overestimated frequencies.

## 4.2. Seismic analysis of an industrial water tank

### 4.2.1. Model
As an example of a strongly-coupled model with sloshing effects, we use the described superelements and experimental procedure to study an industrial water tank. Compressibility effects are not significant in this model, however we use the general vibroacoustic formulation considered in this development. The tank is of radius 7 m and height 12 m. Its steel structure of thickness 1.5 cm is modeled with quadratic quadrangle shell elements. The tank is filled to 80% of its height with water modeled with quadratic hexahedral acoustic elements. The mesh used is represented by figure 2 and corresponds to 3056 degrees of freedom: 2184 displacement unknowns and 872 pressure unknowns. The FSI boundary $\Sigma_b$ corresponds to 1080 displacement unknowns and 204 pressure unknowns. Because quadratic elements are used, the partitioning requires that some additional pressure unknowns are assigned to the FSI boundary $\Sigma_b$. The master domain $C_m^u$ is defined as the tank bottom face and corresponds to 672 displacement unknowns. An elastic support is added to the tank bottom to model concrete foundations. Possible goals for this analysis setup would be to load the tank from the ground or study directly its interaction with its foundations. In a low-frequency seismic analysis context, we consider the FOM modes below 50 Hz as our reference results.

The coupled ROM bases are truncated to 50 Hz, leading to a total superelement size of 898 unknowns (226 modal and 672 master displacement unknowns). For a fair comparaison with the less complex uncoupled and loaded uncoupled ROM, we truncate the uncoupled bases to reach the same superelement size. This leads to truncating the uncoupled ROM bases at around 110 Hz and the loaded uncoupled ROM bases at around 70 Hz.

### 4.2.2. Results
Figure 3 shows the ROM modal frequencies plotted versus the FOM frequencies they are supposed to correspond to. Deviation from the reference FOM line is explained either
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Tank mesh.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Spectrum deviations between ROM and FOM.</figcaption></figure>

by errors on the frequencies, or missing or spurious modes shifting the mode numbers. We can identify below 2 Hz a weakly-coupled sloshing regime where the dynamics is dominated by fluid sloshing modes statically coupled with the structure. The three ROM accurately describe the FOM spectrum in that regime. However, only the coupled ROM we propose is accurate beyond this regime into strongly-coupled solutions. The coupled ROM spectrum can be considered accurate until around 35 Hz.

Figure 4 shows the subspace similarity between the spaces spanned by the reduction bases and the FOM eigensubspaces corresponding to truncations of the FOM eigenbasis at progressively larger cutoff frequencies. The accuracy of the uncoupled reduction basis declines rapidly despite its high frequency truncation of 110 Hz relative to the bandwidth of interest. This suggests that the reduction basis is inadequate for this problem. Both the loaded uncoupled and coupled reduction bases are satisfactory representations of the FOM eigenspace across the studied bandwidth. However, the proposed coupled method demonstrates an appreciable improvement over the loaded uncoupled method below 30 Hz.

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Reduction bases and FOM eigenspace subspace similarity.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** ROM and FOM eigenspaces subspace similarity.</figcaption></figure>

Figure 5 shows the modal subspace similarity of the ROM eigenspaces and the FOM eigenspace, both truncated at a progressively larger cutoff frequency. The loaded uncoupled ROM eigenspace becomes dissimilar beyond the weakly-coupled sloshing regime due to numerous missing modes at higher frequencies. The coupled ROM eigenspace maintains a high level of
similarity up to around 35 Hz.

**5. Discussion**

This study shows that superelements built from uncoupled bases cannot describe a strongly-coupled vibroacoustic cavity beyond the weakly-coupled sloshing regime. We proposed a reduction method centered around loaded FSI boundary modes. In a realistic analysis setup, the proposed superelement demonstrates accurate results with a significant model size reduction. Therefore, this superelement can be used for substructuring of strongly-coupled vibroacoustic heavy-fluid cavities with compressibility and sloshing effects.

The proposed reduction procedure is designed to scale efficiently to most large problems by exploiting the subdimensional property of the fluid-structure boundary. In the present work, we have not quantified the computational cost of our method. We will investigate the optimization of an eigensolver for the computation of the loaded FSI boundary modes. We will then characterize the computational efficiency of our reduction procedure and compare it to a fully multiphysics approach.

It would also be interesting to compare the accuracy and computational cost of the proposed superelement to methods based on enrichment of uncoupled bases with multiphysics information [1][2].

**Acknowledgments**

We wish to thank Jean-Daniel Beley from Ansys France for his insights during the research leading to this study, as well as his help reviewing this publication.

**References**

[1] Tabak U and Rixen D J 2010 Reduced iterative correction algorithm for coupled vibroacoustic problems
*Proc. of ISMA 2010 - Int. Conf. on noise and vibration engineering* (Leuven) pp 4685-95

[2] Tabak U and Rixen D J 2011 Globally enriched substructuring techniques for vibro-acoustic simulation
*Linking Models and Experiments Proc. of the 29th IMAC* (Jacksonville) vol 2 ed T Proulx (New-York:
Springer) pp 263-80

[3] Maess M and Gaul L 2006 Substructuring and model reduction of pipe components interacting with acoustic
fluids *Mech. Syst. Signal. Process.* **20** 45-64

[4] ANSYS, Inc. 2023 ANSYS 2023R1 *Mechanical APDL Theory reference* pp 268-70

[5] Chen J-T, Huang W-S, Lee J-W and Tu Y-C 2014 A self-regularized approach for deriving the free-free
flexibility and stiffness matrices *Comput. Struct.* **145** 12-22

[6] R R Craig and M C C Bampton 1968 Coupling of substructures for dynamic analysis *AIAA J.* **6** 1313-19

[7] Kim S M, Kim J-G and Park K C 2019 A strongly coupled model reduction of vibro-acoustic interaction
*Comput. Methods Appl. Mech. Eng.* **347** 495-516

[8] Sigrist J-F, Broc D and Lainé C 2006 Dynamic analysis of a nuclear reactor with fluid-structure interaction
part I: seismic loading, fluid added mass and added stiffness effects *Nucl. Eng. Des.* **236** 2431-43

[9] Morand H J-P and Ohayon R J 1995 *Fluid-structure interaction: applied numerical methods* (Chichester
New York Brisbane: Wiley, Paris Milan Barcelona: Masson)

[10] Yao Y, Huang X, Yang X and Jézéquel L 2023 Symmetric formulations combined with component mode
synthesis for analyzing coupled fluid-structure problems. *Reviews and extensions Appl. Math. Model.* **115**
645-60

[11] Campanile L F, Kirmse S and Hasse A 2021 A measure for the similarity of vector spaces *Preprint*
10.13140/RG.2.2.30164.73606