Solving inverse electromagnetic scattering problems via domain derivatives

We employ domain derivatives to solve inverse electromagnetic scattering problems for perfect conducting or for penetrable obstacles. Using a variational approach, the derivative of the scattered field with respect to boundary variations is characterized as the solution of a boundary value problem of the same type as the original scattering problem. The inverse scattering problem of reconstructing the scatterer from far field measurements for a single incident field can thus be solved via a regularized iterative Newton scheme. Both the original forward problem and the problem characterizing the domain derivative are formulated as boundary integral equations and we carefully describe how these formulations are obtained in the case of Lipschitz domains. The integral equations are solved using the boundary element library Bempp. A number of numerical examples of shape reconstructions are presented.


Introduction
The use of iterative methods based on derivatives with respect to variations of the boundary has been a standard tool for solving inverse problems of shape reconstruction for many years. For many types of boundary value problems, characterizations of such shape derivatives have been given, either based on variational approaches or using boundary integral equations, and these characterizations have been employed with great success for computing derivatives in actual implementations of reconstruction algorithms. See e.g. chapters 4 and 5 in [5] and [10] and the many references given by these authors for an overview on what has been accomplished in this direction, and also [7] for examples of actual reconstructions obtained by this approach in three-dimensional acoustic wave scattering.
However, much less has been achieved so far for time-harmonic electromagnetic wave scattering problems. The difficulties originate from the much more complicated regularity theory of solutions to Maxwell's equations as compared to solutions of elliptic PDEs. Only in recent years have characterizations of shape derivatives been found [6,9] that are suitable for implementations.
The situation under consideration in this paper is that of an electromagnetic wave in vacuum being scattered by either a perfectly conducting or penetrable homogeneous obstacle. In this case, the total electromagnetic field is a solution to a boundary value problem consisting of the Maxwell system with constant coefficients, the Silver-Müller radiation condition for the scattered field and a boundary or transmission condition on the boundary of the obstacle. It turns out that the derivative of the far field of the scattered field with respect to variations of the boundary only depends on the domain derivative of the scattered field which solves essentially the same type of boundary or transmission problem, only with a different inhomogeneity in the boundary or transmission condition. Thus, the domain derivative can be computed by essentially the same numerical methods as the scattered field itself. However, the conditions that characterize the domain derivative involve traces and surface derivative operators that are not implemented in most libraries for computational electromagnetism.
One numerical method that is well suited to the problem type we are considering is the Galerkin boundary element method. The numerical library Bempp (https://bempp.com) [16,17] provides the necessary implementations of boundary element spaces, potentials, int egral operators and Calderon-based preconditioners to efficiently solve such scattering problems. Moreover, as we will show below, it is not difficult to implement all the necessary boundary operators for implementing domain derivatives in this framework. Thus, gradient-type iterative methods become a feasible tool to solve inverse electromagnetic scattering problems.
Using these techniques, in this paper we present boundary integral equations of the problems characterizing the domain derivatives and for the first time present an implementation and actual reconstructions of domains using these techniques in full 3D inverse electro magnetic wave scattering. Availability of Bempp has facilitated the implementation of this algorithm tremendously. Nevertheless, the mathematical foundation of the application of boundary int egral equations to electromagnetics and the analysis of corresponding boundary element methods is neither elementary nor straight-forward. For the full theory, we refer to [2][3][4]14]. We will provide a brief synopsis of those parts of the theory that are vital to the characterization of domain derivatives.
In section 2, we state the formulation of the direct problems both in variational form and as boundary integral equations. To obtain these formulations, we provide a brief description of how to define the proper Sobolev spaces and boundary differential operators for the Maxwell system. In section 3, the inverse problem is formulated and gradient-based solution schemes are discussed. A boundary integral formulation of the domain derivative of the boundary-to-far-field operator is provided. Details on how to discretize this scheme by choosing an appropriate space of boundary parameterizations are given in section 4. Section 5 discusses the necessary extensions of the Bempp library to implement the iterative solution algorithm and in section 6, we provide various numerical examples.

Electromagnetic scattering problems and boundary integral formulations
The Maxwell system, with constant wavenumber k = ω √ εµ describes the propagation of time harmonic electromagnetic waves of frequency ω in a linear isotropic homogeneous medium with electric permittivity ε and magnetic permeability µ. In this notation, the time-dependent physical electric field is given by , respectively. Nevertheless, to simplify the nomenclature, we will also refer to E as the electric and H as the magnetic field. Consider a scatterer given by a bounded Lipschitz domain D ⊆ R 3 which may represent either a perfect conducting obstacle or a penetrable inhomogeneity. We assume throughout that R 3 \ D is connected. Outside of D, µ and ε are assumed to be equal to the material constants in vacuum, ε = ε 0 , µ = µ 0 .
Consider an incident field (E i , H i ) which is assumed to be a solution of (1) in all of The presence of the scatterer gives rise to a scattered field E s , H s in R 3 \ D which is also a solution to (1) and additionally satisfies the Silver-Müller radiation condition In the case of a perfect conductor, the total field (E, where ν denotes the outward drawn normal to ∂D. In the case of a penetrable scatterer, we will assume that (ε, µ) are constant in D but different from (ε 0 , µ 0 ). In addition to (E s , H s ), the scatterer gives rise to a transmitted field inside of D which we will also denote by (E, H). The tangential traces of the physical fields outside and inside of D do not jump across the interface ∂D, which in our notation translates to Throughout this paper, we will consider solutions of these problems in the weak sense. For a Lipschitz domain Ω, we will use the Sobolev space H(curl, Ω). Let us briefly outline how the relevant trace and surface differential operators may be defined in this setting. See [2,3] and the references contained therein for a comprehensive presentation of this subject. For any connected boundary component Γ of Ω, define the tangential trace operators which can be extended continuously to γ t , γ T : are in general different from each other. They are Banach spaces with the norms Suppose that u, v ∈ H 1 (Ω, C 3 ) vanish outside of a neighborhood of Γ. As a consequence of the divergence theorem, we obtain the integral identity From (4) we can deduce that γ t : H(curl, Ω) → V T and γ T : for any u ∈ H(curl, Ω), v ∈ H 1 (Ω, C 3 ) which vanish in a neighborhood of any component of ∂Ω different from Γ. [2, proposition 3.4] states that the surface gradient ∇ Γ : Thus, the weak definition of the surface divergence as the adjoint of ∇ Γ , has a well-defined bounded extension Div Γ : V T → W . Using (5), it turns out that for u ∈ H(curl, Ω), there holds Div Γ (γ t (u)) ∈ H −1/2 (Γ) and that the mapping is bounded and surjective. In a similar way, the scalar surface curl operator Curl Γ is defined and it is proved that is also bounded and surjective. It can be shown [2, lemma 5.6] that the right-hand side of (4) extends to a duality ·, ·, Γ between H −1/2 (Div Γ , Γ) and H −1/2 (Curl Γ , Γ) such that for all u, v ∈ H(curl, Ω) which vanish in a neighborhood of any component of ∂Ω different from Γ.
We additionally define the magnetic trace operator γ N : The term magnetic comes from the fact that γ N v is the tangential trace of the magnetic field whenever v is an electric field. With these tools, we are able to state weak formulations of our scattering problems for the scatterer D. Denote by B, an open ball sufficiently large such that D ⊆ B. On the artificial boundary ∂B, we use the Calderon operator Λ : H −1/2 (Div ∂B , ∂B) → H −1/2 (Div ∂B , ∂B) that maps the tangential traces on ∂B of a radiating electric field in R 3 \ B to the tangential traces of the corresponding magnetic field, Λγ t E s = γ N E s (see [14, section 9.4]).

A weak solution of the perfect conductor problem is a field
For a penetrable scatterer, we additionally introduce the space together with the transmission conditions (3) on ∂D.
Here and in what follows, let superscripts ·+ and •-indicate traces taken from R 3 \ D and D, respectively. Then (3) translates to Solutions to these problems can be computed using boundary integral equations. Our derivation follows the presentation in [4,12], but see also [5] for a classical derivation for smooth boundaries. For the perfect conductor problem, introduce the electric potential denotes the fundamental solution of the Helmholtz equation. The corresponding electric boundary operator E k : is defined by averaging the traces of the potential from both sides of the boundary, A solution to the perfect conductor problem (7) is then given by is a solution to the electric field integral equation (EFIE), Note that the EFIE is not uniquely solvable if k 2 0 is an interior electric eigenvalue (see [7, definition 4]). Thus, from here on we will assume that this is not the case.
For the transmission problem, the magnetic potential H k : and the corresponding magnetic boundary operator by The traces of E k and H k are related to E k and H k in the following way: where I denotes the identity operator.
Using the multitrace operator from the Stratton-Chu representation formulas we obtain the equations The operators C ± k = 1 2 I ∓ A k are called Calderon projectors and map pairs of elements in H −1/2 (Div ∂D , ∂D) to admissible Cauchy data of the Maxwell system. We write (8c) as

This equation can be rewritten as
which is a boundary integral formulation of (8). Theorem 12 in [3] establishes that this equation is uniquely solvable for every incident data pair.

Inverse scattering problems and domain derivatives
So far we have discussed formulations of the direct problems of computing the scattered fields from knowledge of the incident field and the shape and physical properties of the scatterer. Now we wish to discuss the inverse problem of reconstructing the shape of the scatterer from knowledge of the incident field and of the scattered field away from the scatterer. We will approach this problem using iterative regularization methods. For this type of method, it is also important to know a priori if the scatterer is penetrable or a perfect conductor. The scattered electric field E s in R 3 \ D satisfies an asymptotic representation It is an analytic tangential vector field on the unit sphere S 2 , and E ∞ uniquely defines E s . However, the reconstruction of E s from E ∞ is a severely ill-posed problem.
We fix an incident field (E i , H i ) and restrict ourselves to an appropriately chosen class of admissible boundaries Y . We define the non-linear (electric) boundary-to-far-field operator The inverse problem can then concisely be formulated as: given an E ∞ ∈ L 2 t (S 2 ), find ∂D ∈ Y that satisfies the equation In order to formulate an iterative solution technique for this non-linear equation, we need to find its linearization. To do this in a mathematically rigorous way requires to formulate variational problems in perturbed domains generated by boundary variations. We consider variations of D, described by sufficiently small η ∈ C 1 (R 3 , R 3 ), compactly supported in a neigborhood of ∂D such that the diffeomorphism ξ defined by ξ(x) = x + η(x) gives a perturbed domain D η with an admissible boundary ∂D η = {y = ξ(x) : x ∈ ∂D}. For such variations of D there exists the so-called domain derivative E [8,9], a radiating solution of Maxwell's equations with far field E ∞ depending linearly on η, such that If we choose a certain type of parameterizations Y in a subset of a normed space X , for example star-like domains, (14) means that the operator F possesses a Frechét derivative for any admissible boundary ∂D with In order to solve (13), we use a regularized iterative Newton scheme as follows. First, we choose a starting guess ∂D 0 . Every iteration i ∈ N consists of the following steps: We stop the iteration, if the residual r falls below a chosen threshold.
Step (iii) needs some further explanation. Assuming that for small η the linearization F(∂D η ) ≈ F(∂D) + F [∂D]η is a good approximation, we consider the equation In general we cannot expect (15) to be solvable. However, applying Tikhonov regularization with some regularization parameter α > 0, transforms (15) to the uniquely solvable equation The regularization in (iii) is needed both for solvability of (15) and for damping in (iv) since we need to enforce that the updated boundary is admissible. The identity in the Tikhonov equation (16) corresponds to a penalty with respect to the norm of X . Stronger norms can be considered, which will be explained below. Consider [11] for details on such iterative regularization schemes and convergence results, which are to our knowledge not known for inverse scattering problems. The implementation of the algorithm above requires the computation of F [∂D] or, equivalently, of the domain derivative E . It is possible to characterize the domain derivative via a boundary value problem [8,9], however, it is necessary to impose additional regularity on the boundary ∂D to do so. Hence, we will from now on assume that D is of class C 1 . Note that this assumption has strong implications for the regularity of the solutions to the scattering problems (7) and (8), respectively; we obtain E, H ∈ H 1 (D), E s , H s ∈ H 1 (B \ D) [1] whenever the incident field is smooth enough.
We will briefly outline how domain derivatives can be obtained for electromagnetic scattering problems. We use the approach from [8,9] via variational methods. Different approaches exist [6,13,15]. The case treated in [9] is that of a penetrable scatterer. Introducing the normal trace operator we conclude that γ ν v ∈ H 1/2 (∂D) where v denotes any of the fields E, H, E s , H s . By [2, proposition 3.6], we also have ∇ ∂D (η) × ν ∈ H −1/2 (Div ∂D , ∂D) for any η ∈ H 1/2 (∂D), as this operator is nothing else than the vectorial surface curl in disguise.
The domain derivative E depends linearly on a variation of the geometry η. In [9] a characterization of E as the weak solution of Maxwell's equations with a certain transmission condition is derived. The jumps across ∂D depend on the normal component of η on ∂D, h = γ ν η ∈ C 1 (∂D). It is proved that E is a solution of (8a) and (8b) with E i = 0 such that the following transmission conditions are met Using the scaling matrix S from section 2, we abbreviate these conditions as Arguing as in the derivation of (12), we obtain the system of integral equations Note that (19) differs from (12) only in the right-hand side, i.e. the same solver can be used to compute E and E . For a perfectly conducting obstacle, the derivation of E follows along the same lines. The final result is that E ∈ H(curl, B \ D) is a solution of (7) for E i = 0 together with the boundary condition

The corresponding boundary integral equation is
which again is an equation with the same integral operator as in (9).

Discretization of the Newton scheme
As mentioned above, the derivation of the Newton scheme requires that the set Y of admissible boundaries is an open subset of a normed space. Thus, let Y be the set of smooth starlike domains with center in the origin. The boundaries can then be identified by positive functions on the unit sphere S 2 via spherical coordinates, i.e.
for some smooth r : S 2 → R >0 . To be more precise, we choose the open set in the normed space X = C ∞ (S 2 ) as domain for the boundary-to-far-field operator F. There are two straight-forward possibilities to discretize (16). One can discretize the full Tikhonov operator F [∂D] F [∂D] + αI or one can discretize every operator involved and multiply them on the discrete level. In general, one expects differences in the results. Here, we choose the second idea.
We start by discretizing X . Let Y m n , n ∈ N 0 , |m| n denote the normalized spherical surface harmonics, i.e. eigenfunctions of the Laplace-Beltrami operator ∆ S 2 with eigenvalue −n(n + 1) ∆ S 2 Y m n + n(n + 1)Y m n = 0, explicitly given in spherical coordinates (x, y, z) = (cos ϕ sin θ, sin ϕ sin θ, cos ϕ) ∈ R 3 by Y m n (θ, ϕ) := (2n + 1)(n − |m|)! 4π(n + |m|)! P |m| n (cos θ)e imϕ , ϕ ∈ (0, 2π), θ ∈ (0, π) with the associated Legendre functions P m n , n ∈ N 0 , m n. The functions Y m n form an orthonormal system in L 2 (S 2 , C). Since we are looking for real-valued functions, we choose as discretization of X the finite dimensional subspace X N , given by which leads to the set of admissable boundaries Y N , given by functions r ∈ X N with r > 0. We furthermore pick M ∈ N evaluation points x 1 , . . . ,x M ∈ S 2 for the far fields. Now, equation (13) reads as where α and β denote the vectors of coefficients α m n , n N , m n and β m n , n N , 1 m n. Using the linearity of the domain derivative, we can write Again, using only finitely many evaluation points we have for fixed n and m: where η k denotes the kth element of B. The product of F [∂D] F [∂D] is a complex quadratic (N + 1) 2 × (N + 1) 2 matrix, given by Note, that the adjoint of the discrete F [∂D] is just given by complex conjugation. Similarly, if we denote by E ∞ (x, ∂D) the far field with respect to ∂D, evaluated in x ∈ S 2 the right-hand side of equation (16) is an element of C (N+1) 2 , given by Instead of the the identity I, we choose the diagonal penalty matrix J, given by (J) kk = 1 + λ(k), k = 1, . . . , (N + 1) 2 . Here λ(k) := n(n + 1), such that η k = ReY m n or ImY m n .
This corresponds to an H 2 (S 2 )-penalty, since the H 2 (S 2 )-norm is equivalent to the graph norm · ∆ S 2 of the Laplace-Beltrami Operator ∆ S 2, given by and since the basis elements η ∈ B are eigenfunctions of ∆ S 2. So, solving the Tikhonov equation (16) after discretization of Y becomes solving a linear system of (N + 1) 2 equations.
is in general complex-valued. For the update in (iv), we discard the imaginary part. The full discretization of equation (23) follows by numerically evaluating F for a boundary ∂D of a given scatterer D and by evaluation of the far field of the domain derivative for a given perturbation η. We realized this with the help of the boundary element software Bempp.

Implementation remarks
The actual implementation of the Newton scheme requires the solution of many integral equations (9), (12), (19) or (20). The library Bempp (https://bempp.com/) provides a sophisticated basis for such calculations. For an introduction on how to solve electromagnetic scattering problems with Bempp, we refer to the overview paper [16] for the case of the perfect conductor and to the commented Jupyter Notebook on electromagnetic scattering from multiple dielectric spheres available on the Bempp homepage (https://bempp.com/documentation).
The application of Bempp is particularly attractive for mathematicians as the formulation of discrete problems is exactly analogous to the corresponding full mathematical problem. Thus, in this section, we will only write down equations for the full problem, in the Sobolev space H(curl, D) or its trace spaces. Once it is known how to discretize each operator, the corre sponding equations in the boundary element spaces are exactly the same.
As mentioned before, the boundary integral equations for the domain derivatives (19) and (20) differ only in the right-hand side from the boundary integral equations for the scattering problems (9) and (12), respectively. Let us recall the right-hand sides for the perfect conductor and for the penetrable scatterer This leads to the following tasks: 1. Calculation of the normal traces γ ± ν of E and H. By solving the direct scattering problem, we gain access to ϕ, the solution of on ∂D in the case of the perfect conductor and to (γ + t E s , γ + N E s ) in the case of the penetrable scatterer. We use see [2], the corresponding equation for w ∈ H(curl, B \ D) and Maxwell's equations to obtain the relations For the case of a penetrable scatterer, we arrive at explicit formulas for the normal traces. In the case of the perfect conductor, we can calculate the unknown Neumann trace γ N E with the help of the jump conditions. If ϕ is the solution of the EFIE (9), i.e. E s = −E k0 ϕ, we have Again using (26), we obtain normal traces using the surface divergence Div ∂D . One can see from (6) that the weak formulations for the surface divergence and the surface gradient are coupled. Assuming further smoothness, we can write: Up to the sign, the left-hand side can be seen as weak formulation of the surface gradient and the right-hand side as weak formulation of the surface divergence. Bempp supports the weak formulation of a number of surface differential operators. The details of the implementation are documented in a Jupyter notebook associated with this paper available through the Bempp homepage (https://bempp.com/publications). In order to calculate the product of two discrete functions f and g and represent the result in a given basis {φ j } we solve the linear system of equations to obtain the L 2 projection of the function product on the basis {φ j }. Here, depending on whether the product is scalar (in the case of γ ν hγ ν E), or vectorial (in the case of γ ν hγ T E) we choose the basis {φ j } to be either scalar or vectorial.
So far, we have explained, how we can fulfill tasks 1.-3. The remaining task is the rotation of the surface gradient ∇ ∂D × ν and the calculation of the traces γ T of E and H. Since the two traces γ t and γ T are formally connected to each other by rotation: we only need to implement the rotation. Considering we observe that ·, · ∂D can be seen as weak formulation for the operator R with Rγ T ϕ = γ t ϕ and − ·, · ∂D for the converse operator. Again, Bempp provides the tools to implement this as demonstrated in the associated Jupyter notebook available at https://bempp.com/publications. Observe that for every tangential vector field ϕ, we have (ν × (ϕ × ν)) = ϕ. Therefore, rotating the surface gradient of ϕ can be seen as applying R, i.e. (∇ ∂D ϕ) × ν = R(∇ ∂D ϕ).

Numerical results
In this section we present the results of some of the numerical experiments we have carried out.
We have successfully run reconstructions for the perfect conductor and for the penetrable scatterer. In both cases we have considered exact and noisy data. Below, we present the results for the penetrable scatterer. Results for the perfect conductor are similar, but require less computational effort.
In order to test our implementation of the reconstruction algorithm, we have picked the following shapes, see figure 1: 1. A rounded cuboid, implicitly given by (x 1 /r 1 ) n + (x 2 /r 2 ) n + (x 3 /r 3 ) n = d n with some exponent n ∈ N, positive radius d and side lengths r 1 , r 2 , r 3 > 0. 2. A peanut-shaped object, implicitly given by with R : [−1, 1] → R, R(z) = 3 5 − 2 5 cos π 2 z . We chose the first object in order to have an object close to the non-smooth cuboid. Of course, the actual rounded cuboid is smooth for every n ∈ N, but does not lie in the span of our shape basis functions. Our implementation requires star-like shaped objects, but no convexity. Therefore we picked the second object as an example for a non convex object. To cancel positive effects due to symmetry, we applied a translation such that the center of the two star-like objects does not coincide with the center of our starlike reconstruction in 0.
First, we generated the exact data E ∞ = F(∂D). We have picked 168 evaluation points x i , i = 1, . . . , 168 on the unit sphere S 2 , so that the discrete version of E ∞ is an element of C 3×168 . In order to avoid an inverse crime, we ran calculations of the exact data with meshes unrelated to those used in the reconstruction and yielding a higher accuracy. In the case of noisy data, we multiplied every component of E ∞ ∈ C 3×168 with some perturbation factor of the form where λ 1 , λ 2 are on (0, 1) uniformly distributed random numbers and the noise level δ 0. We call this noise up to δ. The effective noise level is given by Since the noisy far field E δ ∞ is no longer a (discrete) tangential vector field on the sphere, one might think of cancelling the non tangential parts of E δ ∞ , but since we apply the adjoint of F[∂D i ] on the right-hand side of (16), non tangential parts get canceled automatically. For the calculation of δ eff , we did not see any relevant difference, if we just considered the tangential part of E δ ∞ in (27). As an initial guess, we have 1}, the ball of radius one. We have observed that we have to either increase the regularization parameter α drastically or use some a priori information about the size of the scatterer for successful reconstructions.
We have chosen the regularization parameter α by experience. Using too small parameters, especially in the case of noisy data, leads to updates of the parameterization, where negative radii occur, i.e. degenerated objects. But above some critical level, we have observed robust reconstructions. Using larger than necessary α slows down the reconstruction speed, but the effect is barely noticeable. In the case of exact data, we have used α = 3 and for noisy data α = 7 for the peanut and α = 12 for the rounded cuboid, which seems to be necessary for reasonable noise levels. Reconstructions with lower values of α have failed from time to time. Choosing a fixed number of basis elements for the construction, one can calculate the L 2 (S 2 ) projection of the parameterization onto these elements. The resulting shapes are in this sense the best reconstructions, one can hope for. In figure 2, one can see the best approximation of the rounded cuboid using different numbers of basis elements.
For our reconstructions, we have chosen the material parameters ε r = 2.1, µ r = 1.0, k 0 = 1.0472, κ = 1.5175, which correspond to the scattering of Teflon (C 2 F 4 ) illuminated by VHF radiation with wavelength of 6 m. We have considered one incoming pair of plane waves (E i , H i ), given by with polarization p ∈ C 3 and direction d ∈ S 2 , satisfying p · d = 0. Again, to avoid any positive effects due to symmetry, we have chosen Let us present our reconstructions. In all cases, we have run 21 iterations without stopping rule. We have chosen N = 7, i.e. we have used (N + 1) 2 = 64 basis functions. In figure 3, the normalized residuals are plotted for the reconstruction of the peanut and the rounded cuboid with and without noise. Note the relatively large initial error e 0 with e 0 ≈ 0.4 for the rounded cuboid and e 0 ≈ 1.0 for the peanut. Also observe that after some iterations, the residuals stay at the same level. For noise free data, we have achieved final errors of e 21 ≈ 2 · 10 −3 for the peanut and e 21 ≈ 4 · 10 −3 for the rounded cuboid. Considering noisy data, we achieved e 21 = 4 · 10 −2 for the peanut and e 21 = 6.5 · 10 −2 for the rounded cuboid. In figures 4-7 reconstructions of the peanut and the rounded cuboid, each with exact and noisy data, are represented. The arrow in the picture with the exact shape D indicates the direction d of the incoming plane wave. Observe the indention of the reconstruction of the peanut along d even for noise free data, which is a known phenomena for acoustic scattering problems.
We have applied noise up to δ = 0.3, which lead to the effective noise level δ eff ≈ 0.13 for the rounded cuboid and δ eff ≈ 0.12 for the peanut.
For our calculation, we have used machines with 32 or 64 CPU cores. With the help of some parallelization, we have been able to run one iteration of our algorithm in about 10 to 20 min.
In conclusion, we have shown how iterative regularization schemes that have been used in inverse acoustic scattering problems for some time can also be implemented and applied to electromagnetic scattering problems. The results of the numerical experiments have been obtained with a reasonable computational effort and are very promising. In particular, the reconstructions appear to be of better quality than those obtained for acoustic scattering in [7]. We believe that this is due to the vectorial as opposed to scalar nature of the fields, which provides additional information for the reconstruction. Hence, Newton-type schemes may be particularly attractive for inverse electromagnetic scattering problems.