Abstract
We consider the inverse problem that consists in the determination of characteristic sources, in the modified and classical Helmholtz equations, based on external boundary measurements. We identify the location of the barycenter establishing a simple formula for symmetric shapes, which also holds for the determination of a single source point. We use this for the reconstruction of the characteristic source, based on the Method of Fundamental Solutions (MFS). The MFS is also applied as a solver for the direct problem, using an equivalent formulation as a jump or transmission problem. As a solver for the inverse problem, we may apply minimization using an equivalent reciprocity functional formulation. Numerical experiments with the barycenter and the boundary reconstructions are presented.
1. Introduction
In this paper we focus on the inverse problem that consists in the determination of characteristic sources. Inverse source problems have been addressed since 1938, in the context of gravimetry, with the seminal work of Novikov [1], where he was able to prove uniqueness in the two-dimensional Laplace equation model, for star-shaped characteristic sources. For more general geometries, it has been shown by Sakai [2] that there are connected domains for which the problem of reconstruction from the potential outside the source region presents no uniqueness.
Here we consider the inverse problem for characteristic sources, but in the case of Helmholtz equations. When the center is known, uniqueness was proven for a convex characteristic source (see Isakov [3]) with no boundary restrictions. However the proof of Novikov is no longer available by the lack of minimum-maximum principles, and uniqueness has not been established for star-shaped characteristic sources.
While considering Helmholtz equations in a bounded domain this will change the problem in two different ways. If we consider the modified Helmholtz equation, useful for the transient heat problem, the associated boundary value problem has a unique solution, but the decay is exponential. This compromises the quality of the recovery in the inverse problem. On the other hand, considering the classical Helmholtz equation, the models wave phenomena such as acoustic waves, we must avoid values associated with the resonance frequencies, in which the direct problem has no unique solution.
In this paper we focus on two main issues related with the inverse characteristic source problem, that is, the determination of the barycenter of the characteristic source and the recovery of its geometry from a class of star-shaped characteristic sources. As a solver for the inverse problem, we may apply a minimization algorithm using an equivalent reciprocity functional formulation.
In general, there are restrictions on the number and type of sources that can be identified from boundary data. Several particular classes of problems have been considered: combinations of point sources—see [4–6]; linear/affine classes as in [7, 8]; classes of characteristic sources (e.g., [9–16]); and in particular for the Helmholtz equation we refer to the papers [17, 18], where a full identification result was established, but using instead an interval of frequencies.
In this paper, we only consider boundary measurements given by a single frequency.
The paper is organized as follows. In Section 2, we present the direct and inverse problems for Helmholtz equations. In Section 2.1 the inverse jump problem is also presented and we show the equivalence between this problem and the inverse characteristic source problem, which motivates the adoption of the Method of Fundamental Solutions to solve the direct problem (see Section 2.2). The MFS has been applied to solve several inverse problems in recent years (see [19] for an account of this research).
In Section 3 we obtain a new formula for the determination of the barycenter for characteristic sources with symmetric star-shaped support and show that it is the same as obtained for the point source inverse problem. In Section 3.3 the boundary reconstruction is presented as a minimization problem, where the classical Levenberg-Marquardt method is used but the minimization is not considered in terms of the boundary data directly, but in terms of the reciprocity functional.
Numerical experiments related to the barycenter and boundary reconstruction are presented in Section 4. In simulations, we solve the direct problems for different numerical parameters and generate the synthetic experimental data, avoiding the inverse crimes. Beside this, a barycenter formula presented in (46) is tested in Section 4.1. Section 4.2 is entirely dedicated to numerical experiments related to the boundary reconstruction, for the modified and classical Helmholtz equations.
2. The Helmholtz Inverse Source Problem
Let be an open, simply connected and bounded set with boundary, . Given the source term and the Dirichlet data , consider the following problem:where we are considering the cases (Laplace equation), (modified Helmholtz equation), and the usual Helmholtz equation , with denoting the wave number.
Problem (1) has unique solution, except for some (eigenvalues of the Dirichlet-Laplace operator for ) [20], and we can define the direct problem as the problem of finding , with , from the source term and Dirichlet data , where is the normal outward unity vector. On the other hand, the inverse source problem consists in, given the Cauchy data , finding the source term .
For test functions , we introduce the reciprocity functionaland using Green’s second identity, we getfor all . Therefore,with .
In next theorem we establish an equivalence between the source reconstruction from the Cauchy data and from the reciprocity functional.
Theorem 1. The Cauchy data uniquely determines the source if, and only if, is uniquely determined by , for all .
Proof. Consider the sources and for problem (1). Then Imposing the same Dirichlet data on the boundary, for all . Since is homeomorphic to (e.g., [17]), we have Therefore, , in , is equivalent to , in .
2.1. Equivalence to an Inverse Jump Problem
We now show that the inverse source problem for characteristic source is equivalent to an inverse jump problem.
Consider problem (1) with source , where is constant and is an admissible source set. Note that we can rewrite this problem as the transmission problemwhere denotes the jump between the inner and outer traces of the solution, respectively. Let be a particular solution of the equation Notice that if , we can take and if , we can take . So, we can split (8) into two problems andThus, consideringwe have that, up to , problem (8) is equivalent to the inverse jump problemThus, denoting the reciprocity functional of this problem as , we establish the following result.
Theorem 2. If is connected, then . Therefore, the inverse source problem of characteristic source (1) is equivalent to inverse jump problem (13).
Proof. The reciprocity functional for the inverse jump problem, with , is given by, applying the Green formula to the connected set , Notice that the orientation of the normal changes, by taking as the boundary of , instead of as a part of the boundary of . Also, we haveThen, for all , using (2).
Remark 3. If was not connected, this equivalence would not hold, in view of
This jump or transmission problem will be solved using the Method of Fundamental Solutions (MFS) that we present next.
2.2. Solving the Direct Problem with MFS
The Method of Fundamental Solutions (MFS) is a meshless boundary method that has been widely studied in recent years. Its main goal is to approximate the solution of certain boundary value problems by a combination of fundamental solutions of the partial differential equation, in which the singularities lie outside the domain (e.g., [21]).
Furthermore, this method has become a popular tool because of its applicability in several areas where we can cite, as an example, [7, 8, 19, 22, 23], including inverse problems.
Consider the direct problemwhere the source term is , with constant and an open, connected, bounded subset with boundary , .
Using Section 2.1, problem (17) can be rewritten as the transmission problem for (13), whereand is given by (12). In this way, the functions and will be determined by the MFS. For this, consider the following expansions centered at the corresponding source points: with , and with and , where is a fundamental solution of Laplace or Helmholtz equation; that is, where is the Dirac distribution. Beside this, fictitious boundaries , and are chosen such that (i)the curve is outer of , with ;(ii)the curve is outer of , with ;(iii)the curve is outer of , with .
So, is a solution of (13) if it satisfies the following boundary conditions.(i), on ; that is, for , (ii), on ; that is, (iii), on ; that is,In this way, taking collocation points , with , and , with , we can write the above problem as a linear system: where and the matrix M, , of the linear system coefficients, is given by
In general, it is considered that and the system should be considered in the least squares sense. In the 2D case, the fundamental solution can be taken aswhere is a Hänkel function and are Bessel functions.
Therefore, solving the linear system, in the least squares sense, by Tikhonov regularization ( being a regularization parameter), we find a solution, , of (13), in which we may generate an approximation of the direct solution This approximation will be used to generate the synthetic data in Section 3.3.
Remark 4. In the numerical experiments of Section 4, the fictitious boundaries , , and were chosen as circles each one in the correspondent region. The choice of the location of source points is not unique and its best location is a research subject (see, e.g., [24–27]).
3. Barycenter and Point Source Determination
We consider the special class of admissible sets.
Definition 5. A subset is said to be star-shaped with respect to the origin, if there exists a function , such that where is the positive radius function.
We say that a subset is a symmetric star-shaped set if it is star-shaped and the function satisfies When the barycenter is not at the origin, we can use the notation to refer to a shifted domain.
In this section we establish a result for the determination of the barycenter, in the case of symmetric star-shaped characteristic sources for the Helmholtz equation, with no restriction on the sign of , and we bound the reciprocity functional.
In Section 3.2 we see that the same formula also holds for the determination of a point source.
3.1. Determination of the Barycenter
In the case (Laplace equation), since and are functions of , it is well known that the barycenter is determined from where and (for any weight function ).
Consider now , and (therefore ).
Using a test function , where and are arbitrary, the reciprocity functional is On the other hand, by (3), we have, denoting ,
In the following, we will suppose that the origin of the coordinate system, , is the barycenter of an admissible source set , in which we suppose symmetric star-shaped set, whose boundary, , is parametrized by a function .
Therefore, by (35), for all .
Proposition 6. Let be a symmetrical star-shaped subset, whose boundary, , is parametrized by the function . Then
Proof. In fact, let be an arbitrary direction. Consider the sets So, and thenOn the other hand, note thatTherefore, since , , and the star-shaped set is symmetrical, then and we obtain In this way, since is arbitrary, we have , for all .
So, by proposition 6 and by (35), Therefore, that is,
Remark 7. Note that in this case we have used the boundary parametrization associated with the symmetric star-shaped characteristic sources in the proof.
Thereby, we establish the following new result for the barycenter identification of symmetrical star-shaped supports in the inverse source problem posed on the modified Helmholtz operator.
Theorem 8. Let be a symmetrical star-shaped subset and let be the source term for the inverse problem. Then the barycenter, , of is uniquely determined bywhere .
3.2. Determination of Point Sources
In this subsection we will study the inverse source problem for a single point source. Consider the inverse problem with a single point source of type where is the Dirac delta distribution centered at the point and is the intensity of this source. In this way, taking in (4), with , Thus, since , the source intensity is determined in the following way:On the other hand, the source position is obtained fromwhere . The coordinates of the point source are simply given by , with being the canonical basis vector.
Remark 9. Note that, by (50) and (46), we obtain the same expression for the barycenter of symmetric star-shaped sources and for the position of a single point source.
3.3. Determination of the Source Support Boundary from Cauchy Data
Consider the subsets , with same barycenter, whose boundaries, and , are parametrized, respectively, by and by a truncated Fourier serieswith . Here, we suppose that the barycenter of is the origin.
Remark 10. In an arbitrary coordinate system, where the barycenter of is not the origin, the boundary points are given bywith . The barycenter position will be calculated using the reconstruction formula (46); see Section 3.1.
Remark 11. In , is parametrized, in the barycenter coordinate system, by truncated -dimensional spherical harmonics series of the radius , which express the distance between the centroid and the correspondent source boundary point.
Using results from Section 3.1, we can assume that the sources have the same barycenter. Theorem 2 shows us how to implement the boundary reconstruction algorithm using solutions of direct problem (1) with or equivalently solving interface problem (13) with traces at the source boundary of equivalent particular solution.
Remark 12. The normal trace of trial solutions with parameters guess is calculated solving equivalent interface problem (13) with the corresponding trials. As we have pointed out before, the Method of Fundamental Solution is appropriated for the Helmholtz equation model.
Related to these subsets and , we have the associated reciprocity functionals and , for .
Considering, in , test functions and taking different directions , , we define the following discrete functional to be minimized:
Discrete functional (53) represents a possible misfit between the reciprocity functional related to the sources and . Furthermore, note that(i) the given depends only on Cauchy data;(ii) the guessed can be calculated by solving Laplace-Helmholtz direct problem (17) with a trial source term guessed for parameters .We present a numerical optimization method to approximate the boundary support which is based on a Levenberg-Marquardt algorithm implementation for the nonlinear least squares formulation given by (53).
Remark 13. In problems, we can define the functional to be minimized asvalid for test functions , with
In the next sections we present some bidimensional numerical experiments for the minimization of functional (53), defined above, using the Levenberg-Marquardt algorithm.
4. Numerical Experiments
When we study inverse source problems with characteristic sources defined by star-shaped sets, there are two pieces of information that we aim to retrieve: the barycenter of the source support and the boundary parametrization.
4.1. Numerical Experiments for Barycenter Determination
The barycenter reconstruction problem was studied in Section 3.1, and the main objective is to verify the accuracy of the new formula proposed in Theorem 8.
In the numerical Examples 1, 2; and 3 we considered . In Example 4, was considered.
Beside this, relative noise over the collocation points of domain was imposed; that is, the boundary was discretized and the reciprocity functional was determined with the Dirichlet data, in each element of this point list. This value was multiplied by , where was given by a random function with uniform probability distribution.
Example 1. In this bidimensional experiment we consider circular geometries, with . In this case and , using the notation , where stands for the radius and for the center.
The direct problem was solved using the Method of Fundamental Solutions (MFS), where the projection of a solution on the boundary was taken as Dirichlet data, considering collocation points and point sources in the MFS. For the determination of the barycenter, the reciprocity functional was calculated considering the Cauchy data from the solution of the direct problem, . The Dirichlet and Neumann data were given by the trace and the normal trace , on
Thus, the influence of relative noise of , and in the Cauchy data was studied, and the results are shown in Table 1. In this case, it was considered that Dirichlet and Neumann data have the same noise level.
In Table 2, noise in the Dirichlet data was fixed, and the noise in the Neumann data ranges up to . On the other hand, in Table 3, noise in the Dirichlet data was fixed and the noise in the Neumann data ranges up to . The results were worse than when noise was considered in the Dirichlet data.
Example 2. In this experiment we considered the domain as the interior of the ellipse parametrized by , , and the support of the source as the star-shaped set , with barycenter , and boundary parametrized by , , as shown in Figure 1.
In this problem, we consider null Dirichlet data and the direct problem was solved by the MFS (again with 629 collocation points and 252 sources).
The Dirichlet and Neumann data were perturbed with the same amount of relative noise. The results are given in Table 4.
Example 3. In this experiment, the domain was considered as the unitary circle, centered at the origin, and the support of the source as the interior of the square, centered at , with parametrizationas shown in Figure 2.
We again considered null Dirichlet data in the direct problem and solved it by the MFS with only collocation points and source points, and the results are shown in Table 5.
Example 4. In this experiment, the domain was considered as the interior of unitary circle centered at the origin and the support of source as the star-shaped set, with barycenter , whose boundary is parametrized by , as shown in Figure 3.
The main difference between this experiment and the others is the fact that here we consider ; that is, , a complex number. Consequently, the fundamental solution considered in the MFS method is given by where is the Hänkel function of first kind of zero order, and the results include neglectible complex numbers. As in the last experiment, we consider collocation points and source points, with null Dirichlet data in direct problem. Again the obtained results were good, as shown in Table 6.
All these experiments show that the formula proposed for determination of centroid proved to be very efficient, even with considerable amounts of noise level.
4.2. Numerical Experiments for Boundary Reconstruction
In this section we study the boundary reconstruction problem for inverse characteristic source problem. The boundary support reconstruction has been studied by many authors through several numerical methods to estimate this unknown parametrization. In this work, we use the Method of Fundamental Solution (MFS) to estimate the solution of direct problem and the Levenberg-Marquardt algorithm to minimize the functional (53) and estimate the support
4.2.1. Initial Guess for the Parameters
Another feature that is possible to obtain is the effective or mean radius, , the first term in the Fourier series parametrization. A guess for this parameter can be obtained by changing the test function for modified Helmholtz problem by the test function for the Laplace problem identically equal 1 (when is small). In this way, where and are the reciprocity functional for modified Helmholtz inverse characteristic source problem, with test function , and for Laplace inverse point source problem, with test function identically equal 1, respectively. So, where is the Lebesgue measure of the subset . In this way, we obtain the following approximation for the bidimensional case
The initial guess for the other parameters are taken as zero; that is, we adopted , as starting parameter value in the iterative minimization algorithm of least squares functional (53) done using the Levenberg-Marquardt method. The algorithm has been checked for different numbers of parameters with good performance. In these experiments we adopted .
4.2.2. Kinds of Noise Considered
In the following experiments, two kinds of noise were considered: the relative and the absolute.
The relative and absolute noise are imposed over the measurements of the discretized boundary of domain, that is, on the MFS collocation points. In this way, the Cauchy data, in each element of this list, is multiplied by random functions with uniform probability distribution, , and If is the measured data on the boundary, then
In Section 4.2.3 we implemented two experiments considering In experiment 5, was considered the influence of relative noise, but with a variation, domain variation, and support variation. In experiment 6, the square was considered as domain, without the influence of noise. In Section 4.2.4, , that is, complex, and the action of both kinds of noise were considered.
4.2.3. Numerical Experiments for Boundary Reconstruction: Case
In the next experiments, was considered and the damping factor, in Levenberg-Marquardt algorithm, was considered initially equal to and, for each iteration, it is multiplied by a decay factor . The initial values of damping factor and the decay factor were obtained empirically, after several experiments. Beside this, in all experiments directions were taken, to minimize functional (53).
Furthermore, all experiments were implemented in Mathematica software, version 9, in a notebook with processor Intel Core i7, 2.8 GHz, RAM memory 4 GB, 1333 MHz DDR3 with operational system OS X 10.9.5.
Example 5. In this experiment, as a unitary circle and a support were considered, with center in , whose boundary, , is parametrized by , as shown in Figure 4(a), where the Dirichlet data was considered null.
Beside this, the number of parameters in Fourier expansion was considered as 7; that is, . In the original problem the MFS was considered with collocation points and source points, and for the MFS reconstruction collocation points and source points were considered. The initial parameter vector in reconstruction was , with mean radius , calculated by (59).
After 7 iterations (of approximately seconds each), the reconstruction was obtained without noise shown in Figure 4(b), where the figure in red is the original and in green is the reconstructed one.
Beside this, relative noise of , and was also considered in the Neumann data. In all cases good reconstructions were obtained in the fifth iteration, as shown in Figure 5 (again, the figure in red is the original support, in green is the best reconstruction, and in grey is the reconstructions for each minimization iteration).
Another interesting feature is the influence of In experiments related to Figure 5, was considered. In Figure 6 we present the reconstructions with , , , and , and as the value of increases, the value of mean radius also increases, as it can be seen in Figure 6(d). The reconstructions were again very efficient, even with considerable noise.
An experiment with the same kind of domain and source support was also implemented, but with different size. In Figure 8, the domain is a circle with radius 2, and the source support boundary is centered at with parametrization . The reconstruction, without noise, was very efficient (beyond the third iteration there is no visual difference), as shown in Figure 7(b).
Another experiment was done considering a source support having the parametrization, given by , closer to the domain boundary, as can be seen in Figure 8(a). This proximity did not influence the reconstruction, as it can be seen in Figure 8(b) (without noise) or in Figure 9, with relative noise of , , , and even .
(a) Domain considered
(b) Reconstruction, after 7 iterations, no noise
(a) 1%
(b) 5%
(c) 10%
(d) 30%
(a)
(b)
(c)
(d)
(a) Domain considered
(b) Reconstruction, after 7 iterations, no noise
(a) Domain considered
(b) Reconstruction, after 7 iterations, no noise
(a) 1%
(b) 5%
(c) 10%
(d) 50%
Example 6. In this experiment, the original source support, , was considered as the same square centered in , and the domain, , was considered the unitary circle, as shown in Figure 10(a), and the Dirichlet data was considered null.
Again was considered, but in this case we took . After iterations, with approximately seconds each, the reconstruction obtained, without noise, is shown in Figure 10(b).
(a) Domain considered
(b) Reconstruction, after 10 iterations, no noise
4.2.4. Numerical Experiments for Boundary Reconstruction: Case
Now we considered and .
Example 7. In this experiment, the source support, , was considered with center and boundary parametrization , inside an unitary circle domain , as shown in Figure 11(a), again with null Dirichlet boundary data.
The number of parameters in Fourier expansion was taken to be . In the direct problem, collocation points and source points were taken for the MFS, and collocation points and source points were used in the MFS iterations.
The mean radius, calculated by (59), was . So, we considered the initial parameter vector as , where denotes the real part of .
After 7 iterations, of approximately seconds each, the reconstruction was obtained without noise shown in Figure 11(b), where the figure in red is the original and the figure in green is the reconstructed one. It can be observed that from the third iteration the results are very close.
Although the time in each iteration was larger, the reconstruction was better than with experiments with real.
We had considered both kinds of noise in the Neumann measurements (relative noise of , and and absolute noise of and
In experiments related to relative noise only 3 iterations were considered, due to the fast convergence. The results can be seen in Figure 12 (the figures in red are original and the figures in green are the reconstructed ones).
In experiments with absolute noise, the method also converged quickly, after 3 iterations, but the reconstruction only maintained some meaning for a small amount of absolute noise (with 5% the meaning was lost). The results can be seen in Figure 13.
(a) Domain considered
(b) Reconstruction, after 7 iterations, no noise
(a) 1%
(b) 10%
(c) 30%
(a) 1%
(b) 5%
Remark 14. Note that to avoid problems known as “inverse crimes,” the MFS used for the simulation of the direct problem used different points from the MFS used for the inverse iterations. Furthermore, beside the addition of absolute and random noise, notice that in the square support example the boundary square shape could not be reproduced by a finite Fourier series expansion.
5. Conclusions
In this work, we studied the inverse characteristic source problem for Helmholtz equations. Firstly we proved that the associated direct problem is equivalent to a direct jump problem. This equivalence justifies the use of the Method of Fundamental Solution. Beside this, considering support of a characteristic source as a bounded symmetric star-shaped subset of the domain, we present a formula for the barycenter of and a shape reconstruction method. Numerical experiments related to barycenter determination and boundary reconstruction have shown very accurate results.
Competing Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work was partially financed by Brazilian agencies CNPq, Process no. 141829/2012-5, and CAPES, Process no. 99999.007941/2014-5.