Alexandre Leblanc
email: alexandre.leblanc@univ-artois.fr

Antoine Lavie

A meshless method for the Helmholtz eigenvalue problem based on the Taylor series of the 3-D Green's function

Keywords: numbers: 43, 20, Ks, 43, 55, Ka

The solution of the Helmholtz eigenvalue problem is considered through the use of the method of the fundamental solutions. Taylor series of these solutions are employed to form a polynomial eigenvalue problem. The presented method differs from other methods such as the multiple reciprocity method. Here, the Green's function itself is expanded and no integration is performed. Results on classical geometries (sphere, parallelepiped box and finite cylinder) demonstrate the accuracy of the method for the determination of the eigenvalues with Neumann, Dirichlet and Robin boundary conditions. Furthermore, the center of the Taylor approximation is shown to be adjustable, allowing the method to be theoretically effective for any arbitrarily part of the eigenvalue spectra.

Introduction

The characterization of acoustic cavities still represents a crucial aspect in many applications, from the design phase to the challenges of renovation and acoustic comfort improvement. In the automobile industry, acoustic eigenvalue analysis is routinely employed to have car interior resonances away from the vibration frequencies (engine, gearbox, shaft etc.) [START_REF] Brebbia | Boundary element, an introduction course[END_REF]. The determination of the acoustic modes is usually performed with numerical methods, since no exact solution exists when the studied structure has a complex geometry. The mesh generation of arbitrarily shaped cavities could be time consuming, especially with the finite element method [START_REF] Zienkiewicz | The finite element method[END_REF], which has been employed since decades for determining the acoustic eigenvalues of irregularly shaped cavities [START_REF] Petyt | A finite element method for determining the acoustic modes of irregular shaped cavities[END_REF][START_REF] Marburg | Normal modes in external acoustics. Part I: Investigation of the one-dimensional duct problem[END_REF][START_REF] Marburg | Normal modes in external acoustics. Part II: Eigenvalues and eigenvectors in 2D[END_REF]. The boundary element method (BEM) has been also widely developed in this purpose [START_REF] Nowak | The multiple reciprocity method. A new approach for transforming BEM domain integrals to the boundary[END_REF][START_REF] Chen | A 3D boundary element method for determination of acoustic eigenfrequencies considering boundary conditions[END_REF][START_REF] Steinbach | A boundary element method for the Dirichlet eigenvalue problem of the Laplace operator[END_REF][START_REF] Peters | Modal decomposition of exterior acoustic-structure interaction[END_REF] as this formulation reduces the problem's spatial dimension by one. Recently, new numerical methods have recently spawned a renewed interest in meshless (or called mesh-free) methods such as the indirect Trefftz method [START_REF] Jin | Trefftz method for Kirchhoff plate bending problems[END_REF], the smoothed particle hydrodynamics method [START_REF] Monaghan | Simulating free surface flows with SPH[END_REF], the radial basis function collocation method [START_REF] Amaziane | Radial basis function and genetic algorithms for parameter identification to some groundwater flow problems[END_REF], the method of fundamental solutions (MFS) [START_REF] Kupradze | A method for the approximate solution of limiting problems in mathematical physics[END_REF]. These methods have been successfully applied to the eigenproblems of acoustical cavities. Karageorghis et al. [START_REF] Karageorghis | The method of fundamental solutions for the numerical solution of the biharmonic e q u a t i o n[END_REF] employed the MFS to calculate the eigenvalues of 2-D cavities and Chen et al. [START_REF] Chen | Boundary collocation method for acoustic eigenanalysis of threedimensional cavities using radial basis function[END_REF] applied the boundary collocation method with radial basis function for acoustic eigenanalysis of 3-D cavities. For these two examples the eigenvalue parameter still embedded in the resultant matrix system leads to a nonlinear eigenproblem. As a consequence, the extraction of the eigensolutions is not straightforward and can lead to numerical spurious eigenvalues [START_REF] Leblanc | A Wave Superposition method based on monopole sources with unique solution for all wave numbers[END_REF] or iterative solution schemes (with necessary search interval, solution accuracy depending on the step size).

In this paper, an eigenproblem formulation dedicated to 3-D acoustic cavities, based on the MFS and requiring the polynomial approximation of the free-field Green's function, is proposed. This method is truly meshless as elements are not required for both interpolation and integration. Only a discrete set of nodes are used on the boundary of the domain while a set of virtual sources surrounds the cavity. By employing the Taylor series of the Green's function and applying the Dirichlet, Neumann or Robin boundary conditions, polynomial eigenproblem are obtained. The detailed description of this acoustic eigenvalue problem is shown in Section I, recalling the classical MFS and introducing its formulation with the polynomial approximation of its kernel (eg. Green's function). Three numerical examples are shown in Section II: sphere with Dirichlet or Neumann conditions, rectangular parallelepiped with Dirichlet conditions and finite cylinder with mixed boundary conditions. The last example reuses the sphere with Robin conditions.

I. PROBLEM FORMULATION

A. Equations of the Problem

The harmonic wave propagation problem is governed by the Helmholtz equation [START_REF] Morse | Methods of theoretical physics[END_REF] ∆p(r)+k

2 p(r)=0 r ∈ Ω , ( 1
)
where ∆ is the Laplacian operator, p is the pressure amplitude at the point r, k = ω/c is the wave number with ω the circular frequency and c the speed of sound through the fluid medium. The variables exhibit an implicit e -iωt time dependence. This equation is satisfied in an interior domain Ω and subject to a homogeneous condition at its boundary Γ of the form

a(r)p(r)+b(r) ∂p(r) ∂n =0 r ∈ Γ , (2)
where a(r)a n db(r) are known complex-valued functions of r and n is the outward normal to the boundary at r. The non-trivial solutions k * and p * (r) are expressed as the eigenfrequencies and eigenfunctions. In this paper, four different eigenproblems are considered:

• Dirichlet -a(r)=1andb(r)=0, • Neumann -a(r)=0andb(r)=1, • Robin -a(r)/b(r)=ikρc/Z,
• mixed boundary conditions -{a(r 1 )=0,a(r 2 )=1} and {b(r 1 )=1,b(r 2 )=0}, where r 1 (resp. r 2 )i so nt h es u r f a c eΓ 1 (resp. Γ 2 ) on which Dirichlet (resp. Neumann) conditions are located, with Γ=Γ 1 ∪ Γ 2 . The Robin boundary conditions is set in the context of an acoustic cavity covered by an absorbent, locally reacting at r, with Z being the surface impedance.

B. The Method of Fundamental Solutions

The fundamental solution of the Helmholtz equation ( 1) is defined by [START_REF] Karageorghis | The method of fundamental solutions for the numerical solution of the biharmonic e q u a t i o n[END_REF] ∆+k 2 G k (r, r s )=-δ(r, r s ) ,

where r s are the coordinates of source points. The fundamental solution is then obtained as

G k (r, r s )= e ikx 4πx , (4)
where x = |rr s | is the distance between the collocation point and the source. Using the MFS with a single-layer potential approach, the field solution is assumed to be

p(r i )= j ϕ j G k (r i , r s j ) , (5)
∂p(r i ) ∂n = j ϕ j ∂G k (r i , r s j ) ∂n , (6)
where ϕ j is the intensity of the source point located at r s j . For Dirichlet boundary condition, Eq. ( 5) at r i is rewritten as

j ϕ j G k (r i , r s j )=0, (7)
where r i is located on the surface. This equation can be recast as a N × N linear system U k {ϕ} = {0} where N is the number of source points which is chosen equal to the number of collocation points. The Neumann eigenproblem can be formulated with Eq. ( 6) in the same way. Also, for the mixed or Robin boundary conditions, Eq. ( 5) and Eq. ( 6) (depending on the boundary condition applied at r i ) are used to form a linear system of the same type.

In the context of finding the eigenfrequencies of acoustic cavities, the different linear systems are also nonlinear eigenproblems for k being the eigenvalues as these equations have non-trivial solutions k * p . The eigenfunctions of Helmholtz equation (cf. Eq. ( 1)) could be easily evaluated with Eqs. ( 5) and [START_REF] Nowak | The multiple reciprocity method. A new approach for transforming BEM domain integrals to the boundary[END_REF].

MFS formulations have been successfully applied to the Helmholtz problem in domains with and without interior holes [START_REF] Tsai | The method of fundamental solutions for eigenproblems in domains with and without interior holes[END_REF], connected domains [START_REF] Chen | Eigensolutions of multiply connected membranes using the method of fundamental solutions[END_REF] and for mixed boundary conditions [START_REF] Antunes | A meshfree numerical method for acoustic wave propagation problems in planar domains with corners and cracks[END_REF]. The direct determinant search method and singular value decomposition techniques have been extensively employed to solve the resulting eigenproblem. Recently, Reutskiy [START_REF] Yu | Reutskiy: The method of fundamental solutions for helmholtz eigenvalue problems in simply and multiply connected domains[END_REF] proposed a new algorithm for which the eigenvalue problem is reduced to a sequence of inhomogeneous problems. However, the solution procedures still involves an iteration process, requiring the system matrix to be formed repeatedly and leading to poor results in the case of closely spaced frequencies. Using BEM, many methods have been dedicated to the eigenanalysis of acoustic cavities [START_REF] Ali | Advances in acoustic eigenvalue analysis using boundary element method[END_REF] and some of them allow the direct determination of natural frequencies as the internal cell method [START_REF] Bezine | A mixed boundary integral-finite element approach to plate vibration problems[END_REF], the dual reciprocity method [START_REF] Nardini | A new approach to free vibration analysis using boundary elements[END_REF], the particular integral method [START_REF] Ahmad | Free vibration analysis using a BEM particular integrals[END_REF][START_REF] Leblanc | An acoustic resonance study of complex three-dimensional cavities by a particular integral method[END_REF][START_REF] Leblanc | Numerical analysis of eigenproblem for cavities by a particular integral method with a low frequency approximation of surface admittance[END_REF], the multiple reciprocity method [START_REF] Nowak | The multiple reciprocity method. A new approach for transforming BEM domain integrals to the boundary[END_REF] and its equivalent series expansion method [START_REF] Kirkup | Solution of the Helmholtz eigenvalue problem via the boundary element method[END_REF]. This last is of particular interest here as the series expansion technique is the core of the proposed method. The main difference is that the expansion is performed formally on x and not on the product kx. Furthermore, using MFS kernels, the matrix expansion requires no integral evaluation.

C. Series Expansion of the Fundamental Solutions

First, we recall Taylor's development of the exponential function centered at γ as

f (x)=e ix = e iγ ∞ m=0 i m m! (x -γ) m . (8)
From Eq. ( 8) and using substitution in series for k, the Green function is rewritten as

G k (x)= e iγ 4π ∞ m=0 k m i m m! (x -γ/k) m x . (9)
It is worth noting that using this substitution with γ = 0 divides the original radius of convergence by |k|.F r o m Eq. ( 9), the MFS kernel for Neumann boundary condition is easily computed via the relation

∂G k (x) ∂n = G k (x) x ik - 1 x (r -r s ) • n . (10)
The main idea of our approach for solving the Helmholtz eigenvalue problem is to obtain a polynomial eigenproblem formulation with the approximation of the MFS kernels (Eq. ( 9) and ( 10)). The easiest way to achieve an efficient formulation is to truncate the previous series at order α and set γ = 0 so that the three different boundary conditions linear systems can be immediately recast as

α m=0 k m A m {ϕ} = {0} , (11)
where A m are matrices containing terms relating to Eq. ( 9) or Eq. ( 10), evaluated for each field point r and source location r s .U s i n gγ = 0 is not the only possible choice to obtain a polynomial eigenvalue problem (PEP). Indeed, truncating Eq. ( 9) at order α and using the Newton's binomial series gives

G k (x)= e iγ 4π α m=0 k m x m-1 α-m n=0 i m-n m!n! γ n . ( 12
)
Using polynomial approximations to solve nonlinear eigenvalue problems is not new. In particular, for boundary integral formulations of partial derivative equation eigenvalue problems, the use of Taylor approximations is a wellknown approach [START_REF] Nowak | The multiple reciprocity method. A new approach for transforming BEM domain integrals to the boundary[END_REF][START_REF] Kirkup | Solution of the Helmholtz eigenvalue problem via the boundary element method[END_REF][START_REF] Kamiya | A new complex-valued formulation and eigenvalue analysis of the Helmholtz equation by boundary element method[END_REF]. However, and to the authors' best knowledge, such methodology has not yet been employed with the MFS formulation and for this particular context. Furthermore, the aforementioned Taylor series can be centered ad lib.. This allows the application of this method over the frequency band of interest, while maintaining a relatively low order polynomial.

D. Linearization of the Polynomial Eigenproblem

A usual way of solving PEP is to use a matrix augmentation procedure to transform Eq. ( 11) into an equivalent generalized eigenproblem (GEP). Standard techniques could be applied to solve this linear eigenvalue problem. The linearizations most commonly employed are the companion forms: for

P(k)= m i=0 k i A i , one of these companion forms is [30] L(k)=k ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ A m 0 ... 0 0 I n . . . 0 . . . . . . . . . . . . 0 ... 0 I n ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ + ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ A m-1 A m-2 ... A 0 -I n 0 . . . 0 . . . . . . . . . . . . 0 ... -I n 0 ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ . ( 13
)
An important issue for these linearizations is its conditioning. It is now known that different linearizations for a given PEP can have very different conditioning [START_REF] Tisseur | Backward error and condition of polynomial eigenvalue problems[END_REF], so that numerical methods may produce rather different results for each companion form. Two key indicators must be considered with care to get low forward errors: the matrices condition number and the backward error. The backward error of k * for P and its companion linearization L are only of similar value if

A m ≈ ••• ≈ A 0 ≈1 . (14)
Same constrain applies if considering matrix conditioning. So, in order to get a robust solution [START_REF] Betcke | Optimal scaling of generalized and polynomial eigenvalue problems[END_REF], a basic norm balancing procedure is employed, with

P (µ)=β 1 P (β 2 k) . (15)
The choice of β 2 =( A 0 / A m ) 1/m can be shown [START_REF] Betcke | Optimal scaling of generalized and polynomial eigenvalue problems[END_REF] to minimize the scaling of the eigenproblem defined by Eq. [START_REF] Chen | Boundary collocation method for acoustic eigenanalysis of threedimensional cavities using radial basis function[END_REF]. Furthermore, and to achieve Eq. ( 14), β 1 is chosen as (max

A i ) -1 or i A i / i A i 2 .
Other improvements on the condition number could be performed by diagonal scaling. In this work, no iterative refinement is done to improve backward error so the results shown in the following section are the output of the polyeig function implemented in MATLAB which uses a variant of linearization defined by Eq [START_REF] Kupradze | A method for the approximate solution of limiting problems in mathematical physics[END_REF], on the PEP reformulated by Eq. ( 15). This linearization can be directly derived from the companion form family proposed by Antoniou and Vologiannidis [START_REF] Antoniou | Vologiannidis: A new family of companion forms of polynomial matrices[END_REF]. Following Peters et al. [START_REF] Peters | Modal decomposition of exterior acoustic-structure interaction[END_REF], construction of symmetric matrices is possible if considering energy quantities. In this case, structure preserving linearizations [START_REF] Higham | Symmetric linearizations for matrix polynomials[END_REF][START_REF] Mackey | Good vibrations from good linearizations[END_REF] should be rather used.

II. NUMERICAL INVESTIGATIONS

In this section, some results issued from the above method, called the method of expanded fundamental solutions (MEFS), are presented. The first tested geometry is a sphere. The second is a parallelepiped box. Both bodies have separable geometries for which analytic solutions to the Helmholtz eigenvalue problem are available [START_REF] Morse | Methods of theoretical physics[END_REF].

The full set of solutions of the GEP associated to the companion form (cf. Eq. ( 13)) contains the approximations to the true eigenvalues of the Helmholtz problem. Because of Eq. ( 11), these approximations tend to have small imaginary parts. From the results shown in this section and for sake of clarity, the imaginary part is not shown when below 0.2% of the real component. The other eigenvalues obtained have significant imaginary parts and hence they can be ignored from the true eigenvalues.

MEFS principles reveal that error is driven by four key parameters:

• surface nodes density,

• distances between nodes and MFS sources,

• order α of the Taylor series,

• choice of the center of expansion γ.

Beyond the error induced by the method itself, solving with confidence the PEP is predominant. As stated in previous section, PEP issued from this approach can be treated in different ways. To this end, one may also consider a contour integral method [START_REF] Asakura | A numerical method for polynomial eigenvalue problems using contour integral[END_REF] in order to avoid the linearization. Here, no particular structure of the PEP matrices is assumed. The results shown in this section prove that a matrix scaling procedure is sufficient to obtain accurate PEP solutions for regular geometries, even if a non-optimal linearization is used.

For the numerical investigations performed, the MFS sources are distributed on a sphere of radius r MFS =1.5 × r 0 centered at the geometrical origin of the tested structures (where r 0 is a typical length), with N well-separated sources obtained using a centroidal Voronoi tessellations procedure [START_REF] Du | Tetrahedral mesh generation and optimization based on centroidal voronoi tessellations[END_REF]. Obviously, setting MFS sources location independently of the geometry to be treated is not an optimal choice, as the distances between sources and surface nodes define, with α, the range of validity of the computed eigenvalues. With this arbitrarily placement of the virtual sources, the focus is willingly set to the parameter α for the following results. Further works are planned to give this method its full theoretical background.

For the Neumann boundary condition, normals at nodes are needed. In the following tests, these nodes are computed through a ball pivoting method [START_REF] Bernardini | The ball-pivoting algorithm for surface reconstruction[END_REF] involving two key steps. The first step is the determination of a pseudo-mesh consisting of triangular elements with their outward normals. Then, for each node, normals at connected triangles are weighted (with node distance) and averaged to give the normal estimation at node.

A. Sphere

This test consists on a sphere of radius r 0 . Nodes determination at its boundary are performed through the triangulation of an icosahedron in order to create a geodesic sphere of different densities. Two densities are considered: the first consists of 42 nodes corresponding to an inter-nodes distance d =0 .5646 × r 0 and the second has 92 nodes with d =0 .3625 × r 0 . The Dirichlet (resp. Neumann) eigenvalues for this geometry are the zeros of j l (k)( r e s p . j ′ l (k)), the spherical Bessel function of order l. Results for both boundary conditions and γ =0a r es h o w ni n Tables I and II (zero Hz modes are not shown, being always smaller than 1e -13). For the Dirichlet study, the two results for the α = 60 cases of Table I are displayed in Fig. 1 (only the positive quadrant is displayed). Effects of nodes density is clearly visible for = 42. Higher eigenfrequencies are badly recovered despite the fact that they lie in the radius of convergence of the series. Indication on method accuracy in term of mesh criterion (i.e. node density) is also shown in Fig. 1, where 3 nodes per wavelength ensure a sufficient discretization of the sphere surface: from r 0 ×|k 1 | (approx. 2.7 nodes per wavelength), multiple modes are roughly recovered but above, modes become indistinguishable. Nevertheless, it must be recalled that the increase of surface points density is only useful if α and r MFS are set accordingly. As expected, the order of series truncation α is of prior importance for the range of accuracy of the proposed method. In Table I, when the series is centered at 0, α = 30 gives only the two first eigenvalues while doubling this order gives accurate results for the six first expected values (at least). Also, reliable estimations of out of range eigenvalues can be obtained with the arbitrarily centering of the Green's function expansion. We define q satisfying γ = qπ for the results shown in Table III and Fig. 2 illustrates the use of this centering for the Dirichlet case and N = 42. When compared with higher order results (see N =4 2a n dα = 30 of Table I), third and fourth eigenvalues can be recovered, with time and memory saving costs. Nevertheless, and as the radius convergence of series is translated, lower eigenvalues are attended to be deteriorated to benefit others. Furthermore, and as illustrated by the q = 3 case of Fig. 2, results near the radius of convergence are also linked to the nodes density. This can be harmful when dealing with multiple eigenfrequencies, for example when attempting to determine the fifth eigenvalue with this particular configuration. Indeed, a cluster of points with a significant dispersion around the expected value is observed. However, optimal eigenvalue (lowest imaginary part) in this area is found equivalent to the fifth eigenvalue obtained with N =4 2a n d α = 45 (cf. Table I). The first simulation case involves determination of the eigenmodes of a rectangular parallelepiped with Dirichlet conditions. The dimensions are chosen as (0.9 × 1 × 1.1) × r 0 in order to produce a dense distribution of the first eigenvalues. 172 nodes are used to discretize the parallelepiped. Table IV illustrates the MEFS reliability in this context and its flexibility as eleven modes are well evaluated for q = 2 and with only a 40th order polynomial for the Green's function approximation (thus ending to a GEP of order 6880 through linearization).

α =30 α =45 α =60 α =30 α =45 α =60 1

C. Cylinder with mixed boundary conditions

The problem is now a cylinder of radius 0.4 × r 0 , and height 2 × r 0 with mixed boundary conditions (Neumann condition everywhere except at one end cap where Dirichlet condition is assumed). The surface is discretized with 140 nodes. Various configurations of the series parameters (α and q) are used and the MEFS performance is shown in Table V. Previous considerations obtained from the two first study cases are also observable, and the results testify of the MEFS efficiency to deal with mixed boundary conditions. This last application concerns the sphere cavity described in the previous section. Here, we consider Z = ρc on the entire sphere surface so the boundary conditions is rewritten as

. α =20 α =30 α =40 q =0 q =3 q =0 q =2 1
ikp(r)+ ∂p(r) ∂n =0 r ∈ Γ . (16)
Using this equation with Eqs. ( 9) and ( 10), complex eigenfrequencies must occur. Thereby, and in order to identify the original problem's eigenvalues, using a threshold on the imaginary part is no longer a reliable technique. One possibility is to consider that the expected eigenfrequencies are independent of the MFS radius. Thus, moving this radius should reveal these invariants. Figure 3 illustrates the spotting of the eigenfrequencies of the Helmholtz problem with Robin boundary condition, using two different locations of MFS sources. Although not being the most optimal method in this case, satisfactory results are obtained (references values given by the FEM solver COMSOL). A meshless method for the Helmholtz eigenvalue problem is presented and tested for different cavities. This method is based on the polynomial approximation of the fundamental solutions for the 3D Helmholtz equation. A norm balancing procedure is employed on matrices of the resulting polynomial eigenproblem, then this last is linearized. The results show the reliability of the proposed method with convergence to the exact solution as the number of nodes and the accuracy of the polynomial approximation increase. Furthermore, the present paper demonstrates the possibility for this approximation to be centered on a selected range of sought eigenvalues. Mixed boundary conditions are easily handled by the MEFS and its extension to Robin boundary conditions is possible, although the interpretation of results can be difficult. A sorting scheme to identify physical eigenvalues has been proposed. Nevertheless, removing spurious modes has to be one field of further research. Also, no particular care was provided in this paper on virtual source locations. Indeed, the determination of the number of virtual sources and their location deserve sophisticated numerical procedures [START_REF] Shigeta | Adaptive multilayer method of fundamental solutions using a weighted greedy QR decomposition for the Laplace equation[END_REF]. Here, the virtual sources are set to obtain matrices conditioning good enough to effectively solve the PEP. Additionally, the results show that a basic norm balancing procedure of the matrices is sufficient to ensure the robustness of the MEFS, at least for the investigated problems. Other techniques as diagonal scaling of coefficient matrices or backward error improvement by iterative refinement can be considered in order to enhance the accuracy of the computed eigenvalues. Other eigensolvers, which prevent the inflation of the matrix dimension, must be also investigated. Focusing only on eigenvalues of physical interest, the contour integral method [START_REF] Asakura | A numerical method for polynomial eigenvalue problems using contour integral[END_REF][START_REF] Leblanc | Solving acoustic nonlinear eigenvalue problems with a contour integral method[END_REF] is a potential candidate for the range of applications handled by the MEFS.

FIG. 1 :

 1 FIG.1: Distribution of the Dirichlet eigenvalues for the sphere with α = 60 and for 2 nodes densities

FIG. 2 :

 2 FIG. 2: Distribution of the Dirichlet eigenvalues for the sphere (N =42nodesandα = 20) and for two values of q

5 COMSOLFIG. 3 :

 53 FIG. 3: Distribution of the Robin eigenvalues for the sphere (N =42nodesandα = 40) and for two values of rMFS

TABLE I :

 I Sphere test problem -Dirichlet eigenfrequencies (r0 × k * )

	Mode			MEFS			Theory
	num.		42 nodes			92 nodes
		α =30 α =45 α =60 α =30 α =45 α =60
	1	3.1416 3.1416 3.1416 3.1416 3.1416 3.1416 3.1416
	2	4.4998 4.4946 4.4946 4.5001 4.4934 4.4934 4.4934
	3	-	5.7583 5.7600	-	5.7634 5.7634 5.7634
	4	-	6.2825 6.2828	-	6.2834 6.2831 6.2831
	5	-	6.9630 6.9821	-	6.9650 6.9880 6.9879
	6	-	-	7.7776	-	-	7.7248 7.7252

TABLE II :

 II Sphere test problem -Neumann eigenfrequencies (r0 × k * )

	Mode	MEFS	Theory
	num.	42 nodes	92 nodes

TABLE III :

 III Sphere test problem -Dirichlet eigenfrequencies (r0 × k * )w i t hN = 4 2a n dα =20

	Mode			MEFS		Theory
	num. q =0 q =1 q =2	q =3
	1 3.1455 3.1450 3.1450 3.1603 -i0.07 3.1416
	2	-	4.4937 4.4934	4.4937	4.4934
	3	-	-	5.7428	5.7601	5.7634
	4	-	-	6.2798	6.2826	6.2831

TABLE IV :

 IV Rectangular parallelepiped test problem -Dirichlet eigenfrequencies (r0 × k * )w i t hN = 172 and q =2

	Mode		MEFS	Theory
	num. α =22 α =40
	1	5.4964 5.4964 5.4965
	2	7.3937 7.3947 7.3947
	3	7.7308 7.7345 7.7343
	4	8.1247 8.1712 8.1710
	5	-	9.1811 9.1810
	6	-	9.5512 9.5517
	7	-	9.7706 9.7706
	8	-	9.8170 9.8170
	9	-	10.4452 10.4483
	10	-	10.9962 10.9929
	11	-	11.1890 11.1836
	TABLE V: Cylinder test problem -Computed mixed eigenfrequencies (r0 × k * )w i t hN = 140
	Mode		MEFS	Theory
	num