Mathematical Problems in Engineering

Volume 2013 (2013), Article ID 583079, 8 pages

http://dx.doi.org/10.1155/2013/583079

Research Article

Low-Frequency Acoustic-Structure Analysis Using Coupled FEM-BEM Method

Jinlong Feng,¹ Xiaoping Zheng,¹ Haitao Wang,² HongTao Wang,² Yuanjie Zou,³ Yinghua Liu,¹ and Zhenhan Yao¹

¹School of Aerospace, Tsinghua University, Beijing 100084, China

²Institute of Nuclear and New Energy Technology, Tsinghua University, Beijing 100084, China

³China Academy of Space Technology, Beijing 100094, China

Received 19 July 2013; Accepted 5 September 2013

Academic Editor: Song Cen

Copyright © 2013 Jinlong Feng et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

A numerical algorithm based on finite element method (FEM) and boundary element method (BEM) is proposed for the analyses of acoustic-structure coupled problems. By this algorithm, the structural domain and the acoustic domain are modelled by FEM and BEM, respectively, which are coupled with each other through the consideration of the appropriate compatibility and equilibrium conditions on the interface of the two domains. To improve the computational efficiency, the adaptive cross approximation (ACA) approach is incorporated into the proposed algorithm to deal with the nonsymmetric and fully populated matrices resulting from the coupling of the FEM and BEM. The validity and the high efficiency of the present approach are demonstrated by two examples.

1. Introduction

Acoustic-structure coupling problems can be found in various engineering fields. For example, in space exploration, spacecrafts are subject to heavy acoustic load, particularly during launching, which can impose severe and adverse effect on the astronauts and the structure of the spacecrafts. Therefore, it is necessary to develop numerical methods to predict the structural behaviour in the acoustic-structure coupling problems since these numerical techniques are essential to the improvement of practical engineering design.

The finite element method (FEM) and boundary element method (BEM) are both widely used in various engineering problems. However, FEM and BEM are generally applied in different areas. For instance, the former is especially well suited to the analyses of the problems which involve inhomogeneity or nonlinearity [1, 2], whilst the latter is one of the most powerful and suitable numerical techniques for predicting the noise around a vibrating system [3]. Compared to FEM [4], BEM reduces the mathematical dimension of the problem under analysis by one. Furthermore, in BEM, the Sommerfeld condition is automatically satisfied so that the external domain does not need to be bounded. As a result, BEM is regarded as the most suitable tool to deal with the problems in exterior unbounded domains.

Therefore, for acoustic-structure problems involving submerged bodies, it should be a natural strategy to combine FEM and BEM in a computational model to make use of the respective advantages of the two methods. The combination of FEM and BEM can be implemented by decomposing the concerned domain into several subdomains, each of which, according to its specific physics, is modelled by either FEM or BEM [5–7]. Because of the appropriate compatibility and equilibrium conditions on the respective interface boundaries, the subdomains are coupled to each other.

It should be mentioned that the boundary integral equation has a major defect for exterior problems; that is, it has nonunique solution at a set of fictitious eigenfrequencies associated with the resonant frequencies. In the past, several methods have been proposed to overcome this difficulty. The two most popular methods are the CHIEF method and Burton-Miller approach. The CHIEF method [8] is based on the Helmholtz integral equation on the surface of the radiating body combined with Helmholtz equations for some interior points. The resulting system equation is usually solved by the least-square approach. The Burton-Miller approach [9], inspired by Panich formulation [10], forms a linear combination of the Helmholtz boundary integral equation and its normal derivative, providing also a valid solution at any frequency but leading to hypersingular integrals which can be handled by Guiggiani algorithm [11, 12]. In this study, the Burton-Miller approach is adopted.

Similar to BEM, the coupled FEM-BEM approach leads to a fully populated system matrix. This is thought to be a major drawback of the coupled approach as compared with FEM which normally results in a sparse symmetrical system matrix. If iterative methods like the generalized minimal residual (GMRES) method [13] in combination with suitable preconditioners are used to solve the FEM-BEM system with degrees of freedom, the computational expense is in the order of . If a direct solver is applied, the computational cost is even in the order of . Fortunately, over the last few years, several fast solution methods including fast multipole method (FMM) [14–16] and adaptive cross approximation approach [17, 18] have been developed. These methods reduce the memory requirements and the computational time significantly when dealing with the matrix-vector product. In this study, the ACA approach is adopted since it is much easier to operate compared with FMM which is based on a series expansion of kernel-shape function products.

In this paper, a coupling algorithm combining FEM and BEM is developed for the analysis of acoustic-structure response. This paper is organized as follows. Section 1 is a short introduction to the background and corresponding research status. In Section 2, the FE formulations of the structure domain are constructed. In Section 3, the application of the BEM to Helmholtz problem and its implementation is presented. Then, the coupled equations are built in Section 4. In Section 5, the solving algorithms of coupled equations including the GMRES iterative solver and ACA approach are introduced. In Section 6, several numerical examples are presented.

2. FE Formulation of Structural Domain

A coupling system of acoustic-structure is shown in Figure 1. The structure is assumed to be fully submerged in interior and exterior acoustic fields or one of them. In Figure 1, the structure domain is denoted by . The interior and the exterior acoustic fields are denoted by and , respectively. The inner interface and the outer interface are named as and , respectively. The structure is assumed to be thin and flexible, which is susceptible to acoustic-structure interaction.

Figure 1: Illustration of a coupling system of acoustic-structure.

FEM is predominantly chosen for simulating linear elastodynamic systems. When the structure is discretized by FEM, the resulting system of linear equations reads as

where , , and denote the global mass matrix, the damping matrix, and the stiffness matrix, respectively. The vector is load vector, while the vector resulted from the tractions on the coupling interface, and it will be further discussed in the following section. The global displacement vector is denoted by .

In this paper, Rayleigh damping is considered for the damping matrix, which means that

where and are the damping parameters.

If the load vector and the acoustic pressure vector are both assumed to be harmonic, the displacement vector should be also harmonic. This means that

where the imaginary unit is denoted by and is the circular frequency with the excitation frequency . For the relationship between and , one obtains

Substituting (4) into (1) yields

NASTRAN, a finite element package, is employed to set up , , , and the right-hand side vector . Since , , and are frequency independent, these matrices only have to be calculated once for a given model.

3. The Boundary Element Method for the Helmholtz Equation

Acoustic problems in frequency domain can be described by the Helmholtz equation which has a three-dimensional form given by

where is the Laplace operator, is the sound pressure at point in the fluid, is the acoustic wavenumber, and is the speed of sound in the fluid.

The boundary condition can be expressed as

Additionally, the three-dimensional Sommerfeld radiation condition [19] is

For the exterior problem, (8) has to be satisfied to ensure that the wave is purely outgoing.

A three-dimensional fundamental solution for the Helmholtz equation is given by

where is a field point, is a source point, and is the distance from to . The key idea is to use Green’s second identity in combination with the property of the Dirac distribution:

By this way, the pressure at an arbitrary point within the acoustic domain is given by the integral representation

For scattering problems (11) can be written as

Let the point approach the boundary; then the following conventional boundary integral equation (CBIE) [20] is obtained:

where   if is smooth around . This CBIE can be employed to solve the unknown and on .

It is well known that the CBIE has a major defect for exterior problems; that is, it has nonunique solution at a set of fictitious eigenfrequencies associated with the resonant frequencies [21]. A remedy to this problem is to use the normal derivative BIE in conjunction with the CBIE. Taking the derivative of the integral representation given in (13) with respect to the normal at a typical point on and letting approach , one obtains the following hypersingular boundary integral equation (HBIE):

where if is smooth around . It should be noted that (14) suffers from the same kind of defect as (13). However, if a linear combination of CBIE and HBIE is used, the uniqueness of the results can be ensured for exterior acoustic wave problems. Therefore,

where is the coupling constant. This formulation is called Burton-Miller formulation [9] for acoustic wave problems and has been shown to yield unique solutions at all frequencies, if is a complex number which, for example, can be chosen as with [21]. In (15), , and and denote the density of the fluid and the normal velocity on the interface, respectively.

When the boundary of the acoustic domain is discretized into boundary elements, the resulting linear system equations can be expressed as

where the matrices and are obtained by integrating the fundamental solutions over each boundary element.

Equation (16) is the boundary element equation used to solve the acoustical problem. When there is incident sound, a free term is needed, and (16) becomes

4. Coupled FE-BE Formulations

By taking into account the appropriate compatibility and equilibrium conditions at the respective interface boundaries, the fully coupled FE-BE formulations can be derived. Hereby, the matching grids for the FE and the BE parts are required. The compatibility condition at the respective interface boundaries links the acoustic pressure and the tractions on the structure, which can be written in the following matrix form:

For an element, one has

where .

Substituting (19) into (5) and using uniform symbol for the matrices and vectors yield

From (20), one obtains

Please note that in (21) incorporates the velocity vector of all of the nodes in the structure, and in (17) just denotes the velocity vector of the nodes on the surface of the structures. The relationship between them reads as

where is a transformation matrix from the velocity vector for all nodes in the structure to the velocity vector for the surface nodes [22].

Substituting (21) into (22) yields

Further more, substituting (23) into (17) yields

Equation (24) is the coupled FEM-BEM equation, from which the sound pressure of the nodes on the interfaces can be obtained.

5. Solving Algorithm of the Coupled FEM-BEM Equation

Because of the nonsymmetric and fully populated matrices in the coupled FEM-BEM equations, the ACA approach is used in this study to reduce the memory requirements and the computation cost. During iteration process, the generalized minimal residual (GMRES) method iterative solver is adopted. In the following, the ACA and the GMRES algorithms are introduced briefly. More details can be found in relevant literatures [17, 18, 23].

5.1. Adaptive Cross Approximation

The ACA algorithm produced by M. Bebendorf and S. Rjasanow is an effective technique for solving nonsymmetric and fully populated matrices. The main process is as follows.

(1)Index Octree. The low-rank approximation of the matrix of FEM-BEM is based on asymptotically smooth functions, which happens in a well-separated domain. Thus, a tree structure for all boundary points should be constructed at first. For a three-dimensional domain, an octree is necessary to be built to describe the geometry relationship of sections on surface. The information that a tree should record includes its father, its children, its geometry information, and the points contained. The least nodes are called leaves whose numbers of boundary points contained are less than a given threshold value used to judge whether a new division is necessary.(2)Partitioning of Matrix. Based on the tree, the matrix in the FEM-BEM equation can be divided into some blocks. At first, the relationship between nodes should be established. For some node NP in an octree, its neighbour nodes are defined as the nodes which have the same geometry size and one common point at least. The interaction nodes are those nodes owing the same geometry size, whose fathers are neighbour nodes, but not for them. The relationship is shown as given Figure 2 in a two-dimensional domain. Then for each node of the tree, the elements contained are taken as source points and the elements of its neighbour as field points, and the submatrices coming from the source points and field points are defined as . The submatrices whose entries come from elements of the node and its interaction nodes are defined as . All of the submatrices and form the whole efficient matrix of FEM-BEM equation.(3)Low-Rank Approximation. The low-rank approximation in ACA is defined as

The target of ACA is to achieve

where is the tolerance of ACA.

Figure 2: Relationship between nodes.

The ACA algorithm has been applied to the low-rank blocks achieving approximately for both storage and matrix-vector multiplication [24]. It must be noted that the fully pivoted approach is well known to be much slower than the partially pivoted approach. The main reason behind this is that the fully pivoted approach requires the knowledge of the full matrix, whereas the partially pivoted approach would only require generation of individual matrix entries.

5.2. Generalized Minimal Residual Method

The generalized minimal residual (GMRES) method is one of the most popular iterative solvers for nonsymmetric linear systems. It was proposed by Saad and Schultz [25] and further developed by other researchers [26–28]. It has the property of minimizing at every step the norm of the residual vector over a Krylov subspace. The algorithm is derived from the Arnoldi process for constructing an orthogonal basis of Krylov subspace. Because of existence of error, the gained vectors lose orthogonality gradually. Therefore, preconditioner is needed. In the computer codes for this study, the solution of the diagonal blocks is taken as the initial solution, and the maximum number of iteration is set to 500.

6. Numerical Examples

In this section, the proposed coupling approach is applied to two examples for the simulation of acoustic-structure coupling problems.

6.1. Benchmark Example

The spherical test structure is an elastic sphere shell with a radius of 5 m and thickness of 0.15 m, full of some kind of liquid. A point sound source is located in the center of the sphere emitting harmonic spherical wave. The outer surface of the sphere is free, and the inner surface of the sphere is coupled with the acoustic field. The detailed geometrical and material parameters for the sphere and the fluid are shown in Table 1. In this example, the damping is ignored.

Table 1: Characteristics of the fluid and the sphere.

Materials
Properties
Values

Fluid
Acoustic velocity
1524 m/s

Density
1000 kg/m³

Sphere
Young’s modulus
2.07 × 10¹¹ Pa

Poisson’s ratio
0.3

Density
7669 kg/m³

Thickness
0.15 m

Radius
5 m

This problem is simulated by using FEM and BEM, respectively. The BEM results are obtained using the proposed algorithm. The FEM results are obtained using NASTRAN software. The sound pressure amplitudes obtained by the two methods are plotted in Figure 3. From the two curves, we can see that the results using BEM agree well with the results obtained from FEM. However, from Figure 4 in which computational efficiencies of different solving methods are compared, it can be seen that the ACA approach requires the least computing time.

Figure 3: Sound pressure curves obtained by BEM and FEM.

Figure 4: Comparison of computing time.

6.2. Analysis of the Frequency Response of a Rectangular Box

Further more, the frequency response analysis of a square box is also presented. The rectangular box is a closure chamber made of thin plates. Its inner and outer surfaces of the box are coupled with the interior and the exterior acoustic fields, respectively. In this example, the exterior acoustic field is assumed to be infinite. A point acoustic source located outside the box initiates the coupled response of the structure and its surrounding acoustic fields. The mesh of the box is shown in Figure 5.

Figure 5: Mesh of the square box.

Parameters of the box and the fluid are listed in Table 2.

Table 2: Characteristics of the fluid and the structure.

Materials
Properties
Values

Fluid
Acoustic velocity
344 m/s

Density
1.297 kg/m³

Structure
Young’s modulus
1.176 × 10¹¹ Pa

Poisson’s ratio
0.31

Density
4500 kg/m³

Thickness
20 mm

Length
2.36 m

Width
2.1 m

Height
3.6 m

To model this problem, the origin of the Cartesian coordinate system is set to coincide with the center of the bottom surface of the box, while the coordinate axes are parallel with the edges of the box. The box’s bottom surface is restrained. The point acoustic source generating harmonic spherical wave () is located at the point (4.0 m, 0.0 m, 4.0 m). Damping is ignored.

Two cases are considered. In the first case, the interior acoustic field is neglected, whilst in the second case the interior and the exterior acoustic fields are both considered. In the two cases, the sound pressures at a selected point (1.18 m, 0.0 m, 2.057 m) are calculated. They are compared in Figure 6.

Figure 6: Sound pressure under different frequencies.

From Figure 6, it can be seen that the trends of the two pressure-frequency curves are similar, but the pressures are a bit different at some specific frequencies; for instance, the frequency is around 20. Figure 7 demonstrates the acoustic pressure on the outer surface of the box. According to the graphs, it is found that, when the frequency varies from 20 to 100, the pressure distributions on the outer surface of the box are almost the same for the two cases. Therefore, it can be concluded that, for this example, the neglection of the inner acoustic field in the modelling does not lead to ruinous result.

(a)

(b)

Figure 7: (a) Sound pressure distribution for the coupled system without considering the inner field. (b) Sound pressure distribution for the coupled system considering the inner field.

7. Conclusions

In this study, an algorithm based on the coupling of FEM and BEM is developed for the analysis of acoustic-structure response. The FEM is employed to model the structure part, while the BEM is used to discretize the acoustic domain. The two domains are coupled with each other through the consideration of the appropriate compatibility and equilibrium conditions on the interface of the two domains. To improve the computational efficiency, the ACA approach is incorporated into the proposed algorithm to deal with the nonsymmetric and fully populated matrices resulting from the coupling of the FEM and BEM. By the proposed algorithm, the fictitious frequency problem can also be avoided. The validity and the high accuracy of the present algorithm are demonstrated by a benchmark example. The modeling for the frequency response of a rectangular box is also presented. The numerical examples show that the coupling method has great potential to deal with large-scale complex acoustic-structure problems.

Acknowledgment

This work is jointly supported by the National Natural Science Foundation of China (nos. 11272182 and 11072128).

References

K. J. Bathe, Finite Element Procedures, Prentice Hall, Englewood Cliffs, NJ, USA, 1996.

O. C. Zienkiewicz and R. L. Taylor, The Finite Element Method, Butterworth-Heinensann, Oxford, UK, 5th edition, 2000.

C. A. Brebbia and J. Dominguez, Boundary Elements: An Introductory Course, McGraw-Hill, London, UK, 1989.

O. C. Zienkiewicz and R. L. Taylor, The Finite Element Method, vol. 1-2, McGraw-Hill, London, UK, 4th edition, 1991.

O. C. Zienkiewicz, D. W. Kelly, and P. Bettess, “The coupling of the finite element method and boundary solution procedures,” International Journal for Numerical Methods in Engineering, vol. 11, no. 2, pp. 355–375, 1977.

C. C. Spyrakos and D. E. Beskos, “Dynamic response of flexible strip-foundations by boundary and finite elements,” Soil Dynamics and Earthquake Engineering, vol. 5, no. 2, pp. 84–96, 1986.

O. von Estorff, “Coupling of bem and fem in the time domain: some remarks on its applicability and efficiency,” Computers and Structures, vol. 44, no. 1-2, pp. 325–337, 1992.

H. A. Schenck, “Improved integral formulation for acoustic radiation problems,” Journal of the Acoustical Society of America, vol. 44, pp. 41–58, 1968.

A. J. Burton and G. F. Miller, “The application of integral equation methods to the numerical solution of some exterior boundary-value problems,” Proceedings of the Royal Society A, vol. 323, pp. 201–210, 1971.

I. O. Panich, “On the question of solvability of the external boundary value problem for the wave equation and Maxwell's equation,” Russian Mathematical Surveys, vol. 20, no. 1, pp. 221–226, 1965.

M. Guiggiani, G. Krishnasamy, T. J. Rudolphi, and F. J. Rizzo, “A general algorithm for the numerical solution of hypersingular boundary integral equations,” Journal of Applied Mechanics, vol. 59, no. 3, pp. 604–614, 1992.

J. J. Rgo Silva, H. Power, and L. C. Wrobel, “A numerical implementation of a hypersingular boundary element method applied to 3D time-harmonic acoustic radiation problems,” in Boundary Elements XIV, Computational Mechanics Publications, pp. 271–287, Southampton and Elsevier, London, UK, 1992.

Y. Saad, Iterative Methods for Sparse Linear Systems, SIAM, Philadelphia, Pa, USA, 2nd edition, 2003.

V. Rokhlin, “Rapid solution of integral equations of classical potential theory,” Journal of Computational Physics, vol. 60, no. 2, pp. 187–207, 1985.

L. Greengard and V. Rokhlin, “A fast algorithm for particle simulations,” Journal of Computational Physics, vol. 73, no. 2, pp. 325–348, 1987.

A. P. Peirce and J. A. L. Napier, “A spectral multipole method for efficient solution of large-scale boundary element models in elastostatics,” International Journal for Numerical Methods in Engineering, vol. 38, no. 23, pp. 4009–4034, 1995.

M. Bebendorf, “Approximation of boundary element matrices,” Numerische Mathematik, vol. 86, no. 4, pp. 565–589, 2000.

M. Bebendorf and S. Rjasanow, “Adaptive low-rank approximation of collocation matrices,” Computing, vol. 70, no. 1, pp. 1–24, 2003.

A. Sommerfeld, Partial Differential Equations in Physics, Academic Press, New York, NY, USA, 1949.

A. F. Seybert, B. Soenarko, F. J. Rizzo, and D. J. Shippy, “An advanced computational method for radiation and scattering of acoustic waves in three dimensions,” Journal of the Acoustical Society of America, vol. 77, no. 2, pp. 362–368, 1985.

R. Kress, “Minimizing the condition number of boundary integral operators in acoustic and electromagnetic scattering,” The Quarterly Journal of Mechanics and Applied Mathematics, vol. 38, no. 2, pp. 323–341, 1985.

M. J. Allen and N. Vlahopoulos, “Integration of finite element and boundary element methods for calculating the radiated sound from a randomly excited structure,” Computers and Structures, vol. 77, no. 2, pp. 155–169, 2000.

O. von Estorff, S. Rjasanow, M. Stolper, and O. Zaleski, “Two efficient methods for a multifrequency solution of the Helmholtz equation,” Computing and Visualization in Science, vol. 8, no. 3-4, pp. 159–167, 2005.

I. Benedetti, M. H. Aliabadi, and G. Davì, “A fast 3D dual boundary element method based on hierarchical matrices,” International Journal of Solids and Structures, vol. 45, no. 7-8, pp. 2355–2376, 2008.

Y. Saad and M. H. Schultz, “GMRES: a generalized minimal residual algorithm for solving nonsymmetric linear systems,” SIAM Journal on Scientific and Statistical Computing, vol. 7, no. 3, pp. 856–869, 1986.

C. Y. Leung and S. P. Walker, “Iterative solution of large three-dimensional BEM elastostatic analyses using the GMRES technique,” International Journal for Numerical Methods in Engineering, vol. 40, no. 12, pp. 2227–2236, 1997.

M. Merkel, V. Bulgakov, R. Bialecki, and G. Kuhn, “Iterative solution of large-scale 3D-BEM industrial problems,” Engineering Analysis with Boundary Elements, vol. 22, no. 3, pp. 183–197, 1998.

S. Amini and N. D. Maines, “Preconditioned Krylov subspace methods for boundary element solution of the Helmholtz equation,” International Journal for Numerical Methods in Engineering, vol. 41, no. 5, pp. 875–898, 1998.