A Comprehensive Review of Boundary Integral Formulations of Acoustic Scattering Problems

This is a review presenting an overview of the developments in boundary integral formulations of the acoustic scattering problems. Generally, the problem is formulated in one of two ways viz. Green’s representation formula, and the Layer-theoretic formulation utilizing either a simple-layer or a double-layer potential. The review presents and expounds the major contributions in this area over the last four decades. The need for a robust and improved formulation of the exterior scattering problem (Neumann or Dirichlet) arose due to the fact that the classical formulation failed to yield a unique solution at (acoustic) wave-numbers which correspond to eigenvalues (eigenfrequencies) of the corresponding interior scattering problem. Moreover, this correlation becomes more pronounced as the wave-numbers become larger i.e. as the (acoustic) frequency increases. The robust integral formulations which are discussed here yield Fredholms integral equations of the second kind which are more amenable to computation than the first kind. However, the integral equation involves a hypersingular kernel which creates ill-conditioning in the final matrix representation. This is circumvented by a regularisation technique. An extensive useful list of references is also presented here for researchers in this area.


Introduction
T he main purpose of this review is to give the general mathematical public an overview of the developments that occurred in the field of Boundary integral methods for acoustic scattering problems.Here, we are mainly concerned with the acoustic wave.The study of wave propagation in solids and fluids is very important in many branches of physics and engineering.We concentrate mainly on the area of time-harmonic acoustics.Acoustic problems like assessing and reducing the noise produced by aircrafts, vehicles, household appliances and machines in general have been of major concern and interest to industry.
Scattering problems (or radiation problems) in mathematical physics are concerned with finding the solution of the Helmholtz equation,

∇
(1) 0 2 2 = + φ φ k in an exterior domain.In the process of finding this solution, the integral equation method is very convenient.The Helmholtz equation (or the reduced wave equation) arises when we seek steady-state solutions of the scalar wave equation for monochromatic (single frequency) sources.To give a general idea, consider a typical scattering problem: an incident wave with harmonic time-dependence is intercepted by a bounded scatterer.Consequently, reflected and diffracted waves are produced which travel outwards from the scattering region, the propagation of the waves being described by the wave-equation, where c is the speed of propagation and Φ is as defined in (3) below.It may also be regarded as the excess pressure at in the exterior domain and at time t.Let us suppose that the wave-length of the source radiation is λ which corresponds to a frequency c/λ Hz, and assume that steady state has been reached so that all waves present have harmonic time-dependence of this frequency.We can then write, p where φ is a complex function in the spatial co-ordinates and ω = kc.The constant k is the wave-number, defined as k = 2π/λ is a dimensionless quantity.By substituting (3) in equation ( 2) we see that φ is a solution of the Helmholtz equation given in (1).
The problem as described above is generally speaking a diffraction problem.Such problems occur widely in acoustics, electromagnetic theory, hydrodynamics etc.It is very important to know how the incident wave field is modified by the scatterer.Another type of scattering problem arises when the sources are actually on the scattering surface.In such situations one needs to know how the surface modifies the radiation characteristics of the source.As far as the mathematics is concerned, these two classes of problems are almost similar.For diffraction problems it is convenient to regard the total wave-function φ t at any point as being made up of two parts : a known incident wave φ inc , and a scattered wave φ sc which is to be determined (the time factor e -iωt is omitted), so we have . The incident wave is simply the wave-function that would exist in the absence of the scattering surfaces, and the scattered part is a wave diverging from the scattering region.A radiation condition is necessary for the scattered wave in order to ensure that it represents a diverging wave.At the boundary of the scatterer, conditions are imposed on the total wave φ t , or its normal derivative or possibly a linear combination of the two, and the mathematical problem is to find an exterior wave-function which satisfies these boundary conditions and is such that the scattered wave obeys the radiation conditions.In reality the shape of the scatterer is usually arbitrarily shaped .This demands a numerical method of solution especially if the solution is needed in the midfrequency range.And as such, it is convenient to reformulate the problem as an integral equation.
Integral equations have been used extensively in formulating diffraction problems for analytical solution.The integral equations cannot usually be solved exactly in closed form, although the formulation very often provides a very useful basis for asymptotic analysis.During this century much effort has been concentrated on developing new alternatives, such as numerical techniques, to deal with the mathematical modelling of the problems.Time-harmonic acoustic scattering problems are among the engineering problems for which the boundary element methods (BEM) have already been shown to yield powerful numerical techniques.As we will see later, the fundamental feature of the integral equation method is one of shifting the problem of solving a domain equation to a boundary equation which is usually an integral equation of the second kind.By applying BEM numerical techniques the solution is achieved, at least as an approximation.The advantages of the BEM technique are that: only the boundary of the region (i.e. the surface scatterer) needs discretisation; the explicit calculations of results at domain points i.e. in the region; the possibility of using continuous and discontinuous elements in the mesh; and the capability of modelling infinite domains by discretisation of the internal boundary (surface of the scatterer only).During the latter half of the seventies, the appearance of texts like Jaswon and Symm (1977) and Brebbia (1978), on integral equations and their numerical solutions by discretisation, proliferated the popularity of integral equation methods among the engineering community.Later on, Brebbia et al (1984) broadened the range of applications of BEM.The basic ideas of the BEM are not new to acoustic and structural engineers.Since 1960s several integral equation methods have been used.In acoustics, these include, Shaw and Friedman (1962), Banaugh and Goldsmith (1963), Chen and Schweikert (1963), Copley (1967).The BEM was first applied to general transient elastodynamics in 1968 by Cruse and Rizzo (1968).Since then some more S. I. ZAMAN dynamic applications emerged such as Niwa et al (1976), Domínguez (1978), Manolis and Bescos (1981), and Cole et al (1978).
BEM is best suited when dealing with unbounded regions.In domain-based numerical methods, such as finite element method (FEM) there are well known disadvantages.The use of integral equation methods presents many advantages, notably in that the meshing region is reduced by one dimension i.e. from the infinite domain exterior to the body (scatterer) to its finite surface, and the Sommerfeld radiation condition at infinity is automatically satisfied.However, there is a well-known drawback in that the standard integral equation formulations obtained from the Green's representation formula, for time-harmonic acoustics, fail to have a unique solution at frequencies corresponding to the eigenfrequencies of the adjoint associated interior problem.This is not a physical drawback and nor indeed does it appear when the problem is represented by partial differential equations, but comes up entirely as a result of a deficiency of the integral equation representation.We shall look into this in detail later.Over the last thirty years several techniques have been proposed to circumvent these fictitious eigenfrequencies.Indeed the existence and uniqueness of solutions of the radiation problem, i.e. the exterior Neumann boundary value problem for the Helmholtz wave equation, has been proved for all values of wave-numbers by Leis (1958).However, this matter is still one of the main areas of research in the BEM.For acoustic problems, some of the most promising techniques have been proposed by Panich (1965), Schenck (1968), Kussmaul (1969), Burton and Miller (1971), Ursell (1973) and Jones (1984,85).More recently, Zaman (1994) proposed a technique using a modified hybrid layer potentials, based on the theory of potentials.
Let us first consider the equations governing acoustics.The propagation of small amplitude acoustic waves through a fluid medium is as described by the linear wave equation given in (2), where Φ is a velocity potential and c is the velocity of sound in the medium.When the motion is considered to be timeharmonic, the function Φ may be represented as in (3), where φ is a reduced velocity potential and ω is the angular frequency, and thus equation ( 2) reduces to the Helmholtz equation ( 1), where k is the wave-number as defined before.Hence the excess pressure, after substituting (5) in (4), is And consequently by using (3), the velocity potential φ is related to the acoustic pressure P by ρ o is assumed constant here.Helmholtz equation ( 1) is in fact the governing equation of a steady-state sound wave propagation problem.The acoustic scattering problem seeks the solution of this equation in the exterior region where the solution satisfies some boundary conditions.

Boundary Conditions
There are basically two types of scattering surfaces which yield two different types of boundary conditions.One condition is for perfectly soft scatterers and the other for perfectly hard scatterers.At the surface of a soft scatterer the excess pressure p is zero which in view of ( 7 p The condition ( 8) is a form of homogeneous Dirichlet boundary condition and the condition ( 9) is a form of homogeneous Neumann condition.In order to determine the exterior solution uniquely, it is not sufficient to have the boundary condition on the surface only.But from the physical point of view it is essential that the scattered wave be entirely outgoing.In order to ensure this, we have to impose a radiation condition: In two dimensions, where r is the length of the vector measured from a fixed origin in the scattering region to a general field point.This condition was first expressed formally by Sommerfeld (1949), and later by Kupradze (1965), Atkinson (1949), Rellich (1943) andWilcox (1956).

p
Let us now collate all the above and consider the mathematical statement of the general scattering problem.We have seen how time-harmonic wave problems in acoustics and indeed in diffraction, can be reduced to finding solutions of Helmholtz equation, and we have also considered some of the more appropriate boundary conditions, (8) and (9).In scattering problems, wave-functions are sought in exterior domains, which satisfy not only boundary conditions at the surface of the scatterer but also the radiation condition at infinity.The surface boundary conditions are usually of three main types: on the boundary (iii) with φ sc satisfying the radiation condition in (10).
In a radiation-only situation, the problem is similar.Here the solution behaves like a scattered wave satisfying the radiation condition, hence the problem to can be specified in the form (11), with the incident wave terms replaced by given functions describing the surface properties of the sources and the scatterer.Analytical solutions to typical scattering problems are seldom obtainable.Surface shapes and specified normal velocities are quite general.Even asymptotic approximation methods (large or small k) are usually inappropriate because the solutions of interest are often those for which the wavelength of the source radiation is comparable with the size of the scatterer.As such, numerical approach is the only logical option.Although for an elliptic equation (i.e.Helmholtz equation) it is natural to think in terms of either finite difference method or finite element method (FEM), there are some basic difficulties, one of which is that the region of interest is at infinity.On the other hand, the boundary element method (BEM) is appropriate for this problem.

Background Theory
Before considering different integral equation formulations for the interior as well as the exterior Helmholtz problems (acoustic scattering problems for instance), it will be worthwhile looking into the definitions of simple-layer and double-layer Helmholtz potentials, and also investigating their properties.

Helmholtz Surface Potentials
Consider a domain − B bounded by a closed surface ∂B, and let σ be a density function defined on ∂B.The simple-layer and double-layer potentials are defined as follows: Simple-layer : where L k and M k are the simple and double-layer potential operators, x and y are positions in three dimension.The integration is performed at y∈∂B, n y is the outward normal at y and x is any other point either in the interior , exterior B − B + or on the boundary ∂B.g k (x, y) is the free-space Green's function for Helmholtz equation, i.e., ( ) where δ is the Dirac delta function and g k satisfies the radiation condition in (10).When the time factor e -iωt is included, g k represents the effect experienced at x of a point y (or vice-versa) radiating into free space in either two or three dimensions, more precisely, it represents the potential at x generated by a unit point source at y (and vice-versa).In three dimension, this function is The free space Green's function is closely related to the corresponding functions for Laplace's equation (k = 0), in that they possess the same singularity for | x -y |=0.This may readily be seen by expanding the functions for small | x -y |< ε, where ε is very small.So we find that in three dimensions, as It is this weak singularity which is responsible for the continuity properties of the surface potentials, particularly the jump properties as the point x crosses the surface ∂B.Because the singularity is the same for Helmholtz potentials as for the more familiar potentials associated with Laplace's equation, the properties can be obtained as in classical potential theory (Kellogg, 1929;Günter, 1967).

Surface
It is clear from the definitions (12) and (13) of the potentials that whether or not the integrals exist and indeed whether or not they define the operatorss L k and M k , depends on three interrelated factors: The development of classical potential theory followed two paths.Russian workers (Günter, 1967) built up the theory based on the concept of the Lyapunov surface, whereas Kellogg (1929) in his classic work introduced the idea of a regular surface.For these two types of surface slightly different conditions are required on σ in order to establish the existence of the potentials.In the present review we will not delve into details of these two types of surfaces.However let us look into the definitions of these surfaces.

S. I. ZAMAN
A Lyapunov surface S satisfies the following three conditions (Smirnov, 1964;Günter, 1967): (i) a tangent plane exists at every point of the surface S.
(ii) ∀x ∈ S, ∃ε > 0 such that S ∩{y | |x -y| < ε} intersects lines parallel to the normal at x in at most one point.It is self-evident that if this property holds for any given ε it also holds for any smaller ε.

(iii) if
, where are normals at where A, α are constants, and 0 < α ≤ 1.As a result of condition (ii), the portion of the surface within the sphere may be represented in terms of a tangent-normal system of coordinates centred around the point .

x
Kellogg's regular surface : the definition is in several parts (Müller, 1969), (i) a curve is regular if it is continuous, has finite number of continuously differentiable arcs, has no double points and has finite length.
(ii) a regular surface element is a point-set such that, for at least one cartesian coordinate system (x, y), it can be represented in the form z = F(x, y) , (x, y) ∈ H where F is continuously differentiable in H, and H is a finite closed area of the x-y plane which is bounded by a regular curve.The x-y plane would usually be the tangent plane at some representative point in the surface element.
(iii) a regular surface is a point set which can be divided into a finite number of regular surface elements.
(iv) a surface is closed if each of the edges of every regular surface element belongs to two surface elements.
(v) A point is regular relative to a regular surface if there exists a decomposition of the surface into regular surface elements so that the point does not belong to an edge.It is regular relative to a curve if the curve is continuously differentiable at that point.
(vi) A curve or surface is called 'smooth' if all its interior points are regular points.
(vii) A regular region is a compact point set bounded by a finite number of closed regular surfaces such that no two of the surfaces have a point in common.
If ∂B= U and each i is a Lyapunov surface, then provided x is not on an edge or vertex of ∂B, it can be shown that L k and M k exist for all other x∈ℜ 3 provided σ is Hölder continuous over ∂B.The existence and continuity of the normal derivative of M k requires that the source density σ ∈ C (2) , (Burton, 1973;Colton and Kress, 1983).

Properties of Helmholtz Potentials and Their Derivatives
Let the potentials be defined on closed surfaces bounding the regular regions.The surfaces are assumed to be regular in the sense of Kellogg.And also assumed to be 'smooth' in the sense described above.+ .However, they define different solutions in the interior and exterior regions.At the boundary the potentials exhibit various irregularities due to the singularity of g k .The existence everywhere of L k and M k is assured if σ is merely continuous (Günter, 1967).Kellogg (1929) showed that the assumption σ ∈ C ( 2) is sufficient to ensure the existence (and continuity) of the derivative, and this is the condition we are imposing.The i B ∂ immediate consequence of ( 16) is that the jump and continuity properties of the Helmholtz potential at the boundary parallel to those of the Laplace potential at the boundary.For the smoothness condition on the boundary ∂B= , we merely require that each ∂B n i=1

The
U i be a Lyapunov surface.We also require that σ be Hölder continuous over the boundary i.e. σ satisfies the following inequality for any y 1 ≠ y 2 on the boundary.
Before we can review the continuity properties of the potentials at the boundary, we must understand the difference between the limiting value and the direct value of the integrals defining the potentials i.e.

Limiting values :
, are on the normal line at x∈∂B Subject to the conditions described above, the simple-layer potential L k has the following properties: A. It exists and, is continuous and differentiable everywhere in B + ∪ and satisfies the Helmholtz equation i.e.
σ and it also satisfies the radiation condition (10).
B. It exists and is continuous on the boundary despite the singularity ( 16), as this is essentially a weak singularity.Also, its value at x in ∂B is continuous with its neighbouring values in and in i.e.
C. There exists two distinct normal derivatives at x in ∂B, one on each side of ∂B.
We use the convention that n -both increase moving away from ∂B.At any point x ± in B + ∪ − B on the normal line through x in ∂B we have, where g k ′ (x ± , y) denotes the derivative of g k in the direction of the normal passing through x ± keeping y fixed (Jaswon and Symm, 1977).But at the initial point x, For the double-layer potential M k , the treatment is similar, subject to the same conditions.Here are the properties: A. It exists, and is continuous and differentiable everywhere in B + ∪ − B and satisfies the Helmholtz equation ( ) and in B + it satisfies the radiation condition (10).
B. It exists on ∂B even though there is a singularity as |x -y|→ 0 , which is essentially a weak singularity and consequently for x ∈ ∂B However, its value at x in ∂B is not continuous with its neighbouring values in B + ∪B , i.e.
Here are some more useful notations, which is the normal derivative at x of the simple-layer potential ; and which is the normal derivative at x of the double-layer potential.Note that is the transpose of , hence the symbol.The derivative with respect to n consequently the kernel of N k becomes hyper-singular i.e. non-integrable.We shall look into some integral equation formulation of scattering problems which involve such hyper-singular kernels.This poses numerical complications which require special treatments.Rêgo Silva et al (1992,1993) used Hadamard's finite part integration in their treatment of the hyper-singular kernel.

Helmholtz Formulae
These are integral representations of the scattering problem derived from a variant of Green's second theorem.They also utilize the Helmholtz potentials defined above.When the field point of the representation is allowed to approach the surface, we get integral equations connecting the boundary values of the wave-function and its normal derivative.The representation utilizes Green's second theorem: Let φ 1, φ 2 be any scalar functions with continuous second derivatives in a domain S ∪ Ω, where S is the boundary of the domain.Then we have,

[
]dV If both φ 1, φ 2 satisfy Helmholtz the equation in Ω, then (17) takes the form, 0 We get the Helmholtz formulae by replacing φ 1 by any solution of the Helmholtz equation φ(x), regular in Ω, and φ 2 by g k (x, y).So, applying this to an exterior wave-function φ + in B + which satisfies the radiation condition (10), and to the free space Green's function g k , we get, And, as X → w ∈ , due to further 'jump' of the double-layer potential we have, The boundary equation ( 19b) is known in the literature as the Surface Helmholtz Equation (SHE), which in operator form becomes For the Dirichlet boundary condition we have which is Fredholm integral equation of the first kind with φ + given on the boundary.
For the Neumann condition we have which is Fredholm integral equation of the second kind with ∂φ + ⁄∂n y given on the boundary.A similar application of Green's second theorem to an interior wave-function φ -in − B gives the following:

Differentiated Helmholtz Formula
We can construct integral representations of the derivatives of φ in either the interior or exterior domains by differentiating the Helmholtz representations for φ presented above.For the wave function in the interior, For an exterior wave-function,

Layer-Theoretic Formula
Since both L k σ and M k σ are radiating wave functions in B + , it would be convenient to express the exterior solution of the Helmholtz equation by means of layer potentials.

Exterior Dirichlet Problem (EDP):
We seek a solution in the form of a simple-layer potential (22a) and by continuity across the boundary the corresponding boundary equation becomes This is a Fredholm integral equation of the first kind for the surface density σ in terms of the prescribed boundary value of φ .Alternatively, we seek a solution in the form of a double-layer potential (23a) and by using the jump properties the boundary equation becomes This is a Fredholm integral equation of the second kind for µ in terms of the prescribed boundary value of φ.

Exterior Neumann Problem (ENP):
We seek a solution in the form (22a) and then differentiate in the direction of the normal at x, taking the limit towards x ∈ ∂B, to get the boundary equation This is a Fredholm integral equation of the second kind for σ in terms of the prescribed boundary value of ∂φ ⁄∂n.Alternatively, seeking a solution in the form (23a), differentiate in the direction of the normal at x, taking limit towards x ∈ ∂B, gives the boundary equation This is a (non-Fredholm) integral equation of the first kind for µ in terms of the prescribed boundary value of ∂φ ⁄∂n.As noted before, the N k operator has a hypersingular kernel.Hadamard's finite-part integration (Hadamard (1923); Rêgo Silva et al 1992), treatment is required here.A formulation along this line has been proposed by Jin (1993).
Interior Dirichlet Problem (IDP): Simple-layer form gives The boundary equation takes the form Alternatively, the double-layer formulation gives The boundary equation takes the form Interior Neumann Problem (INP) : Applying the Neumann condition to the simple-layer formulation above yields, Also, by applying the Neumann condition to the double-layer formulation above yields Note that (24a) which refers to ENP is the transpose of (25b) which refers to IDP.Also note that (23b) which refers to EDP is the transpose of (26) which refers to INP.The discussion on this connection follows.

Connection between the interior and exterior problems
Integral equation formulations for Helmholtz problems (acoustic) in the exterior region repeatedly make reference to the related interior problem.As suggested in the introduction that there may be no physical connection between the two regions (e.g. the surface may be perfectly rigid), nevertheless the integral equation formulations for the two regions often have similar integral operators.Burton (1973) gives a very clear analysis of this.When the wave-number k takes certain discrete values, the interior problems with homogeneous boundary conditions have non-trivial solutions.These values of k are the eigenvalues of the interior problem and the corresponding non-trivial solutions are the eigenfunctions.It may be shown (Burton, (1973) that these eigenvalues must be real.At these values the corresponding integral equation governing the exterior case breaks down i.e. it may either fail to yield a unique solution or it may be insoluble for the given inhomogeneous term.Zaman (1994)  writing the resulting equation in operator form, we find that ∂φ ⁄∂n satisfies Now inserting the same boundary condition into the (differentiated) boundary equation ( 21b) and writing the resulting equation in operator form, we get, ( ) where I is the identity operator.It is clear that the boundary values of ∂φ ⁄∂n satisfy both (28) and ( 29 Therefore it follows from (20c) that φ ≠ 0 on ∂B, suggesting that (30), (31) have non-trivial solutions.This result has an important implication to the exterior Dirichlet formulation.
Let us now consider how the above phenomena affect different formulations of the exterior problems.We look at first the layer-theoretic formulations.
Connection with exterior problems: Suppose we seek a solution of the exterior Neumann problem using a simple-layer formulation Taking the normal derivative and inserting the Neumann condition with ∂φ ⁄∂n x prescribed on ∂B, we get, This boundary equation has a unique solution if the homogeneous equation has only the trivial solution (Zaman, 1994).However, since (37) is mathematically identical to (29), it has non-trivial solutions for wave-numbers k∈K D i.e.ENF fails to yield unique solution for these wavenumbers.Arguing similarly it can be shown that EDF fails to give unique solution if k∈K N .Also, by using the Helmholtz formula we arrive at the same conclusion.However a distinction can be drawn between the Layer-theoretic formulation and the Helmholtz formula.The inhomogeneous equation ( 36) is insoluble because the free term ∂φ ⁄∂n x is not orthogonal to the eigenfunctions of the homogeneous adjoint equation.By contrast, the Helmholtz formula which yields at least non-unique solutions since the free term satisfies the compatibility condition [cf.Fredholm Alternative].Zaman (1994) gives an interesting illustration of this connection between the interior and the exterior problems with reference to a spherical surface.Jaswon (1963) is one of the early significant works on integral equation methods in Potential theory.It presents the ground work for most of the subsequent works in this area.
This peculiar defect and hence a drawback in the classical formulation led several workers to propose alternative formulations in order to circumvent this.This phenomenon has been known for at least fifty years (Lamb, 1945).In the next chapter we will look at those proposals which made significant contribution to the development of the integral equation formulations of Helmholtz problem.

Improved Integral Equation Formulations
The formulations discussed at the end of the previous chapter are known as the classical integral equation formulations.They were defective in that those formulations fail to yield a unique solution at certain critical values of the wave-number k, which are the eigenvalues of the corresponding homogeneous interior problem.We discuss below some of the significant propositions which improve this situation.Brundrit (1965) in his paper presented a proposal, utilizing a layer potential i.e. indirect formulation, and developed a numerical method, which did not suffer from the effect of critical wave-numbers.Considering the problem of hard acoustic scattering by a fixed obstacle of boundary ∂B, he proposed that the scattered wave be represented as ( 38)

Brundrit
This is a classical formulation.Applying the Neumann boundary condition, it is noted that the source density σ satisfies the boundary equation However, the corresponding homogeneous equation as a result of using direct formulation, which yield Brundrit further noted that by continuity of a simple-layer potential and by uniqueness, the wave-function generated by the eigensolutions of ( 40) must vanish identically in the exterior i.e.
Therefore he argued that component of a general solution of (39) would make no contribution in the exterior.
* j σ Alternatively, we can utilise the non-trivial eigensolutions to generate an exterior potential, where (+) signifies that the normal points are into the exterior.Since, ∂V/∂n + = 0 for each j , it follows, by uniqueness that V vanishes identically in the exterior i.e.
and by continuity for each j.Therefore as noted by Brundrit, the component of the general solution does not make any contribution in the exterior domain.
* j σ Nevertheless, Brundrit overlooks the fact that the boundary integral equation ( 39) has no solution whenever k coincides with any of the critical values, since the eigenfunctions of the corresponding adjoint equation of ( 40) are not orthogonal to the inhomogeneous term in (39).Note that Brundrit's approach can also be formulated by Helmholtz formula which yields a consistent non-unique solution.

Panich
By using layer potentials, Panich (1965) where µ is a complex constant, chosen to be +i if ℜ(k) ≥ 0, see Burton (1974Burton ( , 1976)).Accordingly, σ satisfies the following boundary equation, Since ( 42) contains a hypersingular operator N k , which is non-compact, Panich proposed replacing (41) by the hybrid layer-potential ( 43) where the subscript o denotes putting k = 0 in the definition of L k .Therefore σ satisfies Here, all the operators are weakly singular and hence compact.This is one of the earliest examples of the regularization technique for weakening the singularity of the kernel.Using Fredholm theory, Colton and Kress (1983) deduced that (45) has a unique solution for all non-negative real k provided ℑ(µ) > 0 .
However, numerical implementation is very expensive and cumbersome due to the application of the product of matrix approximations to the operators.

Copley
The proposal by Copley (1967) utilised the Helmholtz formula at points in the interior region i.e.
Thus avoiding the singularity problem.The idea was originally proposed by Kupradze (1965).In operator form, which is a functional relation between φ and ∂φ ⁄∂n on ∂B, involving the non-singular kernel g k .This becomes an integral equation of the first kind for φ in terms of ∂φ ⁄∂n.Copley calls it the interior Helmholtz integral relation.He proves two uniqueness theorems: , then the solution of the functional equation is unique.
if the surface ∂B is axisymmetric and (46) holds for all x along the axis of symmetry in , then the solution is unique.
Copley proposed a scheme for enforcing (46) at a finite number of points {x 1 , x 2 , x 3 ,..} denoted by x, in .However a bad choice of x∈ − B could cause an interference with the interior eigenfunctions , A typical bad choice is when x lies on a nodal surface of eigenfunctions .These nodal surfaces are such that (x) = 0 i.e. ( 47) is the same as ( 46).If this is true for the choice X then the eigenfunctions will also satisfy (46), thereby causing an interference.Cunfare et al (1989) suggested a linear combination of (46) and its normal derivative, but as the points in X do not lie on ∂B the normal derivative operation cannot be well defined.It is also noted that any numerical implementation of this fails as the number of such points increases, which is detailed in the next proposition.Further it can be shown that any numerical scheme based on (46) has undesirable characteristics as the number of points increases as we shall see in the next proposition.Copley solved some axisymmetric radiation problems but noted that his method fails to give accurate results for bodies which are not sufficiently smooth.Fenelon (1969) obtained results for acoustic radiation by a finite cylinder, using Copley's approach.[ ] Schenck showed that for any k there is only one solution to (48a) which also satisfies (48b) for all points in the interior region.Discretisation of (48a) and of (48b) for only a few points in gives a non-square system of linear equations, which is then solved by a least-square method.This method is better known in literature as the Combined Helmholtz Integral Equation Formulation (CHIEF).For lower ranges of values of k, successful applications of CHIEF have been reported (Amini and Harris, 1990;Seybert et al, 1984Seybert et al, , 1985Seybert et al, , 1986Seybert et al, , 1987)).Seybert et al have used CHIEF successfully to model radiation and scattering problems using quadratic isoparametric elements.However, the major drawback of the Schenck method is that there is no way of knowing how many interior points are needed and where best to locate them.If a point is placed on a nodal surface of an interior eigenfunction, the method will still fail for the reason given before.No criterion has yet been reported which enables one to isolate these "bad" interior points.However, since this method is relatively simple to implement, it is often successfully used in practice for modest values of k (Amini and Harris, 1990;Seybert and Rangarajan, 1987).

Kussmaul
This proposal utilizes hybrid layer-potentials.Kussmaul (1969) proposed for a scattered wave, ( 49) It is a superposition of a simple-layer and a double-layer distribution over the boundary.Zaman (1994) proposed an adaptation of this method and successfully obtained some numerical results for higher k.The hard boundary condition on (49) yields This proposition is analogous to the mixed potential used by Brakhage and Werner (1965) and Leis (1958Leis ( ,1962Leis ( ,1965) ) in their treatment of the Dirichlet problem.The constant µ is chosen to be +i if ℜ(k) ≥ 0 and -i if ℜ(k) < 0. This choice ensures the uniqueness of the solution of the equation (50).Since the operator N k is non-compact, so the operator on the right-hand side of (50) becomes non-compact, and as a result, Fredholm theory cannot be applied to deduce whether or not a solution exists.Kussmaul uses a regularization technique, which, by eliminating the non-compactness, enables him to prove the existence of solution.Due to the complexity of the regularization technique few numerical implementations have been reported.Kirkup and Henwood (1992) used this method to compute solutions of some acoustic radiation problems.

Burton and Miller
This proposition is the direct analogue of the indirect formulation of Panich (1965).Burton and Miller (1971) proposed a linear combination of the Surface Helmholtz Equation (19b) and its normal derivative (21b), which in operator notation takes the form, where α is an arbitrary coupling parameter.It is proved (Burton and Miller, 1971;Tzu-Chu, 1984) that if α is chosen to be a complex constant with ℑ(α) > 0 , then the homogeneous equation has only the trivial solution φ = 0 for real values of k, ensuring that (51) has a unique solution for all wavenumbers.Harris (1989), in a similar way to Tzu-Chu (1984), extends the uniqueness proof to the case where α is taken to be a function of x, i.e. for the condition ℑ(α(x)) > 0, ∀x∈∂B it is proved that (52) has only the trivial solution.This result is used in the study of the conditioning of the formulation.Analytically, the proposal is undoubtedly elegant.Nevertheless, it suffers from a numerical drawback due to the involvement of the non-compact (i.e.hypersingular) integral operator N k .To resolve this problem Burton (1974), following Panich (1965), regularizes N k to obtain the Regularized Burton and Miller formulation (RBM): where all the integral operators are weakly singular and compact.The choice of α(x) is such that ℑ(α(x)) > 0 for k > 0. The main drawback of RBM is having to compute the product of operators (by multiplying their matrix approximations) which proves to be very expensive and cumbersome.The reason for the regularized version is twofold.Firstly to establish the uniqueness of solution using classical Fredholm theory, and secondly to make it amenable to numerical treatment despite its complexities.However it seems possible to compute a direct approximation to N k .Terai (1980) gives an approximation to N k for a flat boundary element, and Stallybrass (1967) introduced a pointwise variational principle to approximate N k , which Meyer et al (1978Meyer et al ( , 1979) ) The second integral on the right-hand-side is weakly singular and the first integral can be interpreted in the sense of Cauchy principal value.Amini and Harris (1990) gave a comparison between the RBM and the direct Burton and Miller form. Reut (1985) has developed a Burton and Miller formulation which is called Composite Normal Derivative Overlapping Relation (CONDOR), based in part upon Terai (1980).Ursell (1973Ursell ( , 1978) ) noted that g k (x, y) can be replaced by any fundamental solution to the Helmholtz equation in the exterior domain B + which satisfies the radiation condition.So we have a modified Green's function

Ursell and Jones
where Γ is an analytic wave function in B + .Using this in the classical formulation (Layer-theoretic or Helmholtz formula) with (53) replacing the classical free-space Green's function, and applying the hard boundary condition, the boundary equation takes the form, [Layer-theoretic]: where γ is a complex constant such that ℑ(γ) > 0 and S R is a sphere enclosed entirely within ∂B, then Ursell (1973) shows that the solutions to ( 54) and ( 55) are unique.Ursell noted that Γ(x, y) may be chosen to be an infinite series of spherical wave-functions, which converges quickly for small k but slowly as k increases.Jones (1974) proposed a modification to this by replacing the infinite series with a finite one, i.e.  54) and ( 55) have unique solutions provided that the b mn 's are real and non-zero.Kleinmann and Kress (1983) suggested that b mn be chosen so as to minimize the condition number of the resulting integral operator.The major drawback of this method is that for an arbitrary surface one does not know how large to take M since K is usually unknown.Also it is known (Burton, 1974) that the number of elements in K less than a given k is proportional to k 3 , and so for a moderately large k we have to take a considerably large number of terms.Consequently Γ(x, y) becomes very costly.As such there has been little or no report of its implementation in practical problems.

Piaszczyk and Klosner
This approach (Piaszczyk and Klosner, 1984) is similar to CHIEF, but it uses SHE (19b) together with the exterior Helmholtz relation Naturally φ is not known for points in B + , and for hard scattering ∂φ ⁄∂n only is known on ∂B.Therefore it is supposed that there exists some function Z(x) which gives a simple impedence relationship Using ( 57) and ( 58), one can compute φ for a few points in B + .Utilizing these values and a suitable discretisation of SHE, an overdetermined system for φ in ∂B can be formed.Using a least-squares procedure one can solve this system of linear equations and find φ on ∂B.Then one can recompute φ in B + and repeat the process until a convergence criterion is satisfied.The proposal is analytically sound since a uniqueness proof is presented.However this method also suffers from a major drawback in that it requires the solution of an over-determined system at each iteration, which is very costly from a computational point of view.Also the rate of convergence depends on the choice Z(x) and it is doubtful whether or not the iteration process will converge at all (Cunefare et al, 1989).
Although, the problems which the classical formulations faced are circumvented analytically by the improved approaches, the implementations have been very difficult.Rêgo Silva et al (1993) presented a numerical implementation of a hypersingular boundary element formulation for acoustic radiation problems.Here, Burton and Miller formulation in (51) is used.It is argued that in the hypersingular integral operator N k the behaviour of the integral is determined by conditions at small r = x -yfor which ( ) Clearly the convergence conditions of the integral in N k are the same as those of the normal derivative of a double-layer of a Harmonic (Laplace) potential (with k = 0), given by the Lyapunov-Tauber theorems (Courant and Hilbert, 1989).Following Guiggini et al (1992), the method of Hadamard's (1952) finite part was used here to evaluate the hypersingular integral i.e.Here s ε is the surface exterior to S of a small sphere of radius ε around a point x ∈ S, and e ε is the subregion of S lying inside the small sphere.According to Hadamard's definition, the shapes of s ε and e ε must be preserved while ε → 0.More so, the shape of e ε must be consistent with the shape of s ε throughout the limiting process.

Zaman
Another method (Zaman, 1994;Jaswon and Zaman, 1993) proposed forms an adaptation of the Kussmaul (1969) formulation.Analogous to Kussmaul formulation (48), a superposition of a simple and a double-layer potential is used.The formulation (48) is improved by introducing an internal fictitious surface ∂B * , similar and similarly situated to the given surface ∂B, i.e. for each y ∈∂B there is a corresponding y * ∈∂B * such that y * = ϑy, where ϑ is a real constant such that 0 < ϑ ≤ 1.The solution of the Helmholtz equation in the infinite region exterior to the surface ∂B is represented by a superposition of layer potentials with a coupling parameter η, [ ] is a continuous distribution of modified dipole sources on ∂B * for each ϑ.In (59) the parameter η is chosen to be iϑ -1 so that the coupling parameter eventually takes the form +i. This formulation is called Modified Layer Formulation (MLF).This is a hybrid potential distribution except that the dipole components are located away from the surface.It may be regarded as a classical simple-layer formulation with the freespace Green's function g k replaced by the modified Green's function  ) . Note that G k has the same singularity behaviour as that of g k when |x -y| ≈ 0. The hard acoustic boundary condition (Neumann) applied to (59) yields 60) is a Fredholm integral equation of the second kind for σ in terms of the given term ∂φ ⁄∂n on the left hand side.Also, significant in ( 60) is that the only singularity in the kernel is the weak singularity of ∂g k ⁄∂n x , which presents no computational havoc.Zaman (1996Zaman ( , 1997) ) presents arguments concerning the uniqueness of the solution of (60) and also presents some illustrative works on potential theory applied to a spherical surface.The formulation (59) and the corresponding boundary equation ( 60) are utilized numerically to find some test results comparing the analytical and numerical results for a sphere, cube and cylinder.The radiation problem of a pulsating and oscillating sphere is also performed numerically.Finally the scattering problem from a sphere is performed.It is found that ϑ ≅ 0.5 yields the best agreement with the analytical result.Detailed numerical treatment of the problem is discussed in Zaman (1994).

Conclusion
The purpose of this review has been to investigate different improved boundary integral formulations of scattering / radiation problems designed to overcome the shortcomings of classical formulations.A host of Boundary Element Methods (BEM) has emerged over the last couple of decades for the numerical implementations of the formulations.Numerical aspects of the applications are discussed in many publications.For numerical implementation, one way is the discretisation (Zaman, 1994) process with the collocation method applied to the surface divided into geometrical shaped-elements?with constant on each element, which resulted in a set of several linear equations.The other way is the surface interpolation (Harris, 1989) process which involved spline-interpolation with parametric rectangular or triangular mappings.Burton (1976) suggested that the surface should first be divided into a few surface regions, each of which is mapped parametrically on to a rectangle using bicubic B-spline interpolation.This procedure allows the surface integrals to be expressed as sums over unit square elements and thus leads to a large reduction in the number of quadrature weights needed.The set of linear equations are then solved by a Gaussian-elimination process (Zaman, 1994) or otherwise.The MLF formulation (Zaman, 1993(Zaman, , 1994) ) has been applied to regular axisymmetric bodies as test problems with promising numerical results and its engineering feasibility is still being looked at.There are many papers (Burton and Miller, 1971;Wilton, 1974;Wait, 1973;Van Buren, 1970;Brebbia et al, 1991;Bannerjee and Butterfield, 1981;Ingber and Hickox, 1992) dealing with the feasibilities of other formulations.Atkinson (1976) gives a very useful survey of numerical methods for the solution of boundary integral formulations.Brebbia (1978) gives a good text for BEM.For getting round the problem of fictitious interior eigenvalues, which we come across when wave-numbers are large, Shaw and Huang (1989) suggested an embedding approach.
The mathematical formulations of the problem is reviewed.The factors governing the choice of integral equation formulation were threefold: 1.The integral equation must be numerically stable over the required range of wave-number, preferably without user intervention.2. The integral equation reliable for all surface shapes including long slender bodies.3. The formulation efficient in the use of computer time.
Acoustic shields, screens or barriers are sometimes used as a means for reduction of noise from machines, engines or roads.For determining mathematical model for an acoustic shield, the main property to take into account is that the shield is a thin structure and so its vibratory motion can become coupled to the acoustic field.Thus a shield is equivalent to a discontinuity in the sound pressure field.Due to its topology, the shield is termed a shell.Kirkup (1990Kirkup ( , 1991Kirkup ( , 1992) ) and Kirkup and Henwood (1992) give a detailed study of this.The purpose of their study was to extend the BEM for the solution of 3-D acoustic radiation problems to acoustic-vibratory systems involving vibratory bodies and thin acoustically coupled shells.The collocation method is applied resulting in a method termed Boundary and Shell Element Method (BSEM) which is noted to be a useful generalisation of the traditional BEM.The method was implemented for general 3-D problems in a Fortran subroutine ABSEMGEN.The subroutine is applied to test problems.
Due to numerous recent applications of boundary element methods involving hypersingular integrals, some research has been carried out on developing new sophisticated elements in order to satisfy their existence, as alternatives to the simple but time-consuming discontinuous elements (Rêgo Silva, et al 1993).
Splines C 1,α elements, and B-Splines C 2 elements, are among the most promising options.Nevertheless, further developments regarding their efficient computational implementation are still required.
Problems such as fluid/structure or soil/structure interaction are of importance to many physical and engineering applications, and have been receiving considerable attention as to their numerical modelling.Some approaches discretise the structure using FEM.However, BEM can also be an efficient alternative for these modellings.

Acknowledgement
The facilities provided by the Department of Mathematics and Statistics, Sultan Qaboos University, are gratefully acknowledged by the author.Also, a warm gratitude from the author goes to Professor Ibrahim A. Eltayeb, Sultan Qaboos University, for his encouragement and advice throughout this project.

284 u
Let ρ o, p o be the density and pressure, respectively, of the fluid at rest.At time t, the excess pressure and the particle velocity, denoted by p and , are related by Newton's equation of motion, (I) the behaviour of the layer or the surface density function σ, (II) the singularities of the kernel functions, viz.g k and the nature of the surface ∂B.
Helmholtz potentials are analytic functions of x∈ − B or x∈ B + , the behaviour for x∈ ∂B will depend on the local nature of the surface in the neighbourhood of x.The simple-layer potential, L k , and the double-layer potential, M k , are clearly the solutions of Helmholtz equation, both in and B − B gives a detailed treatment of this.It is shown that the exterior Neumann (Dirichlet) problem breaks down whenever the wave-number k equals the eigenvalue of the interior Dirichlet (Neumann) problem.Let us look at the situation briefly here: Let, K D = { wave-numbers k  interior Dirichlet problem has non-trivial solutions}, K N = { wave-numbers k  interior Neumann problem has non-trivial solutions}.Homogeneous Dirichlet case : Inserting the boundary condition φ = 0 into the boundary equation (20b) and ) simultaneously.Suppose k∈K D .Then the interior homogeneous Dirichlet problem has non-trivial solutions, i.e. φ ≠ 0 in , even though φ = 0 on the boundary.So it follows from (20c) that ∂φ ⁄∂n ≠ 0 on ∂B implying that (28) and (29) have non-trivial solutions.This result has an important implication on the exterior Neumann formulations.− B Homogeneous Neumann case : Inserting the boundary condition ∂φ ⁄∂n = 0 into the boundary equation (20b) and then writing the result in operator form we get the same boundary condition in the (differentiated) boundary equation (21b) and writing the result in the operator form we get values of φ satisfy both (30) and (31) simultaneously.Let k∈K N , the interior homogeneous Neumann problem has non-trivial solutions i.e. φ ≠ 0 in − B even though ∂φ ⁄∂n = 0 on ∂B.
trivial solutions at any failing (characteristic) values of k which are also the eigenvalues of the corresponding interior Dirichlet problem.Therefore the general solution of (39) may be written in the form of k, where σ o is a particular solution, and c j 's are arbitrary constants.Brundrit noted that the non-trivial eigensolutions of (40) are the boundary-values of ∂φ ⁄∂n, for the interior Dirichlet eigensolutions, which satisfy simultaneously This method(Schenck, 1968) utilizes the Surface Helmholtz Equation [see eq.(19b) of previous section], supplemented by the interior Helmholtz relation i.e.

(
utilized in order to show that