COUPLED COMPLEX BOUNDARY METHOD FOR A GEOMETRIC INVERSE SOURCE PROBLEM

. This work deals with a geometric inverse source problem. It consists in recovering the characteristic function of an unknown inclusion based on boundary measurements. We propose a new reconstruction method based on the CCBM and the shape gradient method, the inverse problem is formulated as a shape optimization one, corresponding to a coupled complex boundary state problem. Well posedness and existence results are presented. A computed expression for the shape gradient is used to implement a gradient algorithm. The efficiency and accuracy of the reconstruction algorithm are illustrated by some numerical results, and a comparison between CCBM, Least-squares and Kohn-Vogeluis methods is presented


Introduction
Inverse source problems is a class of inverse problems that aims to find (,  ) solution of the equation  =  in Ω, using a pair of additional data on the boundary Ω of an open bounded set Ω, where  is an elliptic linear differential operator and  is the unknown source term.It is well known that this problem is one of the highly ill-posed problems in Hadamard sense [10], therefore, a general source function cannot be identified uniquely using boundary measurements, see [22] and [13].
A particular source problem where the unknown source  is of the form   ( is the characteristic function) has been studied by Afraites et al. [5] using a shape optimization reformulation, the unknown source support was reconstructed by a gradient algorithm using the shape gradient and the adjoint method, then the stability study was presented.This problem has been studied also by Hrizi and Hassine [20] using a topological optimization formulation, the unknown source was reconstructed using a level-set curve of the topological gradient.Kress and Rundell [21] presented an iterative solution method via boundary integral equations, by reformulating the inverse source problem as an inverse boundary value problem with a non-local Robin condition on the boundary of the source domain.In [27] a method for the reconstruction of star-shaped characteristic sources was developped, by reducing the problem to an algebraic system of equations.The paper [6] traited the inverse characteristic source problem, in the case of Helmholtz equations, from the determination of the barycenter of the characteristic source and the recovery of its geometry from a class of star-shaped characteristic sources, using an algorithm based on an equivalent reciprocity functional formulation.El Badia and Nara [14] have also investigated the inverse source problem of the Helmholtz equation, where the source consists of multiple point sources, an algebraic algorithm was proposed to identify the number, locations and intensities of the point sources from boundary measurements.This inverse source problem was studied also by a coupled complex boundary method (CCBM), originally proposed by Cheng et al. in [9], it consists in recovering the source term  in the equation  =   , when the source support  is known, using some additional Dirichlet and Neumann boundary conditions.The method consists in coupling the two conditions in a single complex Robin boundary condition, and then the boundary data fitting is recast into the whole domain data fitting.The authors have shown that the new CCBM method makes the inverse source problem more robust and more efficient in computations.Our aim here is to apply this CCBM method to a geometric Inverse source problem, we should reconstruct the source support  rather than the function , and we consider the more general case  = ℎ 1   + ℎ 2  Ω∖ω , for some given source functions ℎ 1 and ℎ 2 .
Recently, the idea of using the coupled complex boundary method for solving inverse geometric problem was proposed in [1] for solving inverse obstacle problem, the author have shown that the proposed method is feasible and effective for such problems.
To fix ideas, let Ω be an open, bounded, and connected subset of R  ( = 2 or 3) with  1 boundary Ω.Let  be a positive constant, we define the set of admissible domains denoted by Ω  as the set of all open simply connected subdomains  of Ω with a  2,1 boundary, such that (, Ω) >  for all  ∈ .The notation   denotes the characteristic function of .For some given functions ℎ 1 , ℎ 2 ∈  2 (Ω), and based on some knowledge on the boundary conditions, namely the voltage  ∈  1 2 (Ω) and the current measurement  ∈  − 1 2 (Ω), we try to find  and  solution of the following overdetermined problem where   stands for the outward normal derivative, and  is a non-negative constant.
To be more precise, the inverse source problem is reformulated as follows: Problem 1.1.find  ∈ Ω  and  which satisfy the overdetmined system (1). (2) Note that the right-hand side of the first equation in (1) can be re-expressed as ℎ 2 + (ℎ 1 − ℎ 2 )  .Among the existing methods found in the literature, a robust reconstruction of the unknown inclusion consists in reformulating the inverse problem (1) into the following least-squares minimization problem (see for example [4,15]) where  is the boundary measurement and   is the state function that solves Another well known optimization approach is based on minimizing a Kohn-Vogelius functional (see [3,4,8]), which illustrates a more robust optimization performances, through the following optimization problem and by splitting the over-determined boundary value problem (1) in two auxiliary problems.We denote by   the solution of the first one (4) associated to the Neumann data , and   is the solution of the second one associated to the Dirichlet data  given by We note from the above statements that the minimization problems (3) and (5)  The outline of this paper is organized as follows.In Section 2, we present the new coupled complex boundary method appropriate to our geometric inverse source problem and its reformulation at the shape optimization by introducing the regularized Least Squares fitting for the imaginary part of the complex PDE's solution .In Section 3, we present the well-posedness result of the direct problem.Section 4 is devoted to the existence of an optimal solution of our minimization problem.In Section 5 we establish the shape gradient calculus of the cost function.In the Section 6, we give an algorithm based on the gradient method and we solve an elliptic problem in order to find the steepest descend direction in the space of  1 velocity vector field that satisfies certain boundary conditions.Then, we present detailed numerical results.Finally, we present the conclusion in Section 7 and elaborate the calculation of the shape state derivative in the appendix in Section 8.
From the above discussion, the Problem 1.1 is equivalent to the following inverse problem: where  2 is the imaginary part of the solution  =  1 +  2 of the complex boundary problem (7).
To solve the inverse shape Problem 2.2, we transform it into the following shape optimization functional: and introduce the following minimization problem, find  * ∈ Ω  an admissible shape, such that This problems is unstable under data perturbations.Indeed, the choice of low frequency in measurements data can often be interpreted as a regularization method.But additional regularization is usually recommended to stabilize the numerical algorithm and obtain a satisfactory reconstruction.Thus, we consider the following regularized shape functional where  Ω () is the perimeter of  relative to Ω ( [25], page 48), and  is a regularization parameter (see [2,24] for the choice of the latter).
Remark 2.3.The new method allows us to define the cost function  in the whole domain Ω which brings advantages of robustness in the reconstruction such as the Kohn-Vogelius cost function   compared to the Least Squares fitting   which is defined only on the boundary Ω (see [4,5]).Compared to the Kohn-Vogelius method, the latter requires two problems to be solved at each iteration, however the new method (CCBM), needs a single complex problem to be solved.

Well-posedness result
Before discussing the existence of solution to our problem, we first study the well-posedness of the forward complex problem (i.e., existence, uniqueness, and sensitivity to data).

Existence of an optimal solution
The existence of optimal shape requires some continuity or at least lower semicontinuity of the functional to be minimized.This also implies having continuity of the solution of the associated partial differential equation with respect to the variations of the domain in an adequate topology.In this section, we analyze in detail the continuity of the shape functional and the mapping  →   ∈  1 (Ω), where   is the solution of the complex problem (7) corresponding to a variable open subset  of the fixed open set Ω  .
It is necessary to add extra conditions on the variable domains to expect good convergence of the solutions.We state a classical sufficient condition, in terms of uniform regularity of the domains, the −cone property, and we define the family of admissible domains Proposition 4.1.There exists  ∈  such that Proof.Let (  )  ⊂  be a minimizing sequence such that (  ) −→ inf ().Since  is compact for the topology of characteristic functions (see [25], page 59), then there exist an open set  ⊂  and a subsequence denoted also (  )  ⊂  that converges to  in the sense of characteristic functions.We denote by   ,  the solutions of the complex problem (7) corresponding respectively to   , .From ( 14) we have, This shows that the sequence   is bounded in  1 (Ω).Up to a subsequence, it may be assumed that   converges weakly in  1 (Ω) and strongly in  2 (Ω) to a function  * ∈  1 (Ω).Let us show that  * satisfies the variational formulation of problem (7) corresponding to .By definition of   , for all  ∈  1 (Ω), We know that   ⇀  * in  1 (Ω), then ∇  ⇀ ∇ * in  2 (Ω), we get By compacity of the trace from  1 (Ω) to  2 (Ω), we can extract a subsequence   →  * in  2 (Ω), we deduce that Therefore, for all  ∈  1 (Ω), This proves that  * = ().
Let us denote by   2 the imaginary part of   , and  2 the imaginary part of .We have ‖  ‖ 0,Ω −→ ‖‖ 0,Ω , and therefore Since the perimeter is lower semicontinuous with respect to the convergence of characteristic functions (see [25], page 51), We can infer that () ≤ lim inf (  ).
Which implies that  is a minimizer of (11).

Shape derivative calculus
In order to use a descent method of gradient type, it is necessary to differentiate the shape functional  with respect to the shape variable .To define the shape derivative, we will use the so called velocity method, which is, for instance, introduced in [23].To this end, we consider the variation of a given smooth reference shape  according to the displacement  , defined by   () =  +  (), where  a smooth vector field with compact support in Ω  , and we denote by  the space of admissible deformations  .It is well known that  is invertible, for sufficiently small .
We say that the functional () has an Eulerian derivative at  in the direction  if the limit lim →0 (  ()) − ()  := ().exists.Furthermore, if the mapping  ↦ −→ ().= ⟨∇,  ⟩ is linear and continuous, we say that  is shape differentiable at .When  has an Eulerian derivative, we say that ∇ is the shape gradient of  at .For more details concerning the differentiation with respect to the domain, we refer to the books [12,18,25].

Shape gradient of the cost function
Computing the shape gradient of the cost functional requires to calculate also the shape derivative of the state (7).Numericaly, it implies that we need to solve as many boundary problems as discrete shape variables.To avoid this extra computational cost, we use the classical adjoint state method, which requires to solve only one extra boundary value problem.
By introducing the suitable adjoint problem, we give the shape gradient of the cost function in the following proposition.
Proposition 5.1.Let  ∈ , the cost functional  is differentiable with respect to the shape  in the direction of  and its shape derivative is given by where  =  Γ (), and  =  1 +  2 solves the following adjoint problem Proof.First, we derive the following expression using Hadamard and divergence formulas Where  Γ denotes the tangential divergence on , and  is the mean curvature of .
On the other hand, one can show that the shape derivative of the state solves the following boundary value problem (see Appendix) Then, the weak formulation of ( 21) with  the solution of (20) as a test function, and the weak formulation of (20) with  ′ as a test function, are given by (23) implies that ∫︁ By ( 22) and ( 24), it follows Finally, we conclude that

Algorithm and numerical results
In this section, we present some numerical simulations in order to confirm and complete our previous theoretical results with a comparison between CCBM, Least-squares and Kohn-Vogelius Methods.In order to solve numerically the optimization problem (12), we opt for the classical shape variation descent algorithm.First of all, we describe the algorithm and the framework used, then, we present the numerical simulations and some comparisons.

Algorithm
The shape derivative of the cost function  along a deformation field  can be expressed as where  =  2 (ℎ 2 − ℎ 1 ) +  where  2 is the imaginary part of  the solution of the adjoint problem (20).The deformation field  is chosen to provide a descent direction of the cost function (), thus  = −.on  is a descent direction.In addition, it is well known that the shape gradient is defined on the boundary of the moving shape [28], using this approach, the direction of descent must be defined only on .However, if the boundary measurements (, ) is not sufficiently smooth, the surface expression of the shape gradient may not exist or the direction of descent  may have a poor regularity.Therefore, it is interesting to derive a direction of descent  on Ω from the volumetric expression of the shape gradient.Which requires solving another additional variational problem.Let  be the Riesz representative of −∇(), i.e. (see [7] and [16]) where and <, > is the inner product on  2 defined by The equation ( 25) is the week formulation for the following system We give in the following algorithm the gradient method of our problem.

Numerical results
For the numerical simulations, we consider the dimension two and we use the finite elements Software Freefem++ (see [19]).The exterior boundary Ω is assumed to be the square [−1, 1] × [−1, 1].We construct the synthetic data on Ω, by fixing the shape  and choosing the Neumann boundary condition () = sin(),  ∈ [0, 2], then, we compute the trace of state  solution of (4) to extract the measurement  =   Ω.For the latter equation, we use a  2 finite elements discretization to solve the direct problem.The examples with noisy data are generated by perturbing the Dirichlet data  using a fixed amplitude of  Gaussian noise.The synthetic data has been chosen, in a way to avoid the so-called inverse crime phenomena.To this end, the size of discretization used to obtain the synthetic data is different from the one used for solving the inverse problem.We precise that in all tests, the exterior boundary is represented by the black line, the initial shape by the green line, the exact shape to identify by the blue line and the reconstructed shape by the dotted red line.

Results without noise
we present in Figures 1 and 3, the reconstruction by three cost functionals, in the left, the identification by the Kohn-Vogelius cost function, in the middle the Least-squares, and in the right, the result obtained by CCBM.we notice that the results obtained are effective and similar.
In Figures 2 and 4, we plot the variation of the cost functionals according to the number of iterations and the evolution of their associated gradients also.Clearly, The CCBM converges much quicker than the two other methods.
In Figures 5 and 7, we observe that the results obtained by the CCBM and Kohn-Vogelius is more robust than those obtained by Least-squares.In Figures 6 and 8, the convergence of the CCBM dominates that of the Kohn-Vogelius and Least-squares.

Results with Gaussian noise
In this subsection, we present the results obtained from the noisy data as follows where  is a uniformly distributed random variable in [−1, 1] and  dictates the level of noise.We present in Figure 9, the results with different levels of noise.We observe that for different geometries the CCBM gives good results.

Results with impulse noise
Having shown the stability of the proposed method with noised data infected with Gaussian noise, we present in this section the robustness of this method through a more complicated type of noise.We suppose that the given data is infected by the impulse noise.Impulse noise consists of relatively short duration "on/off" noise pulses, caused by a variety of sources, dropouts, or surface degradation.The impulse noise is not trivial due to its complex statistical nature.To construct noised data, we select randomly a number of nodes in the boundary Ω with a percentage .Then we add the impulsive noise to data as follows: where  is the magnitude of corruption, ‖  • ‖ ∞ is the maximum norm.While  is a vector of impulsive noise which have the following probability density function where (  ()) is the Kronecker delta function and D  (  ()) is a zero-mean Gaussian probability density function.While  1 and  2 are parameters that control the mixture of a discrete probability mass function (  ()) and a continuous probability density function.
In Figure 10, we present the plot of different obtained boundaries by the algorithm 1, using data with different rates of impulsive noise.For the case of 1% rate of impulsive noise, the percentage of corrupted data points is  = 30% and the magnitude of corruption is  = 0.1.In the case of 2% rate of impulsive noise, we have taken  = 50% as percentage of corrupted data points and  = 0.2 as magnitude of corruption.While in the last case, all nodes of boundary Ω was affected by noise which corresponds to  = 100% and the magnitude of corruption was  = 0.5.As we can see in this case of more complicated noise, the proposed approach still present a good approximation which prove numerically its stability.

Real case without exact solution
In this section we consider a case with data from real life situation.We consider an electrostatic problem.We consider a rectangular domain Ω =] − 1, 1[ 2 .Denote by Ψ the source term of electrostatic charges which supposed discontinued through the inclusion  and let (,  ) be given data on the boundary.Our goal is to reconstruct surface charges  that creating the potential field  from Cauchy data (,  ).The Cauchy data are taken as follows : the measured potential  on the boundary is presented in the Figure 11. the output flux  = 0.While the source ℎ 1 = 12 and ℎ 2 = 230.The mesh of the optimal domain is presented in the Figure 12.In the Figure 13, we present the evolution of the cost functional  and the norm of the gradient |∇| according to iterations.In the Figure 14, we present the potential in optimal domain for 2 and 3 dimensions.
As we can be seen, the obtained results in the real-life case show the performance of the proposed approach.

Conclusion
We have presented a regularized Complex coupled formulation for the identification of the geometric source problem.The inverse problem is reformulated as a shape optimization one.The existence of a minimizer is investigated.A numerical algorithm for solving the proposed optimization problem is developed.We also derived an exact computation of the gradient of the cost function of our optimization problem.The numerical experiments    We denote by   the solution of ( 7) with inclusion   =   ().For the computation of the shape gradient of the state, we shall use the classical technique which consists of transporting the quantity   defined in the variable domain   back onto the reference domain  using the following transformation   =   ∘   , the usual methods of differential calculus can now be applied since both functionals   and  are now defined in the fixed domain .The material derivative (or Lagrangian derivative) of the state is then given by u := lim The shape derivative (or Eulerian derivative) of the state is defined by  ′ := u − ∇..
Proof.The proof is divided into four parts.Firstly, we transfer the perturbed problem to the fixed domain, we then prove weak convergence of a given sequence to the material derivative of the state, next, we show its stong convergence, and finally, we deduce the shape derivative.

Figure 2 .
Figure 2. Evolution of the three cost functions and their gradients according to the number of iterations.(a) Cost functions.(b) Gradients.

Figure 4 .
Figure 4. Evolution of the three cost functions, and their associated gradients according to the number of iterations.(a) Cost functions.(b) Gradients.

Figure 6 .
Figure 6.Evolution of the three cost functions and their gradients according to the number of iterations.(a) Cost functions.(b) Gradients.

Figure 8 .
Figure 8. Evolution of the three cost functions and their gradients according to the number of iterations.(a) Cost functions.(b) Gradients.

Figure 9 .
Figure 9. Reconstructions of the unknown domain with noise.

Figure 10 .
Figure 10.Plot of different boundaries obtained with noised data with different rates of impulsive noise.

Figure 11 .
Figure 11.The measured potential on different parts of the boundary Ω.

Figure 12 .
Figure 12.The mesh of the optimal configuration of domain showing the inclusion .

Figure 13 .
Figure 13.The evolution of the cost function and the norm of the gradient.

Figure 14 .
Figure 14.The plot of the optimal potential in 2 dimensions and 3 dimensions.
[26]the Neumann data  and the Dirichlet data  sequentially.In this paper, we propose a new coupled complex boundary method (CCBM) that uses both  and  data in a single PDE.The idea of the CCBM is to couple the Neumann data and Dirichlet data in a Robin boundary condition in such a way that the Neumann data and Dirichlet data are the real part and imaginary part of the Robin boundary condition, respectively.As a result, the data needed to fit defined on the boundary Ω are transferred to the volume problem defined on Ω.The coupled complex boundary method (CCBM) was first proposed byCheng et al. in ([9]) for solving an inverse source problem, Rongfang et al. in([17]applied it to an inverse conductivity problem with one measurement.More recently, this method was also applied in solving inverse obstacle problems by Afraites in[1]and used for solving stationary free boundary problems by Rabago[26].To the best of our knowledge, in the literature, this is the first time that the idea of the coupled complex boundary condition has been explored for solving geometric inverse source problem.Unlike the classical methods, the CCBM method allows us to define the cost function  in the whole domain Ω which brings advantages of robustness in the reconstruction compared to the Least-Squares type functions which are defined only on the boundary.From a theoretical perspective, it is difficult to prove the superiority of the proposed Complex coupled formulation over the Kohn-Vogeluis one.However, numerical examples in Section 6 indicate that it leads to more robust reconstruction results for this geometric inverse source problem.