Port-Hamiltonian Fluid-Structure Interaction Modeling and Structure-Preserving Model Order Reduction of a Classical Guitar

,


Introduction
The modeling and simulation of complex technical systems is indispensable for developing and operating new products.Different physical domains are suitably combined to create entirely new functionalities and obtain a more realistic model and a better understanding of the overall system.Considering multiple physical phenomena poses additional challenges in the system description as these systems need to be reasonably coupled.Their interaction must be taken into account for meaningful simulation results and a deeper understanding of the underlying system behavior.The port-Hamiltonian (pH) representation is based on the idea of interconnecting Hamiltonian systems and including dissipation in the formulation, for this reason the energy is used as the lingua franca at the coupling ports [38].This approach is attracting increasing attention due to its ability to describe complex multi-physics systems in the context of automated modeling and provides the basis  for a unified description of modeling, analysis, control, optimization, and model order reduction.Furthermore, under mild assumptions pH-systems implicitly ensure important system properties such as passivity, stability and energy conservation [13,31].
Fluid-structure interaction (FSI) characterizes the interaction of a moving or deformable mechanical structure with an internal or surrounding fluid [11].The consideration of FSI plays a crucial role in different applications, e.g. the bending of aircraft wings, the blood flow in large veins or the lubricant flow in ball-bearings [32].In the current study we want to analyze an acoustic system in which FSI is also of great importance for its proper technical function.
A classical guitar is a popular instrument consisting of various components whose most important ones are briefly described in order to understand the physical context of the mathematical problem, see Fig. 1.The musician plucks the guitar strings which thereupon start to vibrate.Those vibrations are transferred through the bridge into the soundboard which is the top plate of the guitar body.The guitar body coupled with the enclosed air acts as a resonator for the strings' vibration.The guitar's composition and material properties influence the volume and timbre of the air escaping through the sound hole.The classical guitar represents an impressive example of a multi-physics system where the guitar body is strongly coupled with the enclosed air and the coupling effects cannot be neglected for a realistic representation of the behavior.The modeling of a guitar without strings can be found e.g. in [8,14] and the modeling of strings in [7,15].Port-Hamiltonian representations of acoustic models are presented for instance in [34] for the vocal apparatus and for a Rhodes piano in [19].The modeling procedure of this work follows closely the approach presented in [8].
In the current study this FSI problem is derived from the partial differential equation of linear elasticity for the mechanical structure [35] and the wave equation for the fluid [6].The system is spatially discretized with the finite element (FE) method in order to obtain a coupled second-order ordinary differential equation (ODE) system [42].Structural damping obtained from measurements is included into the model [8].From this basis, a port-Hamiltonian representation is attained by first transforming the system into a first-order system and then applying an adapted state transformation.The coupling effects of the multi-domain system are illustrated in a transient time simulation which gives further insight to classical frequency-based analyses in acoustics.
The spatial discretization via FE methods typically leads to high-dimensional systems that are ill-suited for multi-query simulations, e.g. for usage in optimization loops.Hence, model order reduction (MOR) methods are essential to approximate the system in a subspace of a much lower dimension for a faster and more efficient simulation [3].One major novelty is the presentation of a comprehensive sensitivity analysis for a variety of combinations of projection methods and basis generation methods, i.e. the calculation of projection matrices, for reducing the FSI problem.The projection methods focus on basis-independent structure-preserving1 reduction [40,16] to preserve the important system properties in comparison to non-structure-preserving variants.The basis generation comprises approaches that are based on the system matrices, such as modal truncation [26] or moment-matching via Krylov subspaces [21], and data-based techniques based on timeresponse snapshots, e.g.Proper Orthogonal Decomposition (POD) [39] and Proper Symplectic Decomposition (PSD) [29,10].For more information on recent developments in structure-preserving MOR procedures for both Hamiltonian and port-Hamiltonian systems, we refer to [5,24,2,23] and the references therein.
The paper is organized as follows.In Section 2, the dynamic equations for a fluid-structure interaction problem are derived and adapted on second-order level.Additionally, practical aspects from the enhanced modeling process as an ensemble between the commercial software Abaqus and in-house code in Matlab are described.In Section 3, the necessary background on port-Hamiltonian systems is given and the FSI model of the guitar in a port-Hamiltonian input-output state representation is derived.Based on this, further pH formulations are deduced and the transient time simulation of the full system is performed.Various non-common projection methods, both structurepreserving and non-structure-preserving, are presented in Section 4. In Section 5, we recall classical and introduce modern basis generation methods with adapted snapshot matrices and non-orthogonal bases.The comprehensive sensitivity analysis is given in Section 6.Finally, we conclude the paper.An overview of the entire workflow is given in Fig. 2. In the remainder this workflow diagram will be referenced by a circled number which highlights the workflow area the current section is referring to, e.g. 1 .

Modeling of a classical guitar
When producing a tone with a guitar, many different guitar components interact with each other, with the interaction of the guitar body and enclosed air being of particular importance.An enhanced three mass model focuses on this interaction and serves as a comprehensive system on which the investigated methodologies are presented.The present chapter first describes the mathematical background of this fluid-structure interaction and then also addresses the practical aspects of its implementation using Abaqus and Matlab.

Modal Krylov
Projection structure-preserving non-preserving

Derivation of the dynamic equations for the fluid-structure interaction
The vibration of the mechanical structure can be described as an elastic motion of the wooden plates of the guitar that can be mathematically explained by the methods of continuum mechanics, cf. 1 .The motion of a body is described by the motion of all material points P S for this body, where the current position of each point can be classified by its position R 0 in the initial configuration and the displacement z(R 0 , t).Assuming small strains leads to the elastodynamic problem for the time interval strain-displacement relation: constitutive equation (Hooke's law): known as the strong formulation.The fundamental laws and material-dependence are considered in the force equilibrium (1a) and the constitutive equation (1e).The density of the material ρ S , the second temporal derivative of the displacements z, the divergence of the stress tensor σ and the volume forces b 0 form the force equilibrium (1a), where ∇ describes the nabla operator.The constitutive equation (1e) involves the engineering strain ε s and elastic modulus C EM .We restrict ourselves to the case of orthotropic, linear elastic material behavior which can be parametrized with of nine parameters: three Young's moduli, three Poisson ratios and three shear moduli to consider the behavior in each spatial direction [43].Γ z describe the parts of the body surface where Dirichlet boundary conditions for the displacement z b apply.Likewise Γ F SI describes the contact surface of the mechanical structure and the encased air on which the stresses are given by the surface traction t = pn F SI whose magnitude scales proportionally with the pressure p of the air at the contact surface.Here n F SI denotes the normal vector of unit length pointing into the direction of the encased air.The initial conditions at time t = 0 are given for the displacements z 0 and velocities ż0 [35,36,42].For a numerical simulation the continuum is spatially discretized with the finite element method (FEM) resulting in the symmetric mass matrix M S ∈ R N S ×N S and the symmetric stiffness matrix K S ∈ R N S ×N S .To model dissipation effects, i.e. the damping parameters, we employ the so-called Rayleigh damping with the coefficients α and β.Measurements of the damping ratio were conducted on a real guitar at different frequencies ω D , see [8], and the overdetermined equation system is solved by using a least square fit for the coefficients that results in α = 15.958 1 s and β = 2.821 • 10 −6 s.This complements the dynamic equation of the mechanical structure resulting in the secondorder space discretized system, c.f. 3 with the damping matrix D S ∈ R N S ×N S and the excitation forces f S , f p ∈ R N S for the mechanical structure, whereas the latter is the specific force stemming from the coupling term (1c).
The geometry and division of the guitar into the different domains can be seen in Fig. 3.The model consists of three different parts: the top plate or soundboard with a sound hole, the back plate and the air inside the guitar body in the shape and of the size of a classical guitar.Modeling the guitar behavior with these parts (in lumped parameter form) is known as a three mass model [14].The approach of [8] is used as a reference to create an enhanced three mass model, which will be discussed further in Section 2.2.top plate (solid) enclosed air (fluid) back plate (solid) The fluid behavior, c.f. 1 , can be described with acoustic wave equation: boundary conditions: where the acoustic wave equation is derived from the Euler equations for ideal, compressible fluids with Dirichlet boundary conditions p b on Γ p , zero Neumann boundary conditions, which model natural boundary conditions at the parts where the fluid is not encased by the structure, i.e.Γ n and non-zero Neumann boundary conditions at the interface, which describe the coupling relationship between the fluid and the mechanical structure.The force applied to the mechanical structure along the normal direction, i.e.
∂p ∂n F SI is equal in magnitude but of opposite direction compared to the force applied by the mechanical structure along the normal direction, i.e. −ρ F ∂ 2 z ∂t 2 T n F SI .The acoustic wave equation (4a) contains the speed of sound c 0 and the pressure p on the domain P F [6,27].
The wave equation is spatially discretized by using the Galerkin method with linear shape functions.This leads, c.f. 3 , to the dynamic equations of the fluid in matrix form [25] with the symmetric mass matrix M F ∈ R N F ×N F , the symmetric stiffness matrix K F ∈ R N F ×N F , the excitation vector for the fluid f F ∈ R N F and the excitation vector f z ∈ R N F stemming from the interaction with the mechanical structure, i.e. the coupling term (4d).FSI is an area that describes the multi-physics coupling between the domains of fluid dynamics and structural mechanics.The coupling is considered a strong FSI if the motion of the solid influences the fluid and vice versa [27].The classical guitar as the studied system allows the investigation of these strong couplings.In our case this coupling is expressed via the coupling terms (1c) and (4d).
The effects of the mechanical structure on the fluid and vice versa, c.f. 4 , result in the external forces f p and f z , respectively.The former relates to the quantity whereas the latter corresponds to the quantity Hence, both f p and f z can be represented as a matrix vector product where R describes the coupling matrix [25] representing the term However, in the present case of non-matching meshes, the nodal values of the acoustic and structural domains, which are part of the interface, are related to each other with linear interpolation.Consequently, the separate equations for the mechanical structure (3) and the fluid ( 5) can be combined into the following coupled system: Note that the mass and stiffness matrices, M and K, are unsymmetric which leads to numerical issues when solving large sparse linear systems and to problems for the derivation of a port-Hamiltonian formulation of the system.Therefore, an adapted formulation of the system is derived in the following.
Instead of using the z-p-formulation (10) a state transformation is applied, where instead of the pressure p, the time derivative of the variable q is used which is in close relationship to the velocity potential.This approach is frequently used in the literature to symmetrize the system matrices [37,17].
A minor adaptation is used to make the pH matrices accomplish the mandatory skew-symmetry property after the transformation from a first-order to a second-order system: Inserting (11) into (10), taking the integral of the second line of (10) and dividing by ρ F (instead of −ρ F which is typically considered in literature) leads to the formulation with

Practical modeling aspects
The FE discretization is realized with the commercial software Abaqus, see [1], which allows us to design the complex geometry of the guitar, define material parameters and boundary conditions, and automatically mesh the geometry with sophisticated meshing tools, see 2 .
The main goal of the current study is not a very accurate prediction of the guitar's behavior in accordance with experimental data but to rather investigate the performance of various formulations on an illustrative example for FSI.For this reason, the guitar is not modeled in full detail but still includes the most important aspects.
The top and back plate are modeled with shell elements S4R5, see Fig. 4. A summary of the finite element mesh properties used for the guitar model is presented in Tab. 1.We use S4R5 element as it is accompanied by the advantage that the system size is reduced since there are less DOFs compared to the more common six DOFs S4 element.Furthermore, for the mechanical structure a six degree of freedom (DOF) formulation in Abaqus would lead to a singular matrix M S which would cause additional numerical issues and would lead to a differential-algebraic pH-system.The air model is enhanced by adding length correction to the fluid domain to account for effects of the surrounding air [18], see Fig. 5.As stated in the mathematical derivation of the mechanical structure and the fluid, some boundary conditions need to be considered that are in accordance with the physical behavior of the guitar.Therefore, the mechanical structure is constrained at the edges with homogeneous Dirichlet boundary conditions for the displacement variables.The rotational DOFs are not affected by this condition.This kind of boundary condition is considered as simply supported.The fluid is mainly characterized by the coupling with the solid and the reflection from the sidewalls, due to homogeneous Neumann boundary conditions on Γ n .However, the boundary condition at the top of the sound hole needs to be defined different since there is a free surface.This can be achieved by simply using the assumption of zero pressure [43].The boundary conditions are marked in red and purple in Figs. 4 and   The outcome of the modeling process in Abaqus are the mass and stiffness matrices for the mechanical structure (3) and the fluid (5).The further steps of modeling and simulation are performed in Matlab.This change in the software tool is motivated by various advantages: (a) The modeling in Matlab gives a better insight into the model especially due to the independent coupling procedure.In addition, (b) the obtained flexibility enables the system to be modified in the first place and thus to be reformulated into a pH-system.Furthermore, (c) this choice also gives us the freedom to decide on the integrator and the associated simulation properties.
In a real scenario, the guitarist plucks the string which introduces vibrations that are transferred through the bridge into the soundboard.Since the string and bridge are not modeled, the influence of these components is integrated as a force input at the position where the bridge would normally be located.A standard tuning of a guitar will create sounds in the range of approximately 82 Hz on the lowest pitch and 320 Hz at the highest pitch for the unfretted strings.The frequency range is chosen this way, since the objective is to demonstrate the methods on unfretted strings and the harmonics are neglected.For this reason, there will be a sine wave excitation force with with ω = 2πf in the frequency range f = [82, 320] Hz and with an amplitude of û = 1 N that acts on one node of the top plate.As outputs, different important nodes of the mesh are considered.Those nodes include the aforementioned excitation node where the force acts on, compare Fig. 4, as well as the volume-weighted integral of the pressure over the fluid nodes close to the sound hole of the guitar, compare Fig. 5, and a central mechanical structure node of the back plate.

Port-Hamiltonian formulation
Port-Hamiltonian (pH) systems were initially developed to describe a unified approach for systems from different physical domains that allow network modeling, see references in [9,31].The pH formulation describes an energy-based network modeling approach with the energy as common quantity of various physical systems that interact through ports [38].
The pH structure can imply the underlying physical principles such as conservation laws [4].The input-output state form of a pH-system with neglected feed-through terms appears as where the function H(x(t)) describes the Hamiltonian as an energy function of the system.The Hamiltonian depends on the state vector x(t) ∈ R N with the system order N = 2 N .The matrix J ∈ R N ×N reflects the interconnection of the internal energy storage elements while the dissipation matrix D ∈ R N ×N describes the energy losses in the system.The matrix B ∈ R N ×m is the port or input matrix that specifies the manner in which energy enters or leaves the system.The Hamiltonian and the port matrix motivate the term pH-system.The pH-system matrices need to satisfy the requirement that J = −J T is skew-symmetric and D = D T ≥ 0 is symmetric positive semidefinite.Furthermore, the variables u(t), y(t) ∈ R m describe the time-dependent input and output vectors, respectively [38].
Besides the modeling aspects, the pH formulations come with additional advantages as they implicitly incorporate important system properties.The pH-system ( 15) is a generalization of a classical Hamiltonian system as the conservation of energy is generalized to a dissipation inequality which in combination with the assumption that the Hamiltonian is quadratic and strictly positive H(x) > 0 for all x leads to the properties that the pH-system is both passive and stable [13].Note, that the input and output entries of ( 15) belong to the same variables defined by the port matrix B.
One can choose arbitrary quantities of interest (QoI) as outputs with appropriately chosen matrices but these output quantities may not satisfy the dissipation inequality (16) and only serve the purpose of data observation.Depending on the specific application and the coordinate system for which the problem is formulated, a generalized linear port-Hamiltonian system in the so-called descriptor form [4], with constant matrices E ∈ R N ×N and Q ∈ R N ×N that satisfy the condition E T Q = Q T E in addition to the previously mentioned J , D and B may also be considered.
In this case, the Hamiltonian of the system can be described as a quadratic energy function For our setting, such a descriptor system is obtained by transforming the second-order system (12) into first-order a pH-system (17) via the new state vector x = z T q T żT qT T , which leads to the pH-system matrices with E describing the kinetic energy components and Q containing potential energy terms, c.f. 5 .
The matrices I N F/S represent the identity matrix of size N F/S .The pH matrices satisfy the mandatory properties The pH-system (19) will further be considered as the velocity formulation.

Numerical simulation of full-order results
We discretize the system with respect to time, using the implicit midpoint rule (IMR), cf.6 , which is the Gauss-Legendre collocation method with 1 stage and of order 2 for ODEs.The IMR is a symplectic integrator and hence, preserves the symplectic structure and the quadratic Hamiltonian through the time integration [28].
We use a step size of ∆t = 10 −4 s and simulate for T = 0.1 s to include several oscillation periods for the whole frequency range.All calculations were carried out with Matlab.An advantage of a fixed time step is the better comparability of the trajectories.Since the numerical error of an integrator depends on the step size, we avoid mixing the integration and model reduction error.
The simulation of a full-order model (FOM) of the guitar is important for different reasons.On the one hand, these results will give us an insight about the approximate guitar behavior and the physical plausibility of the modeling.On the other hand, the FOM results form the basis for the calculation of the approximation error in the energy norm.Additionally, the FOM results are required for the snapshot generation, cf.7 , for the data-based basis generation techniques.
We consider a parameter dependent excitation u = u(µ) of the top plate that is modeled by a sinusoidal input where µ ∈ [82, 320] Hz denotes the excitation frequency f .The time simulation of the FOM for an excitation frequency of f = 100 Hz is illustrated in Fig. 6.It can be seen that as a result of the excitation, the top plate at the excitation node starts to vibrate at the same frequency since the guitar describes a linear time-invariant (LTI) system.The vibrations are transferred into the fluid whose pressure values start oscillating which corresponds to a sound of the guitar.The back plate vibrates at a lower magnitude than the top plate which is justifiable by the transmission of energy through the fluid.It shall be emphasized that the back plate is only connected to the excited top plate via the fluid and hence, the strong coupling effects of the fluid can be demonstrated impressively.

Reformulations of the pH-system
We derived our pH-FSI-system via a transformation of the second-order system and using the displacements, the integrated pressure and their derivatives as pH variables.In classical Hamiltonian systems the variables are defined in terms of the canonical position and momentum coordinates to account for additional symmetries [13].For this reason, a system in momentum formulation shall also be considered in the analysis.All of the following reformulations of the pH-system still satisfy the mandatory pH properties (pH1) to (pH3).The coordinate transformation yields such a pH-system in momentum formulation where x m contains the momentum instead of the velocities.
Since the PSD methods presented in Section 6 are based on a system formulation with canonical J , i.e.
the canonical variants of the velocity formulation and momentum formulation will be considered in the sensitivity analysis, which yield the canonical system in velocity formulation and in momentum formulation The transformation matrices P , P −1 depend on the coupling matrix R and are given by The calculation of the FOM is computationally too expensive, especially if it comes to multiquery simulations.Therefore, MOR techniques are required and will be discussed in the following.

Model order reduction -Projection methods
The size of a dynamical system model can be reduced by approximating its behavior in a subspace of a lower dimension.This model reduction comes at the cost of approximation errors.Projectionbased reduction methods can be performed in various ways with different outcomes with respect to the preservation of different structural2 pH properties, cf. 10 .
All reduction methods used in our study are based on a projection-based reduction.The solution for x(t) is approximated in a subspace V of dimension n N which is described by a basis matrix V ∈ R N ×n with colsp(V ) = V and leads to the approximation with the reduced state vector x r ∈ R n .Inserting (28) into the pH-system (17) and using the Petrov-Galerkin condition for W ∈ R N ×n [3] yields the reduced system of size n N , where the matrix W determines the orthogonal projection direction.A special case of the Petrov-Galerkin approach is the Galerkin projection with W = V and takes the form This projection does in general not preserve any of the underlying pH structure properties (pH1)-(pH3) of the system in our formulation.
The quasi-Galerkin projection is derived by inserting the term V V T under the assumption that (J − D)V V T Q ≈ (J − D)Q to the standard Galerkin projection which yields This projection preserves (pH1) and (pH2) but it does not guarantee (pH3) and hence, the Hamiltonian of the reduced system H(x r ) = H(V x r ) changes in the reduction.
Adapting the pH structure-preserving approach [40] with W = QV to the case of a descriptor system leads to a reduced system with the reduced matrices This preserves the pH properties (pH1)-(pH3) and does not change the Hamiltonian in the reduced system [40].The author in [16] presents a general framework for the numerical approximation of evolution problems which is also suitable for pH-systems and preserves the underlying Hamiltonian structure.We will call this approach energy-stable.Applying this approach to a pH-system leads to a reduced system In our case the inverse term can be easily calculated with which allows for a computational efficient implementation.In order to assess the quality of the different projection methods, a direct comparison with the best approximation in the considered subspace colsp(V ) should be examined.Here, we use the best approximation with respect to the energy norm • 2 H induced by the energy inner product where H = E T Q denotes the so-called energy matrix which is positive definite, since both E T and Q are positive definite and their product is symmetric due to the third pH-property (pH3).The best approximation of x in the subspace colsp(V ) is then given by the projection Π V ,H of x onto colsp(V ) with respect to the energy inner product.In terms of the energy matrix this can be expressed via and we have x − x H (38)

Model order reduction -Basis Generation
For all projection-based reduction methods presented in Section 4, a basis matrix V is required.
In the following, different basis generation approaches are presented with the goal to keep the approximation error as small as possible.The bases are generated by using the extensive software packages MatMorembs3 [20] and RBmatlab4 , cf. 8 and 9 .
For the basis generation method based upon modal truncation, the homogeneous solution for (17) with u(t) = 0 is solved with the ansatz function which leads to the generalized eigenvalue problem with A := (J − D)Q for the pH-system.The eigenvectors φ i are also called coupled eigenmodes and describe the deformation of the mechanical structure and pressure distribution in the fluid at the dedicated eigenfrequencies.Modal truncation uses only the most important eigenvectors as the projection basis [3].In the context of the guitar, the most important eigenvectors coincide with the eigenvectors that belong to the lowest eigenfrequencies since those are the crucial eigenfrequencies for the sound emission.Hence, the modal projection basis arises as with n N .In the current study, the eigenmodes belonging to the lowest eigenfrequencies are coupled eigenmodes that allow for dynamics in the mechanical structure and fluid simultaneously.In general, V mod is not real-valued.In this case, the matrix V T mod has to be replaced with V H mod for the projection methods outlined in Section 4, where the subscript H represents the conjugate transpose of a complex-valued matrix.
The goal of Krylov-based reduction is the approximation of the Laplace transformed transfer function of the pH-system with the complex variable s ∈ C which is typically evaluated on the imaginary axis s = iω with the circular frequency ω = 2πf where f denotes the excitation frequency.
The transfer function around an expansion point s 0 can be described with a Taylor series whose coefficients are called moments in this context.The first J b moments of the full and reduced system can be matched by using the block Krylov subspace of a pH-system The approach can be generalized to calculate for multiple expansion points [26].
The direct calculation leads to numerical issues since additional vectors can become linearly dependent.That is why the block Arnoldi algorithm is used for the calculation of Krylov subspaces.This Arnoldi algorithm consists of a LU decomposition, Gram-Schmidt orthogonalization and the iterative calculation of Krylov directions [26].Sometimes the approach is also called tangential interpolation [12].
The Proper Orthogonal Decomposition (POD) approach uses a different idea for the projection basis generation which is based on snapshots instead of system matrices and can therefore even be used for nonlinear systems [39,30].The POD starts with a set of vectors xi ∈ R N with i = 1, . . ., m assembled column-wise into a matrix The goal is to approximate the information contained in the snapshot matrix X by a set of vectors ûi ∈ R N with i = 1, . . ., n, which can be expressed in terms of an optimization problem min û1 ,..., ûn where δ ij describes the Kronecker delta and hence, the vectors ûi form an orthonormal basis [39].
The optimization problem can be solved by solving the eigenvalue problem with the scaled correlation matrix X XT ∈ R N ×N which is computationally inefficient to calculate if m < N due to its size.Therefore, the method of snapshots solves the eigenvalue problem of size m × m with the eigenvectors vi which saves computational time for m N [39,22].The problem (46) then is solved by By only taking the first n POD-modes, the basis matrix can be formed as In the present system, the snapshot matrix X consists of state variables of the pH-system x(t, µ) with µ ∈ [82, 320] Hz.The discrete trajectories that form the snapshot matrix are calculated from the full system with time instances of ∆t = 10 −4 s for a time span of T = 0.1 s.Using the full state vector x = z T q T żT qT T for calculating, the POD basis will further be called as V POD,State .
Another basis will be generated by only using the displacements in the snapshot matrix XDisp = xDisp,1 . . .xDisp,m of the adapted state vector xT Disp,i = z T i q T i and assembling the full projection matrix afterwards as Note that the above is a special case of the so-called tangent lift as presented in [29], where only information in the state is considered.Hence, using the basis for the velocity component might be ill-suited.
A further approach will be the division of the state vector into the four individual parts and calculating the POD-basis for each individual trajectory and assembling the basis as the blockdiagonal matrix with the snapshot matrices Xz , Xq , X ż and X q corresponding to the respective components of the state x.Symplectic MOR is structure-preserving MOR for Hamiltonian systems.One approach to derive a symplectic reduced-order basis is the data-driven Proper Symplectic Decomposition (PSD) which is closely related to the POD [29].
Assuming a suitable basis, a symplectic form ω 2n : R 2n × R 2n → R takes the canonical form where J 2n is the Poisson matrix with the identity matrix I n ∈ R n×n and the matrix of all zeros 0 n ∈ R n×n [33].We call V ∈ R 2 N ×2n a symplectic matrix if If a Petrov-Galerkin projection is used with W T = V + where denotes the so-called symplectic inverse5 , then the Hamiltonian structure and therefore its Hamiltonian is preserved for purely Hamiltonian systems, i.e. in the case of D = 0 in equation ( 15) [29,10].
In analogy to the POD, the minimization problem appears as minimize where the constraint guarantees that the reduced-order basis (ROB) is symplectic [29].In contrast to POD, there is no explicit solution procedure for the PSD optimization problem (55).The PSD Complex SVD approach is the solution of the PSD in the subset of symplectic, orthonormal ROBs [10].It uses an adapted complex snapshot matrix with the imaginary unit i.The minimization problem requires an auxiliary complex matrix U Cs ∈ C N ×n and takes the form minimize and generates the real basis matrix The specialty of the PSD Complex SVD is the choice of the auxiliary complex matrix C s and the computation of E from U Cs .The solution of (57) is the POD for complex matrices based on the left-singular vectors of C s that can be computed with a complex version of the SVD [29].
Most existing basis generation techniques, e.g.Complex SVD, generate a symplectic, orthonormal ROB.In [10] a new symplectic, non-orthogonal basis generation technique based on the so-called SVD-like decomposition is introduced.
There exists a SVD-like decomposition of the snapshot matrix with a symplectic matrix S s ∈ R 2 N ×2 N , a sparse and potentially non-diagonal matrix D s ∈ R 2 N ×m , an orthogonal matrix Q s ∈ R m×m and the symplectic singular values σ i .The rank of the snapshot matrix is rank( X) = 2p + q.
The Frobenius norm of the snapshot matrix can be rewritten as where s i is the i-th column of S s and w i is called the weighted symplectic singular value [10].
The goal is to choose the k indices i ∈ I SVD = {i 1 , . . ., i k } ⊂ {1, . . ., p + q} which have large contributions w i to the Frobenius norm with which yields a ROB V SVD ∈ R 2 N ×2k with k pairs of columns s i ∈ R 2 N from S s leading to reduced model order of 2k = n and the basis matrix The SVD-like decomposition is constructed by computing an eigendecomposition of XT J X, for which the imaginary and real part of eigenvectors corresponding to complex conjugate eigenvalues result in a pair of symplectic vectors.Alternatively, the approach introduced in [41] can be used.However, we chose the former method, as the later one can be quite computationally demanding.

Results
The previous sections featured several alternatives for model order reduction of a classical guitar FSI problem in a pH formulation.The different components that were considered are the 4 system formulations in Sec.3.2, the 4 projection methods from Sec. 4 and the 7 basis generation methods in Sec. 5, cf.Fig. 2.These categories are extended by choosing the reduced system size n, studying 11 sizes between n = 12 and n = 400 and the number of trajectories in the snapshot matrix for the data-based methods, where 6 alternatives were investigated, see Tab. 2. All representatives of each category can be combined with all others, resulting in a number of about 4 • 4 • 7 • 11 • 6 ≈ 7000 possible combinations.All of these combinations have been computed.We will give a focused view onto the most important outcomes in the following sensitivity analysis.

Sensitivity analysis of different reduced-order models
In order to compare the different methods for reducing the FOM to a reduced-order model (ROM), cf.11 , an error measure is needed.The energy norm of the Hamiltonian is used to describe the energy value of a trajectory at time step t.The error ε(t) = x(t)−V x r (t) H between the full and reduced model is measured in the energy norm.In order to obtain a scalar and relative measure, the norm is integrated over the simulation time interval I t = [0, t end ] and related to the full model which yields for the relative error measure.It is important to use the relative error which allows for a comparison between different system inputs.Five randomly chosen frequencies in the standard tuning range (14) and zero initial conditions were used as system parameters for all experiments that were carried out.The mean value of their relative error calculations is displayed in the subsequent results.
The relative error values for the four system formulations over all basis generation techniques and for the pH and energy-stable projections are illustrated in Fig. 7.It can be seen that the displacement variant of POD is not a suitable option in the momentum formulation.This is due to the fact that the state derivatives are not considered in the snapshot matrix and hence, the inertia information gets lost.The canonical alternative for the velocity formulation (25) generally leads to much higher errors and also in the momentum formulation the transformation ( 26) is not increasing the performance.Furthermore, in the transformation process, the block structure of the pH matrices E and Q (19) gets lost.The momentum formulation for both, canonical and non-canonical, leads to the best overall results especially for POD-State and SVD-like.For these reasons, we will focus on the non-canonical momentum formulation in the remainder.
The three best-performing methods from the previous experiment, i.e.POD-State, SVD-like and POD-Individual, are compared more closely in Fig. 8.The plot shows the development of the error of the three methods over the basis size and examines on the one hand the pH and energystable projections, and on the other hand their behavior in relation to the best approximation with respect to the energy norm, i.e. the FOM solution is projected onto the subspace spanned by the corresponding basis matrix respecting the inner product pertaining to the energy norm (38).One can see that the results of pH and energy-stable behave in a similar way for POD-State and SVDlike but for POD-Individual, the energy-stable projection leads to much better error values.The deviation concerning the best approximation does not exceed one order of magnitude for all basis generation methods.The SVD-like method remains the closest to its potential best values for smaller reduced sizes with increasing distance for higher basis size.This is due to numerical reasons as the involved computation of the symplectic SVD-like decomposition is quite sensitive with respect to the condition of the snapshot matrix.A comparison of the basis generation methods uncovers that the SVD-like technique lies in the range of the best approximation error of the other methods and shows an improvement of approximately one order of magnitude for several reduction dimensions compared to the competitors.The overall behavior of all basis generation techniques (even the one not shown in the plot) shows that the errors decay rapidly up to a size of n = 120 and stagnate starting from a size of approximately n = 200.
The heat map in Fig. 9 shows the relative error for all combinations of basis generation techniques and projection methods that were conducted, cf. 13 .The errors are shown for a basis size of n = 120.If the relative error exceeds a value of ε r > 1 then the field entry is left blank.One can see that the Galerkin and Quasi-Galerkin projection perform the worst.These projections do not preserve the pH properties, especially the stability conditions, leading to instability and poor approximations.These projection methods are hence not suitable for the reduction of pH-systems.In general, the pH and energy-stable projection methods lead to the smallest error because of their adaption to the specific underlying pH structure of the system, which shows the importance of adapting the projection methods to the particular structure of the problem.The modal reduction, which is still used in various commercial MOR-packages, leads to high errors compared to the other combinations.The dynamics of the model cannot be described sufficiently accurate by only allowing to move in modal coordinates that belong to the lowest eigenfrequencies.The basis generated with the C-SVD and Krylov algorithms lead to moderate approximation errors which, in the case of the C-SVD, also cannot be improved by the transformation into the canonical structure of the matrix J .The Krylov algorithm focuses on approximating the input-output behavior at certain frequencies, leading to higher error at different frequencies.The overall performance of the remaining data-based methods worked the best for the FSI dynamics of the guitar except for the POD-Displacement variant.The reason for the failure of the POD-Displacement has already been explained above.For the POD-State variant, the information contained in the snapshot matrix and extracted from this matrix via POD is sufficient for reproducing the coupled dynamics of the system.On the other hand, it can also be advantageous to consider the coordinates separately in the basis.The POD-Individual method outperforms the full-state alternative especially for the energy-stable projection.This stems from the fact that the POD-Individual approach accounts for different scalings in the respective components.The SVD-like basis generation technique shows the best overall performance with error values as low as ε r = 3.2 • 10 −3 .The explanation for this lies in the structure of the basis.The SVD-like approach is the only method that builds a non-orthogonal basis and can therefore adapt more flexibly to the system's dynamics.Special attention should be paid here to the fact that the SVD-like method was developed for purely Hamiltonian systems, whereas this study shows that it is also suitable for port-Hamiltonian systems.It is worth mentioning that the information content in the snapshot matrices depends on the number of underlying trajectories which are randomly distributed in the parameter space.In the following study, trajectory numbers between five and 100 are investigated.The results of the associated relative errors are shown in Tab. 2 and are extended by the methods based on system matrices for comparison.It can be seen that the influence of a high number of trajectories is not as significant as for instance the size of the reduced system, see Fig. 8.All methods extract enough information for basis constructions even for small snapshot sizes.A snapshot matrix made of ten trajectories already shows a good trade-off between approximation error and matrix size, which drastically decreases the computational effort in the offline phase due to few evaluations of the full system.
The main goal of MOR is the gain of a speed-up compared to the simulation of the full-order model to make this model suitable for multi-query tasks or real-time scenarios.The speed-up values are determined for the pH projection method and averaged over all basis generation methods since the qualitative behaviour for the other methods was similar.The results are listed in Tab. 3. The greatest speed-up of 540 is obtained with a basis size of n = 12.This choice is not recommendable since the error values are very high for this size.The speed-up values decrease down to 45 for a basis size of n = 400.Depending on the requirements for the ROM, one needs to decide which is the best trade-off between approximation quality and gained speed-up.Since the error decays fast up to a size of n = 120 and still leads to speed-ups of more than 200, this reduction size was focused on in the previous discussion.In Fig. 10 and Fig. 11, the displacement in the excitation node and the averaged pressure in the sound hole, cf.Fig. 5, are depicted.We can see in both cases, that all basis construction methods except for POD-Displacement capture the behaviour quite well.For the POD-Displacement, however, almost no energy seems to be introduced into the system resulting in no visible excitation in the entire FSI model.

Conclusion and Outlook
In this paper, we presented a port-Hamiltonian formulation of a fluid-structure interaction for the case of a classical guitar.Many different model order reduction approaches combined with multiple different reduced basis construction methods as well as various FOM formulations were compared in order to study the effect of structure-preserving model reduction on the quality and behaviour of the reduced model.In particular, we were able to conclude that both structure-preserving model order reduction methods, i.e. the pH-preserving and the energy-stable methods clearly outperform their non-structure-preserving competitors.Furthermore, we can conclude that the reduced bases constructed via a symplectic SVD-like decomposition result in higher quality approximations for lower basis sizes.Future emphasis will be placed on constructing suitable a-posteriori error estimators and extending the above model to allow for further parameter dependency in form of material parameters of the guitar body.

Figure 3 :
Figure 3: sectional view of the FEM multi-physics model of a classical guitar

Figure 6 :
Figure 6: time simulation for one selected node per component at excitation frequency f = 100 Hz

Figure 9 :
Figure 9: relative errors for different basis and projection combinations for n = 120

Figure 10 :Figure 11 :
Figure 10: time simulation of the excitation node displacement for the full and reduced-order models for n = 120 in the momentum formulation and pH projection

Table 2 :
Relative error of various basis generation methods over the number of trajectories from the snapshot generation for the momentum formulation and n = 120.One trajectory refers to an excitation with one specific frequency.

Table 3 :
Averaged speed-up values for different basis sizes for the pH projection