System and Method for Predicting Sound Radiation and Scattering From an Arbitrarily Shaped Object

A. Background

The Kirchhoff integral formulation is one of the most widely used methods for predicting acoustic radiation and scattering from an elastic structure in engineering practice. The advantage of using this integral formulation is a reduction of the dimensionality of the problem by one. The first step of this approach is to determine the acoustic quantities on the surface. For an acoustic radiation problem, the normal component of the particle velocity (or the surface acoustic pressure for an acoustic scattering problem) must be specified. Since the surface is impermeable, the normal component of the particle velocity is equal to that of the surface velocity, which can be measured by an accelerometer. Next, the surface acoustic pressure (or the normal component of the particle velocity for acoustic scattering) is determined by solving an integral equation. Once these quantities are known, the radiated acoustic pressure anywhere can be calculated by the Kirchhoff integral formulation.

It may be difficult to use a conventional accelerometer to measure the vibration response of a flexible or a lightweight structure such as a loudspeaker membrane or a passenger vehicle fuel pump, because the weight of the accelerometer may alter the desired signal. In other cases, it may be unfeasible to use an accelerometer on a structure with a hostile environment such as an engine oil pan, where the temperature on the surface is extremely high. Under these circumstances, we must rely on a non-intrusive measurement technique.

One approach commonly adopted in practice is to use a laser vibrometer to measure the normal component of the particle velocity, which is equal to that of the surface velocity at the interface, and then solve an integral equation for the surface acoustic pressure. The shortcomings of such an approach are well-known: (1) the surface Kirchhoff integral equation may fail to yield a unique solution whenever the frequency is close to one of the eigenfrequencies associated with to the related interior boundary value problem, and (2) the numerical computation may become quite involved. This is because for an arbitrary surface, we must discretize the surface into many segments with several hundreds or even more nodes. Accordingly, we must solve a large number of simultaneous integral equations for the acoustic pressures at these nodes using boundary element method (BEM). Since the central processing unit (CPU) time increases quadratically with the number of the nodes, the computation process may be excessively time-consuming.

Actually, the laser technique can be used to measure the displacement and velocity vectors of a suspended microparticle in an insonified medium. The work in this area, however, has received much less attention than that of measurements of the out-of-plane motion of a vibrating structure. Summarized below are the basic principles and applications of the non-intrusive laser measurement techniques to measurements of the particle displacement and velocity vectors both in fluids and in air.

B. Laser techniques

1. Laser Doppler velocimeter (LDV)

LDV has become a standard tool for non-intrusive measurements of fluid particle velocities. The basic premise in the LDV measurements is that motion of the microparticles in the fluid (either due to natural impurities or due to seeded particles) will scatter the incident light, and produce a Doppler shift in the scattered light which can be detected with appropriate electronics and signal processing. LDV used in the fluid mechanics was extended to the acoustics by measuring the in-air particle velocities associated with steady-state time-harmonic standing waves and travelling waves inside a tube. Laser Doppler anemometry (LDA) has been used for the remote detection of sound. The technique of LDA consists of measurements of the velocity of neutrally buoyant microparticles suspended in an acoustic field by analyzing the spectral content of Doppler-shifted laser light scattered by the microparticles. LDV has been used to measure the acoustic particle velocity in fluids. In particular, the measurements of acoustic particle displacements using different LDV systems has determined that LDV was capable of detecting the particle displacements in the order of a few nanometers with a bandwidth of several kilohertz. The performance and limitations of LDV systems were also analyzed, and the effect of Brownian motion (i.e., thermal agitation in the fluid) on the measured data was shown to produce only negligible broadening of the spectral density of the signal of interest. An equation of motion of microparticles in suspension in an insonified fluid has been derived and it has been determined that the motion of neutrally buoyant microparticles closely emulates the displacement of the surrounding insonified fluid and confirms the basic tenet associated with the laser detection of sound.

2. Differential laser Doppler interferometry (DLDI)

DLDI is evolved from the principle of LDA and used to measure simultaneously the out-of-plane and the in-plane velocities on the surface of a vibrating object. The principle of DLDI is to measure the phase shift of the reflected or scattered light from the surface due to surface vibrational motion. The main component of a DLDI system is a probe head that has three illuminating single mode fibers. Prior to launching, the laser beams are frequency shifted by three acousto-optic Bragg cells by 40.0, 40.1, and 40.3 MHz, respectively, so the interference between the first and second beams occurs at 100 kHz, while those between the second and third and the first and third occur at 200 and 300 kHz, respectively. Geometrically, the first and second beams are positioned symmetrically with respect to the unit normal on the surface at an angle .alpha., and the third beam is aligned with the first and second beams at an angle .beta.(.beta.<.alpha.) with respect to the unit normal. In the differential configuration, the 100 kHz carrier will be modulated by the in-plane motion, and the 200 and 300 kHz carriers will be modulated by both in-plane and out-of-plane motions, respectively.

Mathematically, the surface displacement vector can be written as

where u.sub.in (t) and u.sub.out (t) represent the in-place and out-of-place components of the surface displacement, respectively, and e.sub.in and e.sub.out are the unit vectors in the corresponding directions.

Accordingly, the phase terms .phi..sub.ij, where i, j=1 to 3, in the 100, 200, and 300 kHz carriers can be written as

where k is the optical wavenumber. Hence by measuring the phase shifts .phi..sub.ij in the 100, 200, and 300 kHz carriers, one can determine the displacement vector on the surface. Since there are only two unknowns, one can use any two equations, say, Eqs. (2a) and (2b) to specify u.sub.in and u.sub.out.

In one optical system, the demodulation is done by using a combination of filters and the phase-locked loops (PLL). The PLL demodulates the signal with phase .phi..sub.ij and generates an output which is proportional to the time rate of changes of .phi..sub.ij. From Eqs. (2a) and (2b), one finds ##EQU1##

Therefore, by measuring the instantaneous frequency deviations d.phi..sub.12 /dt and d.phi..sub.13 /dt from the carrier frequencies at 100 and 300 kHz, one can determine simultaneously the in-plane and out-of-plane components of the surface velocity. A three dimensional laser vibrometer was designed based on this principle to measure simultaneously the three components of the velocity on the surface of a vibrating structure.

The DLDI technique can be extended in principle to the measurement of the velocity of a microparticle in the vicinity of a vibrating object. Imagine that an object is surrounded by neutrally buoyant microparticles. As the object vibrates, the acoustic pressure fluctuations will excite the microparticles into oscillations. Suppose that we define a control surface and focus the laser beams on a microparticle lying on that surface. There is no restriction on the formation of the control surface so long as it completely encloses the vibrating object. In the special case in which the control surface coincides with the vibrating surface, the normal component of the displacement of the microparticle will be equal to that of the surface displacement, while the tangible components may be different. In any event, the microparticle displacement in the directions normal and tangential to the control surface will cause a Doppler shift in the phase .phi..sub.ij of the reflected light, which is modulated in the frequency carriers. Once the signal with phase .phi..sub.ij is demodulated, we can calculate the microparticle velocity which is proportional to the time derivative of the phase, d.phi..sub.ij /dt.

3. Electronic speckle pattern interferometry (ESPI)

Alternatively, we can use ESPI to measure the phase term of a microparticle, which is an established optical technique for measuring static and dynamic deformations and surface shapes for more than two decades. Specifically, we can utilize the stroboscopic technique, which "freezes" the dynamic motion of a particle at one position so that during other times of the movement cycle, the particle is not illuminated and therefore is "invisible" to the imaging device. In practice, this technique can be implemented by using a pulsed laser or a light shuttering device with a continuous wave laser. The time interval between two consecutive pulses or shutters is typically in the range of nanoseconds, so ESPI can capture very high frequency oscillations. By using an additive-subtractive speckle pattern interferometry, the accuracy of the phase measurement can be further enhanced.

Suppose that we take five frames of additive speckle patterns of the motion of a microparticle in suspension S.sub.j, j=0, 1 . . . , and 4. Here the speckle pattern S.sub.o is taken with the laser illumination pulsed at the instant when the microparticle reaches its zero amplitude of a harmonic oscillation. The remaining four speckle patterns, S.sub.1 to S.sub.4, are taken with the laser illumination pulsed at the instants when the microparticle reaches its maximum and minimum amplitudes. During the acquisition of S.sub.0 to S.sub.4, the phase of the reference beam is shifted appropriately and is synchronized with the pulses. Accordingly, the displacement-induced phase term .phi. of a microparticle at any surface point x.sub.S can be written as ##EQU2## where F.sub.n, n=0, 1, 2, and 3, are the additive-subtractive fringe patterns which have the same form as that of the Michelson interferometric fringe pattern, except for the randomly distributed modulation term B/cos .psi./ contributed by the speckles,

Once the phase term .phi. is determined, the microparticle velocity which is proportional to the time derivative of the phase can be specified.

The present invention provides an explicit integral formulation for predicting acoustic radiation and scattering from an object of arbitrary shape. The radiated acoustic pressure is shown to be expressible as integrations of normal and tangential components of the particle velocity at the surface of a vibrating object, which are determined by a non-intrusive laser Doppler velocimeter.

The present invention is readily extensible to transient acoustic radiation. Since the acoustic pressure is expressed as an explicit function of the velocity distribution, one can use the Fourier transformation to obtain the acoustic pressure in the frequency domain first, and then take an inverse Fourier transformation to get the time-domain signal.

As illustrative examples only, the present invention can be utilized to predict acoustic radiation from automotive components, vehicle interior noise, the acoustic radiation from a vehicle as perceived by a stationary observer, sound radiation from submarines, etc.

The above, as well as other advantages of the present invention, will become readily apparent to those skilled in the art from the following detailed description of a preferred embodiment when considered in the light of the accompanying drawings in which:

FIG. 1 is a schematic of the system for predicting sound radiation from a vibrating object of the present invention;

FIG. 2 is a flowchart of the method of predicting sound radiation from a vibrating object of the present invention;

FIG. 3 shows the uniqueness of solution obtained by using Eq. (14) for dimensionless surface acoustic pressure for a finite cylinder;

FIG. 4a is a comparison of the dimensionless acoustic pressure given by Eq. (14), real part, and the exact solution;

FIG. 4b is a comparison of the dimensionless acoustic pressure given by Eq. (14), imaginary part, and the exact solution;

FIG. 5 is a comparison of the radiation pattern of the dimensionless field acoustic pressure give by Eq. (15), and the exact solution;

FIG. 6a is a comparison of the dimensionless acoustic pressure given by Eq. (14), real part, and the exact solution;

FIG. 6b is a comparison of the dimensionless acoustic pressure given by Eq. (14), imaginary part, and the exact solution;

FIG. 7 is a comparison of the radiation pattern of the dimensionless field acoustic pressure given by Eq. (15), and the exact solution;

FIG. 8a is a comparison of the dimensionless acoustic pressure given by Eq. (14), real part, and the exact solution;

FIG. 8b is a comparison of the dimensionless acoustic pressure given by Eq. (14), imaginary part, and the exact solution;

FIG. 9 is a comparison of the radiation pattern of the dimensionless field acoustic pressure given by Eq. (15), and the exact solution;

FIG. 10a is a comparison of the dimensionless acoustic pressure given by Eq. (14), real part, and the exact solution;

FIG. 10b is a comparison of the dimensionless acoustic pressure given by Eq. (14), imaginary part, and the exact solution; and

FIG. 11 is a comparison of the radiation pattern of the dimensionless field acoustic pressure given by Eq. (15), and the exact solution.

The present invention provides a system 20 for predicting acoustic radiation and scattering from an object 22 of arbitrary shape. A spherical object 22 is shown only for illustrative purposes. It will be demonstrated that the system 20 of the present invention is useful for objects 22 of arbitrary shape.

The system 20 preferably includes a CPU 26 having a monitor 28 and keyboard 30 or other input device. Further, as will be explained later, the system 20 may further include an accelerometer 32, which is preferably a non-intrusive laser Doppler velocimeter 32.

In operation, the shape, i.e. outer surface 36 of the object 22 is defined mathematically and input into the CPU 26 in step 40. This can be accomplished utilizing input device 30, CAD models, or imaging techniques, all of which are well known in the art. The particle velocity distribution over the surface 36 of the vibrating object 22 is then measured by a non-intrusive laser Doppler velocimeter 32 and input into the CPU 26 in step 42. Alternatively, the particle velocity distribution over the surface of the vibrating object 22 can be described mathematically and input into CPU 26 utilizing the input device 30.

In step 46, the CPU 26 integrates the particle velocity over the surface 36 of the vibrating object 22. In step 48, the CPU 26 determines the acoustic pressure at a point other than that the surface 36 of the vibrating object 22. As will be described below, this point could be interior or exterior to the surface 36 of the vibrating object 22.

The present invention provides an alternate formulation which enables one to predict the radiated acoustic pressure directly, once the particle velocity at the interface of a vibrating surface is determined by a non-intrusive laser Doppler velocimeter. The significance of this alternate formulation is two-fold: (1) solutions thus obtained are unique and (2) the efficiency of numerical computations is high. This is because in this alternate formulation, the radiated acoustic pressure is expressed as a surface integral of the particle velocity that can be implemented numerically using the standard Gaussian quadratures. There is no need to use BEM to solve a set of simultaneous integral equations for the surface acoustic pressures at the discretized nodes, which can be time-consuming for an arbitrarily shaped object. Such an approach can be desirable for analyzing the acoustic and vibration responses of the structure which requires the use of a non-intrusive laser measurement technique.

I. BASIC THEORY

Consider sound radiation from a finite object immersed in an unbounded fluid medium with density .rho..sub.o and sound speed c. Assume that the object vibrates at a constant frequency .omega., so an acoustic quantity can be written as a complex amplitude multiplied by a time dependence of e.sup.-i.omega.t. In this way, the complex amplitude of the radiated acoustic pressure p(x) at any point x can be written as ##EQU3## where G and .differential.G/.differential.n are the free-space Green's function and its normal derivative given, respectively, by ##EQU4## where n.sub.S stands for the outward unit normal on the surface S and R=.vertline.x-x.sub.S .vertline., here x and x.sub.S represent the receiver and source position vectors, respectively.

Equation (6) is an integral representation of the Helmholtz equation in the frequency domain. The normal derivative of the acoustic pressure on the right side of Eq. (6) is related to the normal component of the particle velocity .nu..sub.n (x.sub.S), which is equal to that of the surface velocity at the interface. For an acoustic radiation problem, the normal component of the surface velocity is specified. Hence, one must solve an integral equation obtained by letting the field point x in Eq. (6) approach the surface to determine the surface acoustic pressure p(x.sub.S) ##EQU5## where R is replaced by R.sub.S =.vertline.x'.sub.S -x.sub.S .vertline., here both x'.sub.S and x.sub.S are on the surface S. For an arbitrary surface, there is no closed-form solution to Eq. (8) and p(x.sub.S) must be solved numerically by using BEM. Once p(x.sub.S) and .nu..sub.n (x.sub.S) are all specified, the radiated acoustic pressure p(x) can be calculated by Eq. (6).

II. AN ALTERNATE FORMULATION

A formulation according to the present invention predicts the radiated acoustic pressure directly once the particle velocity at the interface is specified by a non-intrusive laser Doppler velocimeter. For completeness, we give formulations for both exterior and interior regions, respectively.

A. Exterior problems

Derivations of the alternate formulation start from the Euler equation,

Integrating both sides of Eq. (9) along a line which connects one field point x' to another x yields ##EQU6## where e is a unit vector in the direction of the line integral from x' to x.

Obviously, the integral of dp(x) on the left side of Eq. (10) is determined by its upper and lower limits, but independent of the choice of the integration path. Hence, we obtain ##EQU7##

Equation (11) shows that the complex amplitude of the acoustic pressure at any point x can be expressed as the sum of the complex amplitude of the acoustic pressure at another point x' plus a line integral of the apparent force per unit volume over any path that connects these two points. Since there is no restriction on the selection of the integral path, we choose to let it lie on the control surface S.sub.c. Accordingly, we can write ##EQU8## where .nu..sub.n and .nu..sub..eta. represent the normal and tangential components of the particle velocity at the surface S.sub.c, respectively, and dn and d.eta. are the increments in the normal and tangential directions on S.sub.c, respectively. In the special case where the control surface coincides with the real surface, the normal component of the particle velocity .nu..sub.n (x.sub.S) is equal to that of the surface velocity, but its tangential component .nu..sub..eta. may be different.

Without loss of generality, we will omit the subscript c on S.sub.c in the following. Since the integration path remains perpendicular to the unit normal n at all times, dn.ident.0. Hence the first term in the square brackets on the right side of Eq. (12) is identically zero, ##EQU9##

Substitute p(x.sub.S) into Eq. (8) and notice that p(x'.sub.S) represents the acoustic pressure at a fixed surface point x'.sub.S, which is independent of the integration with respect to the unprimed surface coordinates and therefore can be factored out of the integral sign. Combining the coefficients of p(x'.sub.S) on both sides of Eq. (8) yields a solution for p(x'.sub.S), which when substituted back into Eq. (13) gives a solution for p(x.sub.S). ##EQU10##

Substituting Eq. (14) into Eq. (6), we obtain the following integral formulation

where .sub.1,2 represent integral operators operating, respectively, on the normal and tangential components of the particle velocity ##EQU11##

Equation (15) is the main result of this invention. The radiated acoustic pressure p(x) is shown to be expressible in terms of the particle velocity only. This alternative formulation is in contrast with the classical Kirchhoff integral formulation. There is no need to solve an integral equation for the surface acoustic pressure, given the normal component of the surface velocity. Instead, the radiated acoustic pressure can be calculated directly once the particle velocity at the surface is specified. Since the surface integrals in Eq. (16) can be readily implemented by the standard Gaussian quadratures, the efficiently of numerical computations may be significantly enhanced.

B. Interior Problems

Following the same procedures as outlined above and changing the sign of the unit normal derivative .differential./.differential.n, we obtain

where p(X) represents the acoustic pressure at an interior point X enclosed by the surface S, and .sub.1,2.sup.int are defined by ##EQU12## C. Scattering problems

In a similar manner, we can extend this alternate formulation to the acoustic pressure field scattered from a finite object in free or half space. For an object in free space, the solution for the scattered acoustic pressure taken the form of

where .sub.1,2 are defined by Eq. (16), .nu..sub..eta..sup.total and .nu..sub.n.sup.total represent the components of the total particle velocity normal and tangential to the surface of the object, respectively,

where .nu..sup.sca and .nu..sup.int are the scattered and incident components of the particle velocity, respectively, and e.sub..eta. and e.sub.n represent the unit vectors in the tangential and normal directions at the surface of the object, respectively.

For an object in half space bounded by an infinite baffle of certain surface acoustic impedance, the effect of the acoustic pressure reflected from the baffle and that scattered from the object due to this reflected wave must all be taken into account. One way of solving the scattered acoustic pressure in the presence of an infinite baffle is to use the image source method. The resulting formulation can be written as

where .sub.3,4 are given by ##EQU13## where G is the Green's function that accounts for the effect of the image source. For a surface with arbitrary acoustic impedance, there is no closed-form solution for this Green's function. However, for an observer at a point x which is at least one half wavelength away from the surface, G can be approximated by ##EQU14## and its normal derivative is given by ##EQU15## where R and R.sub.I are the distances measured from the observer to the source and to the image, respectively, cos .theta..sub.I =n.sub.b .circle-solid.e.sub.I, here n.sub.b is the unit normal vector on the baffle and e.sub.I is the unit vector in the direction of wave propagation from the image to the observer, .beta. is the acoustic admittance of the baffle surface defined by ##EQU16## where Z (.omega.) is the surface acoustic impedance. For a rigid surface, Z (.omega.).fwdarw..infin., so .beta..fwdarw.0 and G=e.sup.ikR /R+e.sup.ikRI /R.sub.I. On the other hand, for a pressure-release surface, Z (.omega.).fwdarw.0, so .beta..fwdarw..infin. and G=e.sup.ikR /R-e.sup.ikRI /R.sub.I.

III. UNIQUENESS OF SOLUTION

It is well known that in carrying out the numerical computations for the radiated acoustic pressure in the exterior region, the surface Kirchhoff integral equation (8) may fail to yield a unique solution whenever the excitation frequency is close to one of the eigenfrequencies of the interior boundary value problem. This is because Eq. (8) shares the same eigenfrequencies as those of the corresponding integral equation in the interior region, and further, the solution to the adjoint homogeneous equation satisfies the compatibility condition for any velocity distribution .nu..sub.n (x). Since the alternate formulation (14) is derived from the Kirchhoff integral theory, an examination of its uniqueness seems to be in order. In what follows, we follow the procedures outlined by Scheneck to examine the uniqueness of Eq. (14). For brevity, the theorems regarding the characteristic values of an integration kernel, as those regarding trivial and non-trivial solutions to homogeneous and inhomogeneous equations as well as their adjoints are omitted.

The integral formulation governing the acoustic pressure in an interior region is given by ##EQU17## Taking a normal derivative of Eq. (25) at an interior point X and then lettering X.fwdarw.X'.sub.S from the inside, we obtain ##EQU18## Substituting p(X.sub.S) on the right side of Eq. (26) by p(X'.sub.S) plus a line integral as given by Eq. (13) yields ##EQU19## For the homogeneous Dirichlet problem p(X'.sub.S).ident.0. Hence Eq. (27) reduces to ##EQU20##

The integral equation governing the surface acoustic pressure in the exterior region is given by Eq. (14). If we consider the homogeneous Neumann problem in which .nu.(X.sub.S) and take the complex conjugate of this homogeneous equation, we obtain ##EQU21##

Obviously, Eq. (29) is different from Eq. (28), so they do not share the same eigenfrequencies. Consequently, Eq. (14) has a unique solution for the radiated acoustic pressure in the exterior region. However, this uniqueness may break down for the special case of a dilating sphere in which and .nu..ident..nu..sub.n e.sub.n and .nu..sub..eta. .ident.0. Since .nu..sub.n is a constant, it can be factored out of the integral sign and Eq. (28) reduces to ##EQU22## which is identical in form to Eq. (29) for the exterior problem. Therefore p and .nu..sub.n share the same eigenfrequencies .omega.*.

To check if p also satisfies the compatibility condition given by ##EQU23## we return to the integral formulation (25) for the interior region.

Taking the limit as X.fwdarw.X'.sub.S from the inside leads to ##EQU24##

Since .nu..sub..eta. (X.sub.S)=0, the second term on the right side of Eq. (32) vanishes identically. Therefore, for the interior homogeneous Dirichlet problem for which p(X.sub.S).ident.0 and .omega.=.omega.* we have ##EQU25## Because of the equivalence of .nu..sub.n (X.sub.S) and p(X.sub.S), we can rewrite Eq. (33) as ##EQU26##

Interchanging the order of integrations in Eq. (31) yields ##EQU27## Since the square-bracket term in Eq. (35) is identically zero ›see Eq. (34)!, the compatibility condition (31) is satisfied. Consequently, for a dilating sphere Eq. (14) may fail to yield a unique solution for the surface acoustic pressure whenever the excitation frequency is close to one of the corresponding interior Dirichlet eigenfrequencies. However, the ill-conditioning in the numerical computation of Eq. (14) is much less severe than that in the BEM-based Kirchhoff integral equation, as shown below.

For a dilating sphere, the solution for the surface acoustic pressure given by Eq. (14) takes the following form ##EQU28##

The corresponding interior Dirichlet eigenfrequencies are determined by the roots of the first kind of the spherical Bessel function of order zero, T.sub.o (k*a)=0, which yields k*a=m.pi., where m is an integer and a is the radius of the sphere. Using the spherical coordinates, it is easy to show that R.sub.S =2a cos (.theta./2), .differential.R.sub.S /.differential.n.sub.S =cos (.theta./2), and dS=a.sup.2 sin .theta.d.theta.d.phi., where .theta. and .phi. vary from 0 to .pi. and 0 to 2.pi., respectively.

Substituting R.sub.S, .differential.R.sub.S /.differential.n.sub.S, and dS into Eq. (36), we obtain ##EQU29## When ka=m.pi., Eq. (37) reduces to p=-i.nu..sub.0 c.rho..sub.n m.pi.(1-m.pi.).sup.-1 .times.(0/0). However, the computer cannot take the limit to (0/0). Hence any round-off error in the numerical computation may lead to an erroneous result.

In using the BEM-based Kirchhoff integral equation, the surface acoustic pressure is obtained by solving a set of simultaneous integral equations. When the excitation frequency is close to one of the interior eigenfrequencies, the diagonal terms become very small and the matrix becomes ill-conditioned. Table 1 exhibits this trend in detail for a dilating sphere around ka=.pi. and 2.pi.. Here numerical computations are carried out by a general BEM code with the spherical surface discretized into 48 quadratic quadrilaterals and 130 nodes. The results thus obtained are compared with those of Eq. (14).

With Eq. (14) one can carry out surface integrals using Gaussian quadratures directly, without the need of solving a set of simultaneous equations. In evaluating these integrations, we make use of the axisymmetry of the acoustic pressure distribution and discretize the spherical surface into twelve rings along the generator. The integrations with respect to the polar angle .theta. within each ring are carried out by the Gaussian quadrature formula with three interior points. The integration over the azimuthal angle .rho. can be done independently, yielding a factor of 2.pi..

Numerical results in Table 1 show that both BEM and Eq. (14) fail to yield unique solutions at the interior Dirichlet eigenfrequencies. However, the BEM results show the sign of ill-conditioning over a large frequency range around the interior Dirichlet eigenfrequencies, whereas the results obtained by Eq. (14) are more or less correct until ka hits the eigenfrequencies almost directly.

As a second example, we calculate the dimensionless acoustic pressures on the surface of an oscillating sphere. The corresponding interior eigenfrequencies are determined by the roots of the first kind of the spherical Bessel function of order one, T.sub.1 (k*a)=0, which yields k*a=4.493409, 7.725233, . . . . Since in this case the tangential component of the particle velocity is not zero and the normal component of the velocity is not constant, the complex conjugate of the homogeneous equation (29) for the exterior region does not share the same eigenfrequencies as those of the homogeneous equation (28) for the interior problem. Hence, the solution given by Eq. (14) is unique, while those obtained by the BEM-based Kirchhoff integral formulation are not as ka.fwdarw.k*a.

Table 2 lists the dimensionless acoustic pressures obtained by BEM and Eq. (14), respectively, on the surface of an oscillating sphere at .theta.=45.degree.. In using Eq. (14) we divide the spherical surface into eighteen rings along its generator in order to compare the surface acoustic pressure at the same location as those of BEM with 130 nodes. Numerical results demonstrate that the BEM results show the sign of ill-conditioning over a large frequency range around the interior eigenfrequencies, whereas the results of Eq. (14) are accurate and unique.

In the last example, we use Eq. (14) to calculate the acoustic pressure on the surface of a right cylinder of finite length. The aspect ratio of the cylinder is b/a=1, where a and b are the radius and half length of the cylindrical wall, respectively. The eigenfrequencies of the corresponding boundary value problem in the interior region are give by ##EQU30## where m is a positive integer and .alpha..sub.nq is the qth root of the nth Bessel function

Thus, for a=b=1 (m) the first two eigenfrequencies are k*.sub.101 a=2.8724 . . . and k*.sub.210 a=3.9563 . . . .

In the following, we examine the uniqueness of the numerical solutions given by Eq. (14) around ka=2.8724 and 3.9563, respectively. Since there exists no analytical solution for a finite cylinder, comparisons of the numerical result of Eq. (14) are made with respect to the radiated acoustic pressures from a point source of radius a.sub.o (a.sub.o =0.001a) located at the center of the cylinder. The procedures are described as follows. First, we calculate the acoustic pressure distribution on a cylindrical surface due to a point source. Next, we determine the particle velocity .nu.(x.sub.S) on the cylindrical surface via Eq. (9). Once .nu.(x.sub.S) is specified, the surface acoustic pressure is recalculated using Eq. (14). The results thus obtained are compared with the surface acoustic pressure due to the point source.

In carrying out the numerical integrations in Eq. (14), we uniformly divide the surface into 48 rings along its generator, i.e., 12 rings on each of the two flat ends and 24 rings on the side wall. Each of these 48 rings is further divided into 48 equal segments along the circumference. Numerical integrations over each segment are carried out by Gaussian quadratures with nine interior points. Table 3 lists the maximum relative errors of the magnitudes and phases of the dimensionless surface acoustic pressures given by Eq. (14) as compared with those from a point source around ka=2.8724 and 3.9563. It is seen that Eq. (14) shows no signs of the nonuniqueness difficulties around these characteristic frequencies. The relative error in the magnitude of the dimensionless surface acoustic pressure increases monotonically with ka, which is to be expected for a fixed grid size. The relative error in the phase of the surface acoustic pressure remains essentially the same however.

FIG. 3 displays the magnitude of the dimensionless surface acoustic pressure evaluated at the center of a flat end of the cylinder using Eq. (14) with ka varying from 0 to 10 at an increment at .DELTA.ka=0.1. Within this frequency range, there could be many characteristic frequencies for the corresponding interior boundary value problem ›see Eq. (39)!. Nevertheless, the numerical solution given by Eq. (14) is always unique.

In summary, solutions given by Eq. (14) are unique, except for the special case of a dilating sphere. Even under this circumstance, the ill-conditioning in the numerical computation is much less severe than that of the BEM-based Kirchhoff integral formulation.

IV. ANALYTIC VALIDATIONS

In this section, we demonstrate the validations of the alternate formulation derived in this invention on sound radiation from vibrating spheres whose analytic solutions are well-known.

A. A dilating sphere

In the first example, we consider acoustic radiation from a sphere of radius a vibrating radially at a constant frequency .omega. in a free field. Since the sphere dilates uniformly in all direction, the particle velocity .nu..ident..nu..sub.n e, where .nu..sub.n is a constant. Substituting .nu..sub.n into (15) yields ##EQU31##

The surface integrals in the square brackets on the right side of Eq. (40) can be evaluated independently and the results are ##EQU32## Substituting Eq. (41) into (40) leads to ##EQU33## where ##EQU34## where r is the distance from the center of the sphere to the observation point in the filed, .theta. and .phi. are the polar and azimuthal angles of the field point, and .theta..sub.S and .phi..sub.S are the polar and azimuthal angles of a surface point.

Since the radiated acoustic pressure is spherically symmetric, we can set the field point to lie on any axis, say, .theta.=.pi. and .phi.=0. Thus, R and .differential.R/.differential.n can be simplified to ##EQU35##

Substituting Eq. (44) into (42), we obtain ##EQU36## which agrees exactly with the analytic solution. B. An oscillating sphere

The second example concerns sound radiation from a sphere oscillating back and forth along the z axis at a constant frequency, .omega. in an unbounded fluid medium. To check the validity of the alternate formulation, we substitute the analytic solution for the particle velocity ##EQU37## into the alternate formulation and then compare the resulting acoustic pressure with that of the analytic solution. In Eq. (46), .theta. is the angle between the unit outward normal and the z-axis direction, .nu..sub.c is the magnitude of the velocity at the center of the sphere.

Because of the presence of .nu..sub..theta., evaluation of the integrals in Eq. (15) become a bit lengthy. Without loss of generality, we demonstrate the evaluation of the surface acoustic pressure given by Eq. (14). To this end, we first carry out the line integral involved in Eq. (14). Using the axisymmetry, we can set x'.sub.S at (a, .pi.,0) and x.sub.S at (a, .theta., .phi.). Substituting Eq. (46) into the line integral on the right side of Eq. (14) then yields ##EQU38##

Substituting Eqs. (41) and (47) into (14) leads to ##EQU39## which agrees perfectly with the analytic solution. C. A partially vibrating sphere

The third example deals with acoustic radiation from a sphere on which only part of the surface is vibrating at a constant frequency .omega. As before, we substitute the analytical solution for the particle velocity into Eqs. (14) and (15), and then compare the acoustic pressures thus obtained with the known values. For simplicity, we assume that the particle velocity is axisymmetric, and that the normal component of the particle velocity is equal to that of the surface velocity, which is limited to the portion described by the half-vertex angle .theta..sub.o ##EQU40## where H(.theta..sub.o -.theta.) is the Heaviside step function, which is unity when .theta..ltoreq..theta..sub.o and zero when .theta.>.theta..sub.o, .nu..sub.o is the magnitude of the surface velocity, P.sub.m and h.sub.m represent the Legendre and the spherical Hankel functions of order m, respectively, and B.sub.m and .delta..sub.m stand for the amplitudes and phase angles of acoustic radiation of order m, respectively.

Substituting Eq. (49) into (14) and (15), we obtain ##EQU41##

Once again, we use the axisymmetry of the acoustic pressure distribution to set .phi.=0 and replace R and .differential.R/.differential.n in Eq. (50), respectively, by ##EQU42##

In this case, however, an exact solution cannot be found. Hence p(x.sub.S) and p(x) are solved numerically. Here, we discretize the spherical surface into twelve rings along .theta. and twenty-four segments along .phi. within each ring. Numerical integrations over each segment are carried out using Gaussian quadratures with none interior points. For convenience, we set ka=0.1, kr=1, .theta..sub.0 =45.degree., and use sixteen expansion terms in Eq. (49b) to approximate the tangential component of the particle velocity .nu..sub..theta.. Table 4 compares the numerical results of Eq. (50a) at various polar angles .theta. with the analytic solutions. The maximum relative error in the magnitude of the surface acoustic pressure is found to be about six percent, while that in the corresponding phase angles is less than one percent. The magnitudes of the relative errors decay rapidly as the measurement point moves into the field. Table 5 shows the comparison of the numerical result of Eq. (50b) at the same polar angles at kr=1 with the analytic solution. In this case, the maximum relative error in the magnitude of the Field acoustic pressure is less than one percent.

V. NUMERICAL VALIDATIONS

In the preceding section, we demonstrate validations of the alternate formulation on sound radiation from vibrating objects whose analytic solutions are well-known. In this section, we consider cases for which there exist no analytic solutions.

In order for the vibrating objects to be of certain generality, we choose a slender cylinder with two flat ends. The aspect ratio of the cylinder is b/a=10, where a and b are the radius and half length of the cylindrical wall, respectively. Since there are no analytic solutions available, we compare the acoustic pressures given by Eqs. (14) and (15) with those of a monopole and a dipole sources, respectively, located at the center of the cylinder. The radii of the monopole and the dipole are a.sub.0 =0.001a, where a=1(m). The procedures involved in these comparisons are the same as those described in Section III.

In what follows, we first consider the case in which the particle velocity field on the cylindrical surface is generated by a monopole located at the center of the cylinder. Substituting the velocity field into Eqs. (14) and (15) yields the surface and field acoustic pressures, which are then compared with those from the monopole source, respectively. Excellent agreements are obtained for all the cases tested. For brevity, however, we only plot the comparisons of kb=1 and 50 below. The surface integrals in Eqs. (14) and (15) are implemented with the cylindrical surface discretized into 480 rings on the side wall and 24 rings on each end. These rings are further discretized into 48 segments along the azimuth. Numerical integrations over each segment are carried out using Gaussian quadratures with nine interior points.

FIG. 4 demonstrates the comparisons of the real and imaginary parts of the dimensionless acoustic pressure obtained by using Eq. (14) with those of a monopole at kb=1 along the generator of the cylindrical surface. Here the abcisa represents a dimensionless distance s/a along the generator of the cylindrical surface. In particular, s/a=0 indicates the center of the side wall, s/a=.+-.10 the edges that separate the side wall from the flat end, and s/a=.+-.11 the centers of the two ends. FIG. 5 shows the comparison of the radiation pattern of the dimensionless field acoustic pressure given by Eq. (15) at a radial distance of r=20(m) with that of the monopole source.

In a similar manner, we calculate the surface and field acoustic pressure by using Eqs. (14) and (15) at kb=50, and compare the results thus obtained with those of the monopole source (see FIGS. 6 and 7).

Next, we repeat the same procedures and compare the calculated surface and field acoustic pressures with those of a dipole at kb=10 and 100. Because of the increase in the excitation frequency, we double the number of discretization segments. Specifically, we use 960 rings on the side wall and 48 rings on each end, and 96 segments along the azimuth within each of these rings. The comparisons of the real and imaginary parts of the surface acoustic pressure for kb=10 are depicted in FIG. 8, and those of the filed acoustic pressure are displayed in FIG. 9. FIGS. 10 and 11 demonstrate the comparisons of the surface and field acoustic pressures for kb=100.

In accordance with the provisions of the patent statutes, the present invention has been described in what is considered to represent a preferred embodiment. However, it should be noted that the invention can be practiced otherwise than as specifically illustrated and described without departing from its spirit or scope.