G M Vinco

S Theodoulis
email: spilios.theodoulis@isl.eu

O Sename

G Strub

Quasi-LPV Modeling of Guided Projectile Pitch Dynamics through State Transformation Technique

Keywords: Linear parameter-varying systems, Flight dynamics modeling, Model validation

In this paper, we develop a reliable Linear Parameter-Varying (LPV) model of the pitch channel dynamics of a fin-stabilized projectile. Among the available LPV design approaches, the State Transformation is analyzed, being particularly suitable for a class of systems defined as output nonlinear, and compatible with the projectile dynamics formulation. The State Transformation provides a quasi-LPV representation, which corresponds to an exact transformation of the original nonlinear model, preserving relevant couplings that are usually lost through the classical approximation methods. Some important considerations regarding the limitations of this approach are also discussed and verified in the missile dynamics. The accuracy of the obtained quasi-LPV model is assessed by means of open-loop simulations, comparing the performance to the original nonlinear model at selected flight conditions.

INTRODUCTION

In the field of aerospace applications, the Linear Parameter-Varying (LPV) framework has already provided successful results in terms of both system modeling and control design. However, it implies the nontrivial transformation process of the original nonlinear system into the corresponding LPV formulation, which requires an exhaustive investigation. A standard approach, especially for industrial settings, is based on the Jacobian linearization [START_REF] Theodoulis | LPV modeling of guided projectiles for terminal guidance[END_REF], [START_REF] Prempain | A gain scheduled autopilot design for a bank-to-turn missile[END_REF]). The LPV model is obtained as a collection of linearized LTI plants, where the investigated parameters are "frozen" around specific design values. Despite its feasibility, the main drawbacks of this method rely on the inability to fully capture the system dynamics away from the design points, potentially neglecting important nonlinearities and coupling terms that can affect the system behavior.

Several alternative approaches for LPV system derivation have been proposed in the last years with successful applications, such as function substitution [START_REF] Pfifer | Quasi-LPV model of a NDI-controlled missile based on function substitution[END_REF], and Marcos and [START_REF] Balas | Development of linear-parameter-varying models for aircraft[END_REF]), velocity-based techniques [START_REF] Leith | Gain-scheduled and nonlinear systems: dynamic analysis by velocity-based linearization families[END_REF]), and State Transformation [START_REF] Shamma | Gain-scheduled missile autopilot design using linear parameter varying transformations[END_REF], and [START_REF] Carter | Gain-scheduled bank-to-turn autopilot design using linear parameter varying transformations[END_REF]). In this paper, the latter method is investigated to develop a reliable quasi-LPV model for a projectile dynamics pitch channel autopilot. The structural configuration of the projectile, and the Bank-To-Turn flight strategy intended to be implemented, may generate relevant nonlinear couplings in the projectile dynamics. The State Transformation method provides an exact transformation between the original nonlinear system and the obtained quasi-LPV model, thus no approximations are involved in the design, increasing its capability to represent the original dynamics. A simplified version of the original nonlinear model is first developed. This allows to respect the requirements imposed by the State Transformation technique without penalizing the nonlinear coupling terms present in the dynamics. The quasi-LPV formulation derived from the simplified nonlinear model is strongly dependent on the feasible trim map and presents an internal feedback loop that updates in real-time the trimmed input at the current flight conditions. In order to avoid any stability and robustness issues, a second quasi-LPV model is generated by augmenting the system with an integrator at the input. The performance of the two quasi-LPV models is analyzed and compared with both the original nonlinear dynamics and the simplified one to assess the impact of the assumptions on the model reliability.

The investigated nonlinear model dynamics is discussed in Section 2, together with some relevant insight related to the projectile aerodynamic coefficients. Section 3 addresses more in details the structure of the output nonlinear (or output dependent) class of systems, as well as the State Transformation approach employed for the derivation of the quasi-LPV model. Important considerations are also discussed, with the intention of reformulating the nonlinear model in a suitable form that respects the State Transformation assumptions. Finally, in Section 4, openloop simulations are performed, aiming to compare the time response of the nonlinear and quasi-LPV models and to assess the accuracy of the latter.

PROJECTILE NONLINEAR MODEL

The investigated concept consists of a 155 mm finstabilized guided projectile, equipped with a set of four symmetrical rear fins, and two horizontal front canards, in a non-coplanar configuration, as shown in Fig. 1.

Canards Fins

Fig. 1. Projectile concept: details on the control surfaces.

The complete nonlinear reference model has been derived in [START_REF] Vinco | Flight dynamics modeling and simulator design for a new class of long-range guided projectiles[END_REF]. The original equations of motion describing the projectile translational dynamics in the body coordinates (B) have been converted to the more convenient stability coordinates, in terms of the angleof-attack, the sideslip-angle, and the true airspeed (zero relative wind), respectively α, β, and V , as in Zipfel (2014).

Pitch Dynamics

The proposed study focuses on the pitch channel of the projectile dynamics, intended to be analyzed in view of the future autopilot design. The projectile behavior is modeled through the dynamics of α, and the pitch rate, q:

α = - X sin α mV cos β + Z cos α mV cos β + g V cos β
(sin α sin θ + cos α cos θ cos ϕ)

+ q -p tan β cos α -r tan β sin α, q = 1 I yy [M -pr (I xx -I zz )] , (1)
where m stands for the mass of the projectile, g is the acceleration of gravity, I xx , I yy , and I zz refer to the moments of inertia about the x, y, and z axes, respectively. Additionally, p and r, correspond to the roll and yaw rates, respectively, while ϕ and θ indicate the roll and pitch angles. In particular, the air mass is assumed at rest, so no relative wind contributions are modeled and, due to the projectile planar symmetry with respect to the x -z and x -y planes, the products of Inertia I xy , I yz , and I xz result to be negligible.

Aerodynamic Model

Concerning the projectile aerodynamics, the Simple Linear aerodynamic model discussed in [START_REF] Vinco | Flight dynamics modeling and simulator design for a new class of long-range guided projectiles[END_REF] is here employed, in reason of the observed high fidelity level and the lower mathematical complexity. Furthermore, the influence of β on the pitch dynamics is expected to be negligible. The longitudinal and vertical forces, X and Z, and the pitching moment, M , are functions of the dynamic pressure q = 1 2 ρ (h) V 2 , the reference surface S and caliber d, the air density ρ, and the altitude h, as follows:

X = qS C XS + C X δ eff , Z = qS C ZS + d V C ZD q + C Z δq , M = qSd C mS + d V C mD q + C m δq .
(2)

The projectile static coefficients, C XS (M, α), C ZS (M, α), and C mS (M, α) were derived in the Simple Linear model as a function of the Mach number, M, and the angle α.

Similarly, the virtual control coefficients C X δ eff (M, δ eff ), C Z δq (M, δ q ), C m δq (M, δ q ) are expressed with respect to the virtual pitch δ q , and virtual roll δ p control deflections, their nonlinear combination δ eff = δ 2 p + δ 2 q , and the Mach number.

C X δ eff = C X δ 0 (M) + C X δ 2 (M) sin 2 δ eff , C Z δq = C Z δ 1 (M) sin δ q + C Z δ 3 (M) sin 3 δ q , C m δq = C m δ 1 (M) sin δ q + C m δ 3 (M) sin 3 δ q .
(3)

In particular, the virtual deflections δ p and δ q , are generated by the mutual interactions between the individual right and left canard deflections, δ r and δ l , respectively. As shown in Figs. 2a -2b, the virtual roll is defined by the differential deflection δ p = δ l -δr 2 , while the virtual pitch by the concurrent deflection δ q = δ l +δr 2 , respectively.

Right

Left

C Z δr (δr) Coherently, the virtual control coefficients were modeled as an approximation of the global coefficients:

C Z δ l (δ l ) C l δp (δp) B (a) Right Left C Z δr (δr) C Z δ l (δ l ) Cm δq (δq) B (b)
C Z δ = C Z δr (M, δ r ) + C Z δ l (M, δ l ), C m δq = C m δr (M, δ r ) + C Z δ l (M, δ l ), (4)
generated by the linear superposition of the individual canard contributions. The derivation is extensively discussed in [START_REF] Vinco | Flight dynamics modeling and simulator design for a new class of long-range guided projectiles[END_REF]. An interpolation analysis was conducted to define the deflection range in which the canards response results to be linear, motivating the superposition of the individual deflections in a combined contribution.

The results exhibited a reliable level of accuracy for the deflection range δ r , δ l ∈ [-25

• , 25 • ].
Finally, the aerodynamic damping coefficients C ZD (M, α) and C mD (M, α) are also addressed in the model, aiming to provide with an exhaustive characterization of the projectile behavior.

QUASI-LPV MODELING

Linear Parameter-Varying refers to a class of systems whose state space representation is defined as a continuous function of a time-varying vector of scheduling parameters ρ 1 , assumed to be measurable in real-time. For a detailed definition, the reader can refer to [START_REF] Sename | Robust control and linear parameter varying approaches: application to vehicle dynamics[END_REF] and [START_REF] Mohammadpour | Control of linear parameter varying systems with applications[END_REF]. When the state of the system, x, is decomposed into a scheduling subset, z(t), and a nonscheduling subset w(t), the system is specifically indicated as quasi-LPV. Coherently, the scheduling parameters are divided into an exogenous subset, Ω(t), and an endogenous subset, z(t).

x(t) = [z(t) w(t)] T , ρ(t) = [z(t) Ω(t)] T .

State Transformation

The approach was first introduced in [START_REF] Shamma | Gain-scheduled missile autopilot design using linear parameter varying transformations[END_REF] and later further elaborated in [START_REF] Leith | Survey of gainscheduling analysis and design[END_REF]. It consists of a state transformation aimed to remove all the nonlinearities present in the model, that do not depend on the scheduling variables. As extensively discussed in [START_REF] Leith | Survey of gainscheduling analysis and design[END_REF], the State Transformation is suitable for a restricted class of systems, defined as Output Nonlinear or Output Dependent, and expressed in the generic form:

ż ẇ = f 1 (ρ) f 2 (ρ) + A 11 (ρ) A 12 (ρ) A 21 (ρ) A 22 (ρ) z w + B 1 (ρ) B 2 (ρ) u, y = z, (5)
where u(t) ∈ R nu represents the control input, and the model nonlinearities f 1 (ρ) and f 2 (ρ), appear only as a function of the measured output z(t) ∈ R nz , and of a selected set of exogenous parameters Ω(t).

The transformation is achieved by means of trimming functions of the nonscheduling states w(t), and of the control input, assuming that n z = n u . By considering a feasible trim region of the flight envelope, where the trimming functions w eq (ρ) and u eq (ρ), respectively, are continuously differentiable for every z and ρ, such that:

0 0 = f 1 (ρ) f 2 (ρ) + A 11 (ρ) A 12 (ρ) A 21 (ρ) A 22 (ρ) z w eq (ρ) + B 1 (ρ) B 2 (ρ) u eq (ρ), (6)
and by defining the new set of transformed entries:

ξ(t) := w(t) -w eq (ρ(t)), ν(t) := u(t) -u eq (ρ(t)), Ã22 (ρ) := A 22 (ρ) - ∂w eq ∂z A 12 (ρ), B2 (ρ) := B 2 (ρ) - ∂w eq ∂z B 1 (ρ), E(ρ) := - ∂w eq ∂Ω , (7)
1 The scheduling parameters set ρ used in Sections 3 -4 is different from the air density, ρ(h), defined in Section 2.

the original nonlinear system in (5) can be reformulated as a quasi-LPV model, in the form:

ż ξ = 0 A 12 (ρ) 0 Ã22 (ρ) z ξ + B 1 (ρ) B2 (ρ) ν + 0 E(ρ) Ω. ( 8
)
In particular, the additional input term Ω describes the impact of the dynamics of the exogenous variables on the system, through the trimming functions. For the current analysis, this term can be assumed as a disturbance to be rejected and thus neglected in the model, as discussed in [START_REF] Balas | Linear, parameter-varying control and its application to aerospace systems[END_REF]. However, it represents an interesting aspect to be further investigated for future works.

The dependence of this formulation on the feasible trim map can be inferred from ( 7). Indeed, the definition of the deviated input ν depends on the selected equilibrium condition u eq . From a stability and robustness perspective, this internal feedback loop may affect the performance of the system and the following control design. Thus, a straightforward solution would be to require that the trimmed input u eq (ρ) = 0, imposing a strong limitation to the region of the trim map that respects this condition.

An interesting design alternative consists of adding an integrator at the plant input. By defining, u = σ, we can reformulate the system (8), as:

  ż ξ ν  =   0 A 12 (ρ) B 1 (ρ) 0 Ã22 (ρ) B2 (ρ) 0 Ã32 (ρ) B3 (ρ)   z ξ ν + 0 0 I σ, (9) with:
Ã32 (ρ) := -∂u eq ∂z A 12 (ρ); B3 (ρ) := -∂u eq ∂z B 1 (ρ). ( 10)

In this way, the new input is uniformly zero at every equilibrium point and the feedback loop does not affect the system anymore. Furthermore, this solution is motivated by the assumption that the controller intended to be designed contains pure integral action for steady state tracking error elimination. Thus, the controller integral action can be formally included in the system description.

The main advantage of the State Transformation approach relies on the exact transformation between the original nonlinear system and the obtained quasi-LPV model, which avoids any form of approximation, as opposed to the classical Jacobian linearization. Additionally, due to the integration of the control input, the reliability of the quasi-LPV model does not depend on the specific operating point, and the term ẇeq = ∂weq ∂z ż, included in (7), allows to account for any off-equilibrium conditions.

Pitch Channel quasi-LPV Model

As previously mentioned, the State Transformation approach results to be very suitable for aerospace application. It is easy to observe how the pitch channel dynamics, presented in (1), is already formulated as the Output Nonlinear system in (5). In particular, the state variables [z, w] T = [α, q] T represent the measured and the unmeasured output, respectively, while the set of varying parameters intended to be investigated corresponds to [z, Ω] T = [α, (V, h)] T , being α and [V, h] the endogenous and the exogenous variables, respectively.

The reader can notice in the complete projectile nonlinear model in [START_REF] Vinco | Flight dynamics modeling and simulator design for a new class of long-range guided projectiles[END_REF], that the true airspeed V was defined through its components along the main body coordinates, [u B , v B , w B ] T , included in the model state. Thus, V itself could be considered as a state for the pitch dynamics, describing the long period longitudinal dynamics (phugoid) of the projectile. However, because of the lack of any relevant control authority on the projectile longitudinal velocity, u B , V is addressed in this paper as an exogenous scheduling variable. As implied by the Output Nonlinear formulation, all the nonlinearities involved in the model are expressed as a function of the measured output. The only exceptions rely on the control coefficient definitions in (3). Indeed, (C X δ eff , C Z δq , C m δq ) are modeled as nonlinear functions of the virtual control inputs δ q and δ p . Consequently, the system is not affine with respect to the control inputs.

In order to adjust the projectile dynamics into a suitable form, a correction of the control terms is required. To avoid any approximation through a standard Taylor expansion, the aerodynamics data analyzed inVinco et al. ( 2022 By restricting the analysis to the deflections range (normalized) δ r , δ l ∈ [-0.7 • , 0.7 • ], and by employing the first order regression models C Z δr = C Z δ 1 δ r , and C m δr = C m δ 1 δ r , the obtained accuracy results higher than 95% in terms of Coefficient of Determination, for both the coefficients. The same results are observed for the left canard deflection. As discussed in [START_REF] Vinco | Flight dynamics modeling and simulator design for a new class of long-range guided projectiles[END_REF], a second analysis is then performed in order to model the resulting global coefficients as in (4), into the corresponding virtual formulations (C Z δq ,C m δq ).

The interpolation error is expressed in terms of the difference between the global coefficients and the virtual ones, normalized by the value of the original aerodynamic data. As shown in Fig. 4, both the vertical force and the pitching moment exhibit an average error lower than 10% in the investigated deflection range. Thus, the control coefficients can be reformulated in a way to be affine with respect to the pitch dynamics. The longitudinal control force, C X δ eff , nonlinear function of δ p and δ q , is neglected for the purposes of this analysis. The choice is motivated by the intent to operate in a gliding flight condition, a subsonic regime, where the drag generated by the canard deflections does not provide a relevant impact on the projectile aerodynamics.

MODEL VALIDATION

In order to assess the reliability of the obtained quasi-LPV model, a campaign of open-loop simulations has been performed. The scheduling parameters investigated during the analysis are ρ = [α, V ], while the altitude h is assumed as constant. In particular, for each combination of the parameters, the corresponding trimming functions w eq (ρ) = q eq (ρ) and u eq (ρ) = δ q,eq (ρ) are evaluated to define the initial equilibrium condition of the simulation. The analysis shown in Figs. 5a -5b provide a map of the equilibrium points variation as a function of the scheduling parameter. (2022), and the lower equilibrium deflection commands observed in Fig. 5b, the equilibrium conditions considered for the following performance comparison are shown in Tab. 1. Indeed, the low equilibrium deflection, δ q,eq , allows the employment of a larger input command as a form of perturbation to the system, without exceeding the validity range of the aerodynamic model. In the following sections, both the quasi-LPV model formulated as in ( 8), with the internal feedback loop affecting the control input, and the Augmented quasi-LPV model as in ( 9) assuming integral action on the input, are investigated. In particular, their open-loop responses are compared to the ones of the original complete nonlinear model (NL), including the aerodynamics in (3), and the simplified one (NL sim ), presenting the linear aerodynamic model discussed in Section 3.

Quasi-LPV Model

Nonlinear Dynamics α Aerodynamics quasi-LPV Output Deviation q dev α δ q,dev q α δ q δ q α Trimming δ q,eq -+ q dev α δ q,dev δ q,dev Fig. 6. Quasi-LPV model design scheme.

The simulation is performed assuming the initial equilibrium conditions in Tab. 1, and a set of input commands corresponding to a first deflection δ q1 = 20 • at t 1 = 5 s, and second deflection of δ q2 = -15 • at t 2 = 30 s. The results in Fig. 7a -7c show the open-loop time responses related to the transformed states α(t) and q dev (t) = q(t) -q eq (t), and to the transformed input δ q,dev (t) = δ q (t) -δ q,eq (t), respectively. In particular, the q dev (t) and the δ q,dev (t) of the nonlinear models, are obtained in post processing by applying the transformation to the original state q(t), and input δ q (t), of the system, as shown in Fig. 6.

The curves describing the quasi-LPV (q-LPV) and the simplified nonlinear model (NL sim ) are perfectly overlapped since the former has been obtained through an exact transformation of the latter. On the other side, a small mismatch in the convergence values is observed The large oscillations characterizing all the three responses can be associated with the stability of the projectile dynamics since, at the flight conditions in Tab. 1, the system is near the limit of static stability.

Augmented quasi-LPV Model

The same simulation is performed on the Augmented quasi-LPV model. The main difference consists of the input commands (δ q1 , δ q2 ) that are injected into the systems since the input now are integrated for both the quasi-LPV and the nonlinear models2 . In this context, two pulses simulating the canards deflection rates ( δq1 , δq2 ) are employed as perturbations. The deflection rates are selected such that, once integrated, they provide the same control deflection (δ q1 , δ q2 ), previously defined. The curves in Figs. 9a -9c, related to the Augmented q-LPV and the NL sim models, present a mild difference due to the way the integrator is implemented. Indeed, in the former system, the deflection is accounted for as a state, thus the integrator is part of the internal dynamics. Differently, in the two nonlinear the action is applied to the input before being into the system. The advantage of the Augmented formulation consists in the removal of the internal loop updating the trimming values. As previously discussed, it can be critical in terms of performance since any forms of unmodeled dynamics or approximation are injected as input disturbances. Additionally, the RMSE values in Tab. 3 show a better accuracy between the Augmented q-LPV model and the NL one with respect to the results in Tab. 2, confirming the improvement of the augmented formulation. The differences between the NL and the NL sim models, generated by the aerodynamics approximations, appear to be reduced by the integration process. The only exception relies on deviated deflection, δ q,dev , where the integration tends to increase the bias in the convergence values. Finally, it is interesting to observe how the convergence values for the three states are almost the same as the ones in Figs. 7a -7c, but the oscillation amplitudes are way narrower.

Fig. 2 .

 2 Fig. 2. Individual and virtual control deflections: (a) virtual roll, (b) virtual pitch.

Fig. 3 .

 3 Fig. 3. Canards response for deflection δ i , i = r, l: (a) vertical force coefficient, (b) pitching moment coefficient. (All the data have been normalized in reason of confidentiality.)

Fig. 4 .

 4 Fig. 4. Interpolation error: (a) -(b) vertical force coefficient at different Mach values; (c) -(d) pitching moment coefficient at different Mach values. (All the data have been normalized in reason of confidentiality.)

Fig. 5 .

 5 Fig. 5. Trimming functions analysis: (a) pitch rate equilibrium values; (b) control deflection equilibrium values.

Fig. 7 .

 7 Fig. 7. Simulation time responses: (a) angle-of-attack; (b) deviated pitch rate; (c) deviated input deflections.with respect to the full nonlinear model (NL), because of the different aerodynamics. For the same reason, the NL model shows some initial oscillations around the equilibrium conditions, which have been evaluated on the base of the NL sim model. The results in Tab. 2 present the Root Mean Square Error (RMSE) evaluated between the three models for each of the state variables, and normalized by the maximum variation range of each state variable related to the NL model. In particular, the RMSE values obtained between the NL sim and the q-LPV models are very small, while the others are almost the same, coherently to what previously observed in the simulation results.

Fig. 8 .

 8 Fig. 8. Augmented quasi-LPV model design scheme.

Fig. 9 .

 9 Fig. 9. Simulation time responses: (a) angle-of-attack; (b) deviated pitch rate; (c) deviated input deflections.

Table 1 .

 1 Trim point and equilibrium conditions.

	α(deg) V (m/s) h(m) qeq(deg/s) δq,eq(deg)
	14	158.026	6029	-1	0.6

Table 2 .

 2 RMSE results: quasi-LPV model.

		NL -qLPV NL sim -qLPV NL -NL sim
	α	0.1302	0.0018	0.1310
	q dev	0.1332	0.0020	0.1342
	δ q,dev	0.1352	0.0019	0.1361

Table

  RMSE results: Augmented quasi-LPV model.

		NL -qLPV NL sim -qLPV NL -NL sim
	α	0.0968	0.0373	0.0913
	q dev	0.0341	0.0585	0.0805
	δ q,dev	0.1705	0.0655	0.1886

An integrator is added at the nonlinear system input to be consistent with the q-LPV formulation during the comparison.5. CONCLUSIONIn this paper, the pitch dynamics of a guided projectile in investigated. The nonlinear model is converted through the State Transformation approach into two quasi-LPV models. The first model, (quasi-LPV), presents a strong dependence on the trim map due to an internal loop updating the system input. In the second one, (Augmented quasi-LPV), this issue is addressed by augmenting the system with an integrator at the input. Open-loop simulation results and RMSE analysis show higher accuracy in the case of the augmented system. This model will be employed for the development of a LPV-based pitch autopilot.


Olivier Pourret

New Rare Earth Elements Found In Amazon Estuary? Science Trends

New rare earth elements resource large of Amazonia?

Olivier Pourret

The industry is nowadays more and more demanding of many metals, among these many metals, rare earths and more particularly some of them. These elements (e.g., neodymium, Nd) are used in particular to produce powerful permanent magnets that are used in wind turbines, electric motors ... (Figure 1). Others such as europium are used as luminophors. These rare earths are thus part of the materials considered strategic for the economy. The rare earths are in the strict sense the 14 lanthanides, elements of the periodic table (Figure 2) whose atomic numbers go from 57 to 71. The yttrium of atomic number 39 is often added there because it has chemical properties close to these. Their atoms have an external electronic configuration that includes 14 orbitals on the layer 4f, at the origin of this family of 14 lanthanides. Their particular electromagnetic properties come from the fact that the 4f orbitals do not properly mask the nucleus field. Rare earths have similar chemical properties because of their three valence electrons.

Figure 2 Periodic table of the elements

Contrary to what their name suggests (rare earths), they are not so rare. However, they are strategic materials for the global economy (i.e. critical metals; Figure 3).

Figure 3 Criticity of metals

It is important to distinguish between light rare earth elements (LREE) and rare earths (HREE); defines these two groups as follows: -LREE: La, Ce, Pr, Nd, Sm; -HREE: Eu, Gd, Tb, Dy, Ho, Er, Tm, Yb, Lu (plus Y).

There is no single rare earth primary deposit, but mineral deposits containing several rare earths, predominantly either light rare earths or heavy rare earths. However, deposits rich in HREE are much rarer than deposits rich in LREE. For the moment, the main exploited deposits of HREE are located in the south of China, which leads to a monopoly situation (in 2010, respectively 99% and 87% of the supply of HREE and LREE came from China) and creates risks for supply, which are exacerbated by Chinese export quotas that emerged in 2011.

Moreover, in a 100% electric vehicle world, the commodity demand will increase by ~655% (Figure 4, source UBS). In this context, several researchers around the world are working on discovering new potential resources. [START_REF] Pourret | Continental shelves as potential resource of rare earth elements[END_REF] recently proposed to explore continental shelves, especially in Atlantic Ocean front of Amazon estuary (Figure 5).

Figure 1

 1 Figure 1 Rare earth elements used in a car (source Toyota)

Figure 4

 4 Figure 4 Incremental commodity demand in a 100% electric vehicle world (proportion of

Figure 5

 5 Figure 5 Sketch illustrating the REE external cycle and summarizing the processes responsible for REE


When rock meets life

Olivier Pourret All started in 2005, during a trip to the United States, where I went to the famous conference of Geochemistry (Goldschmidt conference, Moscow, Idaho). The following weeks I crossed the Rocky Mountains and visited many national parks including Yellowstone Park. Of the many geological, geothermal, and volcanic curiosities in Yellowstone National Park, and outside of wildlife (including some incredible photographs of bison), I have been fascinated by hot springs and associated termophile organisms as this photograph of the "Morning Glory Pool" north of the "Old Faithful Geyser" (Figure 1). The orange and yellow bacteria thrive in the cooler waters by secreting siderophores (or metallophores) favoring the precipitation of iron (Fe) and manganese (Mn) oxides, while for higher temperatures we find blue or green bacteria. Bacteria grow massively when the temperature of hot springs drops. Siderophores or metallophores are chelating molecules of iron produced by microorganisms when intracellular iron concentrations are low, which allows them to "draw" in its environment the iron necessary for its development. These molecules of low molecular weight have a very strong affinity for trivalent iron Fe(III) and are capable of forming [siderophore-Fe(III)] complexes that will internalize the iron necessary for the functioning of the cell.

Generally, siderophores are classified by their functional groups as catecholates, hydroxamates and hydroxycarboxylates. Other structural chemical and functional group modifications are also observed. The functional groups participate in the iron chelation process when the [siderophore-Fe(III)] complex is formed. Acidophiles, alkaliphiles, halophiles, thermophiles, psychrophiles, piezophiles, extremophiles have particular iron requirements depending on the environmental conditions in which they develop. Most of the research that is concerned with the production of siderophores by extremophiles is based on siderophore concentration and/or genomic studies determining the presence of siderophore synthesis and transport genes. The siderophores produced by extremophiles are poorly known and much remains to be done to elucidate the chemical structures and their role in the survival of microorganisms and the cycles of metals in extreme environments.

In this context, we realized with colleagues from Jacobs University in Bremen, Germany (under the direction of Michael Bau) rare earth elements (REE) sorption and modeling experiments on manganese oxides in the presence of a biogenic compound of desferrioxamine B type (DFOB). This compound belongs to the group of metal-specific organic ligands commonly called siderophores that I just briefly described. These low molecular weight organic molecules are synthesized by bacteria and plants in oxidized environments to cope with the iron deficiency of the medium. In oxidized systems, iron is present in its trivalent form and is therefore immediately fixed in insoluble oxyhydroxide structures and thus removed from the aqueous solution.

Therefore, under modern atmospheric conditions, the essential nutrient Fe is not available to organisms and siderophores are excreted by plants, microbes and fungi to facilitate the dissolution of iron oxides. However, it has recently been shown that siderophores not only form very stable complexes with Fe (III), but also significantly complex other polyvalent metal ions, some of which are considered toxic at higher concentrations, but are essential nutrients at trace concentrations. Several studies have shown that siderophore chelation of metals reduces metal uptake by bacteria, fungi and plants, and siderophores play a vital role in the detoxification of metals.

In the absence of siderophores, trivalent Ce is surface-oxidized to Ce(IV) during sorption onto Mn (hydr)oxides. Such Mn precipitates, often show positive Ce anomalies, whereas the ambient solutions display negative ones. However, REE scavenging in the presence of DFOB produces negative Ce anomalies in the Mn (hydr)oxides and a distinct and characteristic positive Ce anomaly in the siderophore-bearing solution. Furthermore, the heavy REE (Tb to Lu) also remain in solution, creating strongly Ce and HREE-enriched residual solutions during scavenging of REE onto Mn (hydr)oxides. This process efficiently separates light REE (La to Gd excluding Ce) from HREE. Such preferential partitioning of Ce into the residual aqueous phase has been described in the literature for natural Mn (hydr)oxides, but the mechanism for this is still unclear. Our experimental results demonstrate that biogenic organic ligands such as siderophores that form very strong aqueous complexes with Ce 4+ , may produce solutions with positive Ce anomalies and Mn (hydr)oxide minerals with negative Ce anomalies, complicating the use of Ce anomalies as paleoredox proxies.

Our work on the "Development of negative cerium anomaly on manganese hydroxide precipitates in the presence of siderophores" has just been published in the January 2017 issue of the journal Geochimica and Cosmochimica Acta and can be consulted here: https://www.sciencedirect.com/science/article/pii/S0016703716305336 Kraemer, D.; Tepe, N.; Pourret, O.; Bau, M., Negative cerium anomalies in manganese (hydr)oxide precipitates due to cerium oxidation in the presence of dissolved siderophores. Geochimica et Cosmochimica Acta 2017, 196, 197-208.

Figure 1 .

 1 Figure 1. Photograph of the Morning Glory Pool at Yelloswtone (May 2005)


Christian Daveau
email: christian.daveau@u-cergy.fr

Abdessatar Khelifi
email: ab-dessatar.khelifi@fsb.rnu.tn

Soumaya Oueslati
email: ouslati@yahoo.fr

Small pertubations of an interface for Stokes system

Keywords: small perturbations, interface problem, Stokes system, asymptotic expansions, boundary integral method, viscous moment tensors AMS Subject Classification: 35C20, 31B10

We consider the Stokes system for a viscous medium consisting of an impurely (inclusion) merged in a consistent background medium. Appointed on both field expansion and layer potential methods, we strictly derive the asymptotic expansion of the perturbed velocity field because of small perturbations in the interface of the inclusion. We use these techniques to determine a relationship between Stokes solutions measurements and the shape of the object. Moreover, we may prove an asymptotic expansion for the perturbation in the viscosity moment tensors (VMTs) caused by the presence of small changes in the interface of the impurety.

Introduction

In this paper, we suppose that any shape deformation like [START_REF] Ammari | Conductivity interface problems. Part I: small perturbations of an interface[END_REF] occurs inside a bounded domain Ω in a Hölder space. Based on the works [START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF][START_REF] Ammari | Reconstruction of small interface changes of an inclusion from modal measurements II: The elastic case[END_REF][START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF][START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF][START_REF] Boyer | Fabrie Mathematical tools for the study of the incompressible Navier-Stokes equations and related models[END_REF][START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF][START_REF] Boujemaa | Small perturbation of a surface: full Maxwells equations[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF], we may advance the asymptotic expansions for the perturbed solutions of the Stokes system and we may give explicitly both the first term correction and a relationship between Stokes solutions measurements and the shape of the object. Since the Stokes system can be viewed as the incompressible limit (λ, λ → +∞) of the Lamé system [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF][START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF], it is natural to view the viscosity moment tensor VMT as a limit of the elastic moment tensors EMT. We will check formally, in Section 3, this limit when the compressional modulus goes to infinity. Our method in this investigation may be deeply based on the ones developed in [START_REF] Ammari | Reconstruction of small interface changes of an inclusion from modal measurements II: The elastic case[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF] for the elastic case.

Suppose in this section that Ω ⊂ R d (d = 2, 3) is a bounded domain with C 1,1 -boundary. Let D be an open subset of Ω such that dist(∂Ω, D) ≥ d 0 > 0, [START_REF] Ammari | Reconstruction of small interface changes of an inclusion from modal measurements II: The elastic case[END_REF] representing an inclusion made of a different Newtonian fluid material. Assume that the boundary ∂D is of class C 2,1 . We denote by σ 0 and σ 1 the stress tensor fields in Ω\D and D, respectively. We assume that both Ω\D and D are occupied by isotropic and homogeneous Newtonian fluids. Then, the tensors σ 0 and σ 1 may be given by (σ s ) ijlk = µ s (δ ki δ lj + δ kj δ li ) for i, j, k, l = 1, 2 and for s = 0, 1,

where µ 0 and µ 1 are the viscosity constants of the flow in Ω\D and D, respectively. Given two (2 × 2) matrices A and B we denote by A : B the contraction, i.e., A : B = ij a ij b ij . Now, it is useful to introduce the strain rate tensor D for the flow as follows:

D(v) := 1 2 (∇v + (∇v) T ) = 1 2 ( ∂v i ∂x j + ∂v j ∂x i ) 1≤i,j≤d
.

Let σ D := σ 0 χ Ω\D + σ 1 χ D where σ 0 and σ 1 are given by [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF]. Then according to [START_REF] Boyer | Fabrie Mathematical tools for the study of the incompressible Navier-Stokes equations and related models[END_REF], one can rewrite the Stokes problem as the generalized one:

   -div(σD(v) -p Id) = 0 in Ω ∇.v = 0 in Ω v = F on ∂Ω (3)
where v = (v i ) i=1,••• ,d denotes the velocity field while the scalar function p is the pressure. F is a given vector valued function. Moreover, (v, p) ∈ (H 2 (Ω)) d × H 1 (Ω) and Id means the identity. Here we refer the reader to [START_REF] Boyer | Fabrie Mathematical tools for the study of the incompressible Navier-Stokes equations and related models[END_REF]Theorem IV.5.8] for more details about the regularity properties. On the other hand, if S mens the Stokes operator, it is well know that (see for example, [START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF][START_REF] Boyer | Fabrie Mathematical tools for the study of the incompressible Navier-Stokes equations and related models[END_REF]) there exist an orthogonal projection (the Leray projection) P such that:

S(v) = P (-∆v), ∀v ∈ (H 2 (Ω)) 2 ∩ H(Ω).
In consequence, the Stokes operator enters the general framework about spectral properties of the operator -∆. So, througout this paper we suppose that 0 is not an eigenvalue of S.

The conormal derivative ∂v ∂ν on ∂D is defined to be

∂v ∂ν = µ 0 ∂v ∂n -p n ( 5
)
where n is the outward unit normal to ∂D

The problem (3) is equivalent to the following transmission problem:

               -div(σ 0 D(v) -p Id) = 0 in Ω \ D, -div(σ 1 D(v) -p Id) = 0 in D, v | -= v | + on ∂D, ∂v ∂ν | -= ∂v ∂ν | + on ∂D, ∇.v = 0 in Ω, v = F on ∂Ω, (6)
where F is a known vector-valued function. The quantities v | ± on ∂D denote the limits from outside and inside of D, respectively. We will also sometimes use v e for v | + and v i for v | -.

The δ-perturbation, denoted by D δ , of the domain D is given by

∂D δ := {x = x + δh(x)ν(x), x ∈ ∂D} (7)
where here h(x) is assumed to be a real function in C 1,1 (∂D) that satisfies

h C 1,1 (∂D) < c 0 for some positive constant c 0 . (8)
Let σ D δ := σ 0 χ Ω\D δ + σ 1 χ D δ . Therefore, for the same given vector-valued function F , problem (6) can be generalized to the following perturbed Stokes problem:

     -div(σ D δ D(v δ ) -p δ Id) = 0 in Ω ∇.v δ = 0 in Ω v δ = F in ∂Ω (9)
where

(v δ , p δ ) ∈ (H 2 (Ω)) d × H 1 (Ω).
The transmission stokes problem associated to ( 9) is given by:

               -div(σ 0 D(v δ ) -p δ Id) = 0 in Ω \ D δ , -div(σ 1 D(v δ ) -p δ Id) = 0 in D δ , v δ | -= v δ | + on ∂D δ , ∂v δ ∂ν | -= ∂v δ ∂ν | + on ∂D δ , ∇.v δ = 0 in Ω, v δ = F on ∂Ω. (10)
The first main result of this paper is the derivation of the leading-order term in the asymptotic expansion of (v δ -v) Ω as δ → 0, where Ω is a bounded region outside the inclusion D and away from ∂D.

The methods and results developed in this paper can be generalized to higher dimension interface problems and can be extended to other PDEs systems, such as Lamé, Maxwell and also for cases with high oscillating or periodic boundary.

Main Results

Now we are ready to state the main achievements of this paper. The first main result of this paper is the following derivation of the leading-order term in the asymptotic expansion of v δ -v 0 as δ → 0, where Ω is a bounded region outside the inclusion D and away from both ∂D and ∂Ω.

Theorem 2.1. Assume that all hypothesis (1), ( 4), [START_REF] Ammari | Conductivity interface problems. Part I: small perturbations of an interface[END_REF] and (8) are satisfied. Let v and v δ be the solutions to [START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF] and [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF], respectively. Let Ω be a bounded region outside the inclusion D and away from ∂D. For x ∈ Ω, the following pointwise asymptotic expansion holds:

v δ (x) = v(x) + δv 1 (x) + O(δ 2 ), (11)
p δ (x) = p(x) + δp 1 (x) + O(δ 2 ). ( 12
)
The remainder O(δ 2 ) depends only on µ 0 , µ 1 , Ω, D and v 1 is the unique solution to

               -div(σ 0 D(v 1 ) -p 1 Id) = 0 in Ω \ D, -div(σ 1 D(v 1 ) -p 1 Id) = 0 in D, v 1 | --v 1 | + = 0 on ∂D, ∂v1 ∂ν | --∂v1 ∂ν | + = ∂ ∂τ (h([σ 1 -M 0,1 ]D(v i ))τ ) on ∂D, ∇.v 1 = 0 in Ω, v 1 = F on ∂Ω (13)
with τ as the tangential vector to ∂D, and the viscosity moment tensor

M l,k = 2µ k I + 2(µ l -µ k )I ⊗ (τ ⊗ τ ). ( 14
)
Here I means the identity four-tensor, µ l , µ k given by (2), I the identity in R d and τ a tangential vector.

We will discuss the viscosity moment tensor M = (M l,k ) in the next section. As a consequence of the results of Theorem 2.1, we obtain the following relationship between velocity measurements and the deformation h. The scalar product in R d will be denoted by . and sometimes to ease the notations, by < ., . >. Theorem 2.2. Assume that all hypothesis (1), ( 4), [START_REF] Ammari | Conductivity interface problems. Part I: small perturbations of an interface[END_REF] and (8) are satisfied. Let M 0,1 be given by [START_REF] Depauw | Solutions des équations de Navier Stokes incompressibles dans un domaine exterieur[END_REF] and S be a Lipschitz closed curve enclosing D and away from ∂D Let (v, p) and (v δ , p δ ) be the solutions to [START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF] and [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF], respectively, and (u, q) be the solution of the following system:

               -div(σ 0 D(u) -q Id) = 0 in Ω \ D, -div(σ 1 D(u) -q Id) = 0 in D, u | --u | + = 0 on ∂D, ∂u ∂ν | --∂u ∂ν | + = 0 on ∂D, ∇.u = 0 in Ω, u(x) = G(x) on ∂Ω ( 15
)
where G is a given vector-valued function. Then, the following asymptotic expansion holds:

S (v δ -v). ∂G ∂ν dσ - S ( ∂v δ ∂ν - ∂v ∂ν ).Gdσ = δ ∂D h(([M 0,1 -σ 1 ]D(v i ))τ.D(u i )τ )dσ +O(δ 2 ), ( 16
)
where the remainder O(δ 2 ) depends only on µ 0 , µ 1 , the C 1,1 -norm of h, and dist(S, ∂D). Here u i and v i denote the limits from inside D

The following theorem concerning the perturbations of the generalized viscosity moment tensors VMTs can be justified by the same manner as done in two dimensions, see [START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF]. Later, We will give briefly its proof.

Theorem 2.3. Assume that all hypothesis (1), ( 4), ( 7) and ( 8) are satisfied. Let's say that d = 2 and let M 0,1 be given by [START_REF] Depauw | Solutions des équations de Navier Stokes incompressibles dans un domaine exterieur[END_REF]. Suppose that (a α j ) and (b β k ) are constants such that

H = 2 j=1 α∈N 2 a α j x α e j and F = 2 k=1 β∈N 2 b β k x β e k satisfy -div(σ 0 D(•) -q Id) = 0 in R d , d = 2, 3.
Let u and v be the solutions to (15) and ( 6), respectively. Then, the following asymptotic expansion holds:

α,β,j,k a α j b β k M j αβk (D δ ) -M j αβk (D) = -δ ∂D h(([M 0,1 -σ 1 ]D(v i ))τ.D(u i )τ )dσ + O(δ 2 ), ( 17
)
where O(δ 2 ) depends only on µ 0 , µ 1 , the C 1,1 -norm of h, and dist(S, ∂D).

According to previous investigations in others EDPs such that in [START_REF] Ammari | Reconstruction of small interface changes of an inclusion from modal measurements II: The elastic case[END_REF][START_REF] Ammari | Reconstruction of closely spaced small inclusions[END_REF][START_REF] Brühl | A direct impedance tomography algorithm for locating small inhomogeneities[END_REF][START_REF] Khelifi | Asymptotic expansions for the voltage potentials with two and three-dimensional thin interfaces[END_REF][START_REF] Khelifi | Boundary voltage perturbations resulting from small surface changes of a conductivity inclusion[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF], the VMTs associated with an unknown inclusion can be detected from boundary measurements. It is well known that [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF] the model of VMTs has its complements in the conductivity as generalized polarization tensor GPTs and the elasticity as the elastic moment tensor EMTs. Following once again [START_REF] Ammari | Reconstruction of closely spaced small inclusions[END_REF][START_REF] Brühl | A direct impedance tomography algorithm for locating small inhomogeneities[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF], given an arbitrary shape, one can find an equivalent ellipsoid with the same polarization tensor. Therefore, recovering more shape details than the equivalent ellipsoid using a finite number of GPTs is an ambitious question. One of solutions is to extend that optimization approach in [START_REF] Ammari | The generalized polarization tensors for resolved imaging. Part I: Shape reconstruction of a conductivity inclusion[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF] which is based on the shape derivative of α,β a α b β M αβ (k, D), where a α and b β are constants such that H = α a α x α and F = β b β x β are harmonic polynomials. In order to calculate that shape derivative, we use the asymptotic expansion of the GPTs under small perturbations of the boundary of the inclusion D given in [START_REF] Khelifi | Boundary voltage perturbations resulting from small surface changes of a conductivity inclusion[END_REF].

Representation of solutions

We now develop a boundary integral formulation for solving the perturbed Stokes problem [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF]. The components of the fundamental Stokes tensor Γ = (Γ ij ) d i,j=1 and those of the associated pressure vector P = (P ij ) d i,j=1 , which determine the fundamental solution (Γ, P ) of the Stokes system in R d , are given for d = 3 by (see for instance [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF], [START_REF] Ladyzhenskaya | Mathematical theory of the viscous incompressible id[END_REF])

   Γ ij (λ, x) = -1 4π δij e i √ λ|x| |x| -1 4πλ ∂ xi ∂ xj ( e i √ λ|x| -1 |x| ) P i (x) = -1 4π xi |x| 3 , (18)
where ∂x i denotes ∂/∂x i for i = 1, • • • , d and i 2 = -1. We recall that the i th row Γ i of Γ satisfies

-∆Γ i + ∇P i (x) -λΓ i = e i δ(x) in R d ∇.Γ i = 0 in R d ; ( 19
)
in the sense of distributions and where (e

i ; i = 1, • • • , d) is the orthonormal basis of R d .
Note that we used the Einstein convention for the summation notation omitting the summation sign for the indices appearing twice.

The potential theory for the Stokes system

Let us denote by ϕ = (ϕ i ) 1≤i≤3 a complex vector-valued function with class C 0 (∂Ω).

The hydrodynamic single-layer potential with density ϕ ∈ C 0 (∂Ω) d is the vector function S(λ)ϕ(x) defined by

S(λ)ϕ(x) := ∂Ω Γ(λ, |x -y|)ϕ(y) dσ(y), x ∈ R d \∂Ω. (20)
The pressure term Q corresponding to the single layer potential is the function given by Qϕ(x) := ∂Ω P (x, y)ϕ(y) dσ(y), x ∈ R d \∂Ω.

The double layer potential of the density function ϕ(x) on L ( ∂D) associated with the Stokes system is defined by

D(λ)ϕ(x) = ∂D [∇ y Γ(λ, |x -y|)ν(y)] T + (∇ y Γ) T (x -y, λ)ν(y) ϕ(y)dσ(y) (21)
:

= ∂D K(x -y)ϕ(y)dσ(y), x ∈ R d \∂D.
For a careful study of these potentials, one can refer to [START_REF] Kohr | The interior Neumann problem for the Stokes resolvent system in a bounded domain in R n[END_REF], [START_REF] Ladyzhenskaya | Mathematical theory of the viscous incompressible id[END_REF], [START_REF] Varnhorn | The Stokes equations[END_REF].

Taking into account the well known properties of Green function Γ, one obtains the result that the pair (Sϕ, Qϕ) are smooth functions in each of the domains Ω and R d \Ω respectively. Also these functions are classical solutions to the Stokes system [START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF].

The continuity and jump relations of the Stokes surface potentials on the boundary ∂Ω are described in the following proposition (see [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF], [START_REF] Depauw | Solutions des équations de Navier Stokes incompressibles dans un domaine exterieur[END_REF] pp. 41-42 or [START_REF] Varnhorn | The Stokes equations[END_REF] p. 66 ): 2 and let S and D(λ) denote the surface potentials defined in ( 20)- [START_REF] Li | Estimates for elliptic systems from composite material[END_REF]. Then on the boundary ∂D the following continuity and jump relations are satisfied:

Proposition 3.1. For d = 2, let ϕ ∈ C 0 (∂Ω)
(S(λ)ϕ)

+ = (S(λ)ϕ)| -= S(λ)ϕ ∂S(λ)(ϕ) ∂ν (x) ± = ± ϕ(x) 2 + ∂Ω ∂Γ(λ, |x -y|) ∂ν(x) ϕ(y) dσ(y). := (± 1 2 I + K * D )ϕ(x) a.e. x ∈ ∂D, and D(λ)ϕ(x) = (∓ 1 2 I + K D )ϕ(x) a.e. x ∈ ∂D, (22)
where K D is defined by

K D ϕ(x) := p.v.
∂D := {x = X(t), t ∈ [a, b]}.
Then the outward unit normal to D, ν(x), is given by ν

(x) = R -π 2 X (t)
, where R -π 2 is the rotation by -π/2, the tangential vector at x, τ (x) = T (x) = X (t), and X (t) ⊥ X (t). Set the curvature ζ(x) to be defined by

X (t) = ζ(x)ν(x).
We will sometimes use h(t) for h(X(t)) and h (t) for the tangential derivative of h(x). Then, X

(t) = X(t) + δh(t)ν(x) = X(t) + δh(t)R -π 2 X (t) is a parametrization of D δ . By ν(x)
, we denote the outward unit normal to D δ at x. Then, it is proved in [START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF][START_REF] Ammari | Conductivity interface problems. Part I: small perturbations of an interface[END_REF][START_REF] Khelifi | Asymptotic expansions for the voltage potentials with two and three-dimensional thin interfaces[END_REF] that

ν(x) = ν(x) -δh (t)τ (x) + O(δ 2 ), uniformly for x ∈ ∂D. (23)
The following uniformly convergent expansion for the length element dσ δ (ỹ) can be obtained easily (see [START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF][START_REF] Ammari | Conductivity interface problems. Part I: small perturbations of an interface[END_REF][START_REF] Khelifi | Boundary voltage perturbations resulting from small surface changes of a conductivity inclusion[END_REF]):

dσ δ (ỹ) = (1 -δζ(y)h(y) + O(δ 2 )).dσ(y). ( 24
)
Let Ψ δ be the diffeomorphism from ∂D onto ∂D δ given by Ψ δ (x) = x + δh(t)ν(x), where x = X(t). The asymptotic expansion of S * D δ is given by the following lemma, see [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF][START_REF] Khelifi | Asymptotic expansions for the voltage potentials with two and three-dimensional thin interfaces[END_REF] for an idea of the proof. Lemma 3.1. Let N ∈ N. There exists C depending only on N , the C 2,1 -norm of ∂D and the C 1,1 -norm of h such that for any φ ∈ L 2 (∂D δ ),

||(S * D δ φ) • Ψ δ -S * D φ - N -1 n=0 δ n+1 S (n+1) D φ|| L 2 (∂D) ≤ Cδ N +1 ||φ|| L 2 (∂D) ,
where φ := φ • Ψ δ . The operators S

(n) D , n > 0, are compact on L 2 (∂Ω) and can be explicitly computed.

Representation of solutions

The following transmission condition follows. Proposition 3.2. Let v be the solution to [START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF]. Let v i and v e be the inner and the outer components of v. Then, the following identities hold:

(σ 0 D(v e ))τ = (M 0,1 D(v i ))τ, ( 25
) (σ 1 D(v i ))τ = (M 1,0 D(v e ))τ, ( 26
)
where the tensor M 1,0 was introduced previously in ( 14).

Proof 3.1. According to [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF], the solution v of (6) may satisfy the following transmission conditions along the interface ∂D:

v i = v e , (27)
∇v e τ = ∇v i τ, (28)
< D(v i )τ, τ >=< D(v e )τ, τ >, (29)
and so,

p| + + 2µ 1 < D(v i )n, n >= p| -+ 2µ 0 < D(v e )n, n >, (30)
µ 1 < D(v i )n, τ >= µ 0 < D(v e )n, τ > . (31)
Therefore,

(σ 0 D(v e ))τ = 2µ 0 (D(v e )τ ) = 2µ 0 < (D(v e )τ, τ > τ + 2µ 0 < (D(v e )τ, n > τ = 2µ 0 < (D(v i )τ, τ > τ + 2µ 1 < (D(v i )τ, n > τ.
On the other hand, we can expand:

(D(v i )τ =< (D(v i )τ, n > n+ < (D(v i )τ, τ > τ.
Thus,

(σ 0 D(v e ))τ = 2µ 0 < (D(v i )τ, τ > τ + 2µ 1 D(v i )τ -2µ 1 < (D(v i )τ, τ > τ = 2µ 1 D(v i )τ + 2(µ 0 -µ 1 ) < (D(v i )τ, τ > τ = (M 0,1 D(v i ))τ.
By a completely similar manner, one may obtain the desired condition (26).

To find a representation to solutions of ( 10) or ( 6), the reader may see [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF][START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF][START_REF] Ammari | Layer Potential Techniques in Spectral Analysis[END_REF][START_REF] Daveau | Asymptotic formula for the solution of the Stokes problem with a small perturbation of the domain in two and three dimensions[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF], 10] and others. The solution (u δ , q δ ) to [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF] can be represented by

u δ = H(x) + S D δ [ϕ δ ](x) , x ∈ Ω \ D δ , S D δ [ψ δ ](x) x ∈ D δ , (32)
p δ = V D δ [ϕ δ ](x) , x ∈ Ω \ D δ , V D δ [ψ δ ](x) x ∈ D δ , (33)
where the pair (ϕ δ , ψ δ ) is the unique solution in L 2 (∂D δ ) × L 2 (∂D δ ) of

S D δ [ψ δ ] | --S D δ [ψ δ ] | + = H on ∂D δ , (-1 2 I + K * D δ )[ψ δ ] -(-1 2 I + K * D δ )[ϕ δ ] = ∂H ∂ν on ∂D δ . (34)
But, the solution (u, p) to ( 6) can be represented as

u(x) = H(x) + S D [ϕ](x) , x ∈ Ω \ D, S D [ψ](x) x ∈ D, (35)
p(x) = V D [ϕ](x) , x ∈ Ω \ D, V D [ψ](x) x ∈ D, (36)
where the pair (ϕ, ψ) is the unique solution in

L 2 (∂D δ ) × L 2 (∂D δ ) of S D [ψ] | --S D [ϕ] | + = H on ∂D, (-1 2 I + K * ∂D )[ψ] -(-1 2 I + K * ∂D )[ϕ] = ∂H ∂ν on ∂D, (37)
For more details on the previous representations formulae, we may refer to [START_REF] Ammari | Reconstruction of Small Inhomogeneities from Boundary Measurements[END_REF].

Now, we introduce the integral operator K

(1) D , to be defined for any φ ∈ L 2 (∂D) by:

K (1) D [φ](x) := ∂ K 1 (x -y)φ(y)dσ(y), x ∈ ∂D, (38)
where K 1 (x -y) means the first order term in the asymptotic expansion of K(x -ỹ) with K is the kernel of the double layer potential D(λ). The following result may be useful in our next analysis.

Lemma 3.2. Assume that all hypothesis (1), ( 4), ( 7) and ( 8) are satisfied. Then if we introduce (ψ (1) , ϕ (1) ) to be the solution to the following system:

S D [ψ (1) ] | --S D [ϕ (1) ]

| + = h ∂H ∂n -( S (1)
D [ψ] -S (1)
D [ϕ]) on ∂D, (-1 2 I + K * D )[ψ (1) ]) -( 1 2 I + K * D [ϕ (1) ] = κh ∂H ∂ν -∂ ∂τ (h(σ 0 ∇H)τ ) -( K (1)
D [ψ] -K (1)
D [ϕ]) on ∂D, (39
) where (ψ, ϕ) is the solution to (37). Then, we have the following relations:

( S D -S (1) D )[ ψ-ψ-δψ (1) ] | --(S D +δS (1) D )[ ϕ-ϕ-δϕ (1) ] | + = H•φ δ -H-δh ∂H ∂n +O(δ 2 ) on ∂D, (40) and
(- 1 2 I+ K * D +δ K ( 1) D )[ ψ-ψ-δψ (1) ] | --(- 1 2 I+K * D +δK ( 1) D )[ ϕ-ϕ-δϕ (1) ] | + = ∂H ∂ν •φ δ - ∂H ∂ν -δh ∂H ∂ν (41) +δ ∂ ∂τ (h(σ 0 ∇H)τ ) + O(δ 2 ) on ∂D.
Here, the operator S

D is defined for any φ ∈ L 2 (∂D) by

S

(1)

D [φ](x) = -S D [κhφ](x) + (h(x) ∂S D [φ] ∂n (x) + D D [hφ](x)) | ± , x ∈ ∂D (42)
with D D means the standard double layer potential, and K

D is given by (38).

Now we try to check formally that the obtained formula for V M T is the limit of the one obtained in case of elastic medium when the compressional modulus goes to infinity. According to [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF], we define the viscous moment tensor (VMT) denoted by M = (M pq ij ij) i,j,p,q=1,2,3 by

M pq ij = ( µ -µ) D ∇v pq : (∇.(ζ i e j ) + ∇(ζ i e j ) T )dζ, (43)
                   µ∆v pq + ∇p = 0 in R d \ D, µ∆v pq + ∇p = 0 in D, vpq | -= vpq | + on ∂D, (pN + µ ∂ vpq ∂N ) | -= (pN + µ ∂ vpq ∂N ) | + on ∂D, ∇.v pq = 0 in R d , vpq (ζ) -ζ p e q + δ pq d(ζ) → 0 as |ζ| → ∞, p(ζ) → 0 as |ζ| → ∞. (44)
Now, according to [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF] the elastic moment tensor M = (m pq ij ) associated with the domain D and the pairs of Lamé parameters (λ, λ; µ, µ) is given by

m pq ij = D - ∂((ζ i e j )) ∂ν + ∂((ζ i e j )) ∂ ν w pq dσ, i, j, p, q = 1 • • • d. (45)
where w pq means a solution of the Lamé system. Then the following limit holds (formally).

Proposition 3.3.

[3] For i, j, p, q = 1 • • • d, we have:

M pq ij = lim λ, λ→∞ [m pq ij -δ ij /d d k=1 m pq kk -δ pq /d d s=1 m ss ij + δ ij δ pq d 2 d k,s=1 m ss kk ]. Proof 3.2. Following [3, page 2], one can remark that - ∂((ζ i e j )) ∂ν + ∂((ζ i e j )) ∂ ν = -(λ -λ)∇.(ζ i e j )N -(µ -µ) (ζ i e j ) ∂N .
By a closely approach used to prove Theorem 3.1 in [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF], one can prove that

vpq = lim λ, λ→∞ (w pq -δ pq /d k w kk ).
Then, the proof is achieved by using (43) and expending

[m pq ij -δ ij /d d k=1 m pq kk = ( µ -µ) D ∇w pq .[∇.(ζ i e j -δ ij d(ζ)) + ∇.(ζ j e i -δ ij d(ζ)) T ].

Proof of the Theorem 2.1

Regarding Section 3.2, we can easily write that,

∂ ∂t (φ(x)) = ∇φ(x)X (t) = ∂φ ∂t (x). (46)
Then, the following result easily holds.

Lemma 3.3. Suppose that we have all hypothesis of Theorem 2.1. Using Section 3.2, then the following holds:

div(σ 0 D(φ(x))) = 1 2 µ 0 div(∇φ(x) + (∇φ) (x)) = 1 2 µ 0 ∇∇φ(x)n(x)n(x) + ∇∇φ(x)τ (x)τ (x) + ∇(∇φ) (x)n(x)n(x) + ∇(∇φ) (x)τ (x)τ (x) .
Using relation (46), we get that

div(σ 0 D(φ(x))) = 1 2 µ 0 ∇∇φ(x)n(x)n(x) + ∇(∇φ(x)) n(x)n(x) + 1 2 ∂ ∂t σ 0 D(φ(x) τ (x) = 1 2 µ 0 ∇∇φ(x)n(x)n(x) + ∇(∇φ) )(x)n(x)n(x) -κ(x) σ 0 D(φ(x)) n(x) + ∂ ∂t σ 0 D(φ(x))τ (x) x ∈ ∂D.
Next, regarding Section 3.2, for

x = x + δh(x)ν(x) ∈ ∂D δ , ν( x) = ∞ m=0 δ m ν m (x), x ∈ ∂D ν 0 (x) = n(x), n 1 (x) = -h (x)τ (x), x ∈ ∂D.
Consequently, the following Taylor expansion holds:

H(x + δh(x)n(x)) = H(x) + δh(x) ∂H ∂n (x) + o(δ 2 ), x ∈ ∂D. (47)
On the other hand, taking into account relation [START_REF] Ammari | Reconstruction of closely spaced small inclusions[END_REF], we get:

∂H ∂ν ( x) = µ 0 ∂H ∂n (x) -p n( x) = µ 0 ∇H( x) + (∇H) ( x) n( x) -p n( x) = µ 0 ∇H(x)n(x) + δh(x)∇∇H(x)n(x) + δh(x)∇(∇H) (x)n(x) +(∇H) (x)n(x) + δh(x)∇(∇H) (x)n(x) + δh(x)(∇(∇H) ) (x)n(x) n(x) -p n( x) = µ 0 ∇H(x)n(x) + δh(x)∇∇H(x)n(x) + (∇H) (x)n(x) + δh(x)∇(∇H) (x)n(x) n(x) -µ 0 δh (x) ∇H(x)τ (x) + (∇H) (x)τ (x) -p n(x) + δh (x)p τ (x) = ∂H ∂ν (x) -2δh σ 0 D(H)(x)τ (x) + δh (x)p τ (x) + µ 0 δh(x) ∇∇H(x)n(x)n(x) + ∇(∇H) (x)n(x)n(x) .
Then, by Lemma 3.3, we immediately get

∂H ∂ν ( x) = ∂H ∂ν (x) -δh 2σ 0 D(H)(x)τ (x) -p τ (x) + δh(x) 2div(σ 0 D(H(x)) + 2κ(x)(σ 0 D(H(x)))n(x) -2 ∂ ∂t (σ 0 D(φ(x))τ (x))) + o(δ 2 ) = ∂H ∂ν (x) + δh p τ (x) + δh(x) 2div(σ 0 D(H(x)) + 2κ(x)(σ 0 D(H(x)))n(x) -δ ∂ ∂t (2h(x)σ 0 D(φ(x))τ (x))) + o(δ 2 ) = ∂H ∂ν + δh(x)(div(σ 0 D(φ(x)) + κ(x)(σ 0 D(φ(x)))n(x) -δ ∂ ∂t (hσ 0 D(φ(x))τ (x)))) + o(δ 2 ), x ∈ ∂D.
Next, one can easily see that -div(σ 0 D(v 1 )-p 1 Id) = 0 in Ω\D and -div(σ 1 D(v 1 )-p 1 Id) = 0 in D. This may be seen if we introduce the asymptotic expansion of v δ with respect to δ → 0 into the transmission problem [START_REF] Daveau | Asymptotic behaviors for eigenvalues and eigenfunctions associated to Stokes operator in the presence of small boundary perturbations[END_REF].

On the other hand, we may write:

v i 1 -v e 1 = -S D [ϕ (1) ] + S D [κhϕ] -D D [hϕ] + S D [ψ (1) ] -S D [κhψ] + D D [hψ]. (48)
Using relation (39), we get

v i 1 -v e 1 = h ∂H ∂n + S (1)
D [ϕ] -S (1)
D [ψ] + S D [κhϕ] -S D [κhψ] -D D [hϕ] + D D [hψ].
But, by (42), we have

v i 1 -v e 1 = h( ∂H ∂n + ∂S D [ϕ] ∂n | + - ∂ S D [ψ] ∂n | -),
and so, by (35) we find that

v i 1 -v e 1 = h(∇v e n -∇v i n) = 0.
To proceed with the proof, we may expand

∂v 1 ∂ν | -- ∂v 1 ∂ν | + = - ∂S D [ϕ (1) ] ∂ν | + + ∂S D [κhϕ] ∂ν | + - ∂D D [hϕ] ∂ν | + + ∂ S D [ψ (1) ] ∂ν | - = - ∂ S D [κhψ] ∂ν | -+ ∂ D D [hψ] ∂ν | -.
Inserting relation (39) into previous equality, we immediately obtain:

∂v 1 ∂ν | -- ∂v 1 ∂ν | + = κh ∂H ∂ν - ∂ ∂τ (hσ 0 (D(H))τ ) -K (1)
D [ψ] + K (1) D [ϕ] + ∂S D [κhϕ] ∂ν | + - ∂ S D [κhψ] ∂ν | - + ∂ D D [hψ] ∂ν | -- ∂D D [hϕ] ∂ν | + .
Finally, relations (48) and ( 25)-(26) imply that,

∂v 1 ∂ν | -- ∂v 1 ∂ν | + = ∂ ∂τ (hσ 1 (D(v i ))τ ) - ∂ ∂τ (hσ 0 (D(v e ))τ ) = ∂ ∂τ (h[σ 1 -M 0,1 ](D(v i ))τ ),
and the proof can be easily obtained.

Proof of the Theorem 2.2

Let v and v δ be the solutions to ( 6) and ( 10),respectively and let v 1 be the unique solution of (13).

Let Ω be a bounded region outside the inclusion D and away from ∂D. For x ∈ Ω, the following pointwise asymptotic expansion holds:

∂v δ ∂ν (x) = ∂v ∂ν (x) + δ ∂v 1 ∂ν (x) + O(δ 2 ), (49)
where the remainder O(δ 2 ) depends only on λ 0 , λ 1 , µ 0 , µ 1 , the C 2 -norm of X, the C 2 -norm of h, and dist(Ω, ∂D).

Let S be a Lipschitz closed curve enclosing D away from ∂D and let u be the solution to (15). As we done in [START_REF] Boujemaa | Small perturbation of a surface: full Maxwells equations[END_REF], and by integration by parts one may combine both relations 11) and (49) to get that,

S (v δ -v). ∂F ∂ν dσ - S ( ∂v δ ∂ν - ∂v ∂ν ).F dσ = δ S (v 1 . ∂u ∂ν - ∂v 1 ∂ν .u)dσ + O(δ 2 ) = δ S ( ∂u e ∂ν .v e 1 -u e . ∂v e 1 ∂ν )dσ + O(δ 2 ).
Taking into account relation [START_REF] Brühl | A direct impedance tomography algorithm for locating small inhomogeneities[END_REF], we immediately get

S (v δ -v). ∂F ∂ν dσ - S ( ∂v δ ∂ν - ∂v ∂ν ).F dσ = δ ∂D ( ∂u i ∂ ν .v i 1 -u i . ∂v i 1 ∂ ν )dσ =0 + ∂D ∂ ∂τ (h([σ 1 -M 0,1 ]D(v i ))τ ).u i dσ + O(δ 2 ) = -δ ∂D (h([σ 1 -M 0,1 ]D(v i ))τ ).∇u i τ dσ + O(δ 2 ) = -δ ∂D (h([σ 1 -M 0,1 ]D(v i ))τ ).D(u i )τ dσ + O(δ 2 )
which achieves the proof.

Proof of the Theorem 2.3

To prove Theorem 3, we may firstly introduce the notion of the viscosity moment tensors (VMTs) associated with the inclusion D as follows (see for example [START_REF] Ammari | Reconstruction of Small Inhomogeneities from Boundary Measurements[END_REF][START_REF] Ammari | Effective viscosity properties of dilute suspensions of arbitrarily shaped particles[END_REF][START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF]):

For multi-index α ∈ N 2 and j = 1, 2, let the pair (f j α , g j α ) ∈ (L 2 (∂D)) 2 be the unique solution (as seen for ( 49)) to:

S D [f j α ] | --S D [g j α ] | + = x α e j on ∂D, ∂ S D [f j α ] ∂ν | -- ∂S D [g j α ] ∂ν | + = ∂(x α ej ) ∂ν on ∂D. (50)
For some multi-index β ∈ N 2 , the VMTs denoted by M = (M j αβ ) j,αβ are defined by: M j αβ = (m j αβ,1 , m j αβ,2 ) := ∂D y β g j α (y)dσ(y).

(51)

One can see [START_REF] Ammari | A method of biological tissues elasticity reconstruction using magnetic resonance elastography measurements[END_REF] for more details about this definition. Moreover, we recall the readers that similar relations as (51) was introduced before for generalized polarization tensors GPTs in [START_REF] Ammari | The generalized polarization tensors for resolved imaging. Part I: Shape reconstruction of a conductivity inclusion[END_REF][START_REF] Khelifi | Boundary voltage perturbations resulting from small surface changes of a conductivity inclusion[END_REF] and for elastic moments tensors EMTs in [START_REF] Lagha | Small perturbations of an interface for elastostatic problems[END_REF].

where ϕ was defined in (37). But the perturbed m j αβk (D δ ) is given by:

where ϕ δ was defined in (34).

By relation [START_REF] Ammari | Conductivity interface problems. Part I: small perturbations of an interface[END_REF] and by Taylor expansion (as δ → 0), we can write: G(ỹ) = F(y + δh(y)n(y)) = F(y) + δh(y) ∂F ∂n (y) + O(δ 2 ); y ∈ ∂D.

To prove the asymptotic behavior of αβjk a α j b β k m j αβk (D δ ) (as δ → 0), we may use also (53) and the fact that ϕ δ (ỹ) = F(y + δh(y)n(y)) = ϕ(y) + δh(y)ϕ (1) (y) + O(δ 2 ); y ∈ ∂D, where ϕ (1) was defined in (48). Using the asymptotic expansion of dσ δ (ỹ) given in Section 3.2, then by change of variables ỹ = y + δh(y)n(y), we use (53) to immediately obtain: Analyzing each part of the left hand side of (55) one can find by using (15), Proposition 2 and ( 14) that: Then the proof can be deduced easily.

∂DK

  (x -y)ϕ(y)dσ(y) a.e. x ∈ ∂D, and K * D is the adjoint operator of K D . 3.2 Small perturbation of a C 1,1 -interface Let a, b ∈ R, with a < b, and let X(t) : [a, b] → R 2 be the arclength parametrization of ∂D, namely, X is an C 2 -function satisfying |X (t)| = 1 for all t ∈ [a, b] and

  where for ζ ∈ R d , we define d(ζ) := 1/d k ζ k e k and vpq be the solution to

Now let H = 2 j=1 α∈N 2 a α j x α e j and F = 2 k=1 β∈N 2 bF

 222 β k x β e k satisfying -div(σ 0 D(•) -q Id) = 0 in R d .The VMTs m j αβk (D) associated with D satisfy: (y)ϕ(y)dσ(y),

F

  (ỹ)ϕ δ (ỹ)dσ δ (ỹ),

F • ϕ ( 1 )

 1 ϕ(y) + δh(y)ϕ(1) (y) 1 -δτ h dσ + O(δ 2 ).(54) Thus,αβjk a α j b β k m j αβk (D δ )αβjk a α j b β k m j αβk (D) = δ ∂D h(y) ∂F ∂n (y)ϕ(y)dσ(y) +δ ∂D (y) -τ hϕ(y) dσ(y) + O(δ 2 ).On the other hand, using Proposition 1, we can write∂D F • ϕ (1) (y) -τ hϕ(y) dσ(y) = -∂D F • ∂S D [ϕ (1) (y) -τ hϕ(y)] ∂ν | -(55)-∂S D [ϕ(1) (y) -τ hϕ(y)] ∂ν | + dσ(y).

∂DF • ϕ ( 1 )

 1 (y) -τ hϕ(y) dσ(y) = ∂D h(([σ 1 -M 0,1 ]D(v i ))τ.D(u i )τ )dσ.

Acknowledgement

The authors wish to thank the editor and the referees for their valuable criticisms and suggestions, leading to the present improved version of our paper.


Christian Daveau
email: christian.daveau@math.u-cergy.fr

Abdessatar Khelifi
email: abdessatar.khelifi@fsb.rnu.tn

Anton Sushchenko
email: anton.sushchenko@ensea.fr

Reconstruction of closely spaced small inhomogeneities via boundary measurements for the full time-dependent Maxwell's equations

Keywords: Maxwell's equations, inhomogeneities, inverse problem, reconstruction, geometric control 2000 AMS subject classifications. 35R30, 35B40, 35B37, 78M35

We consider for the full time-dependent Maxwell's equations the inverse problem of identifying locations and certain properties of small electromagnetic inhomogeneities in a homogeneous background medium from dynamic boundary measurements on the boundary for a finite time interval.

Introduction

The ultimate objective of the work described in this paper is to determine locations and certain properties of the shapes of small electromagnetic inhomogeneities in a homogeneous background medium from dynamic boundary measurements on part of the boundary and for finite interval in time. Using as weights particular background solutions constructed by a geometrical control method we develop an asymptotic method based on appropriate averaging of the partial dynamic boundary measurements. For stationary Maxwell's equations it has been known that the Dirichlet to Neumann map uniquely determines (smooth) isotropic electromagnetic parameters, see [START_REF] Mcdowall | An electromagnetic inverse problem in chiral media[END_REF], [START_REF] Romanov | Inverse Problems for Maxwell's Equations, Inverse and Ill-posed Problems Series[END_REF], [START_REF] Somersalo | A linearized inverse boundary value problem for Maxwell's equations[END_REF]. We will provide in this paper a rigorous derivation of the inverse Fourier transform of a linear combination of derivatives of point masses, located at the positions z j of the inhomogeneities, as the leading order term of an appropriate averaging of (partial) dynamic boundary measurements of the tangential components of electric fields on part of the boundary. Our formulas may be used to determine properties (location, relative size ) of the small inhomogeneities in case a single, or a few (tangential) boundary electric fields are known. Our approach differs from [START_REF] Ammari | An inverse initial boundary value problem for the wave equation in the presence of imperfections of small volume[END_REF], [START_REF] Ammari | Reconstruction of closely spaced small inclusions[END_REF], [START_REF] Ammari | Boundary integral formulas for the reconstruction of electromagnetic imperfections of small diameter[END_REF], [START_REF] Ammari | Asymptotic formulas for perturbations in the electromagnetic fields due to the presence of inhomogeneities of small diameter II.The full Maxwell equations[END_REF], [START_REF] Volkov | An inverse problem for the time harmonic Maxwell equations[END_REF] and is expected to lead to very effective computational identification algorithms.

For discussions on closely related (stationary) identification problems we refer the reader to [START_REF] Sylvester | A global uniqueness theorem for an inverse boundary value problem[END_REF], [START_REF] Vogelius | Asymptotic formulas for perturbations in the electromagnetic fields due to the presence of inhomogeneities[END_REF], [START_REF] Cedio-Fengya | Identification of conductivity imperfections of small diameter by boundary measurements. Continuous dependence and computional reconstruction[END_REF], and [START_REF] Friedman | Identification of small inhomogeneities of extreme conductivity by boundary measurements: a theorem on continuous dependence[END_REF].

Problem formulation

Let Ω be a bounded C 2 -domain in R d , d = 2, 3. Assume that Ω contains a finite number of inhomogeneities, each of the form z j + αB j , where B j ⊂ R d is a bounded, smooth domain containing the origin. The total collection of inhomogeneities is B α = ∪ m j=1 (z j + αB j ). The points z j ∈ Ω, j = 1, . . . , m, which determine the location of the inhomogeneities, are assumed to satisfy the following inequalities: |z j -z l | ≥ c 0 > 0, ∀ j = l and dist(z j , ∂Ω) ≥ c 0 > 0, ∀ j.

(

) 1
Assume that α > 0, the common order of magnitude of the diameters of the inhomogeneities, is sufficiently small, that these inhomogeneities are disjoint, and that their distance to R d \ Ω is larger than c 0 /2. Let µ 0 and ε 0 denote the permeability and the permittivity of the background medium, and assume that µ 0 > 0 and ε 0 > 0 are positive constants. Let µ j > 0 and ε j > 0 denote the permeability and the permittivity of the j-th inhomogeneity, z j + αB j , these are also assumed to be positive constants.

Introduce the piecewise-constant magnetic permeability

µ α (x) = µ 0 , x ∈ Ω \ Bα , µ j , x ∈ z j + αB j , j = 1 . . . m. ( 2
)
If we allow the degenerate case α = 0, then the function µ 0 (x) equals the constant µ 0 . The electric permittivity is defined by ε α (x) = ε 0 , for all x ∈ Ω. Let n = n(x) denote the outward unit normal vector to Ω at a point on ∂Ω, ∂ t u = ∂u ∂t and ∆ means the Laplace operator defined by ∆u =

d i=1 ∂ 2 u ∂x 2 i .
In this paper, we will denote by bold letters the functional spaces for the vector fields. Thus H s (Ω) denotes the usual Sobolev space on Ω and H s (Ω) denotes (H s (Ω)) d and L 2 (Ω) denotes (L 2 (Ω)) d . As usual for Maxwell equations, we need spaces of fields with square integrable curls:

H( curl ; Ω) = {u ∈ L 2 (Ω), curl u ∈ L 2 (Ω)},
and with square integrable divergences

H(div ; Ω) = {u ∈ L 2 (Ω), div u ∈ L 2 (Ω)}.
We will also need the following functional spaces:

Y (Ω) = {u ∈ L 2 (Ω), div u = 0 in Ω}, X(Ω) = H 1 (Ω) ∩ Y (Ω),
and T L 2 (∂Ω) the space of vector fields on ∂Ω that lie in L 2 (∂Ω). Finally, the "minimal" choice for the electric variational space would be

X N (Ω) = {v ∈ H( curl ; Ω) ∩ H(div ; Ω); v × n = 0 on ∂Ω}.
Now, we introduce the following time-dependent Maxwell equations (associated to the electric field)

               (ε α ∂ 2 t + curl 1 µ α curl )E α = 0 in Ω × (0, T ), div (ε α E α ) = 0 in Ω × (0, T ), E α | t=0 = ϕ, ∂ t E α | t=0 = ψ in Ω, E α × n| ∂Ω×(0,T ) = f, (3)
where E α ∈ R d is the electric field, f the boundary condition for E α × n, and ϕ and ψ the initial data.

Let E be the solution of the Maxwell's equations in the homogeneous domain:

               (ε 0 ∂ 2 t + curl 1 µ 0 curl )E = 0 in Ω × (0, T ), div (ε 0 E) = 0 in Ω × (0, T ), E| t=0 = ϕ, ∂ t E| t=0 = ψ in Ω, E × n| ∂Ω×(0,T ) = f. (4)
Here T > 0 is a final observation time and ϕ, ψ ∈ C ∞ (Ω) and f ∈ C ∞ (0, T ; C ∞ (∂Ω)) are subject to the compatibility conditions

∂ 2l t f | t=0 = (∆ l ϕ) × n| ∂Ω and ∂ 2l+1 t f | t=0 = (∆ l ψ) × n| ∂Ω , l = 1, 2, . . .
it follows that (4) has a unique solution E ∈ C ∞ ([0, T ] × Ω). It is also known (see for example [START_REF] Nicaise | Exact boundary controllability of Maxwell's equations in heteregeneous media and an application to an inverse source problem[END_REF]) that since Ω is smooth (C 2 -regularity would be sufficient) the non homogeneous Maxwell's equations (3) have a unique weak solution E α ∈ C 0 (0, T ; X(Ω)) ∩ C 1 (0, T ; L 2 (Ω)). Indeed, curl E α belongs to C 0 (0, T ; X(Ω)) ∩ C 1 (0, T ; L 2 (Ω)).

Asymptotic formula

We start the derivation of the asymptotic formula for curl E α × n with the following estimate.

Lemma 3.1

The following estimate as α → 0 holds:

||∂ t (E α -E)|| L ∞ (0,T ;L 2 (Ω)) + ||E α -E|| L ∞ (0,T ;X N (Ω)) ≤ Cα, ( 5
)
where the constant C is independent of α and the set of points {z j } m j=1 provided that assumption (1) holds.

Proof. From (3)-(4), it is obvious that E α -E ∈ X N (Ω), then due to the Green formula we have for any v ∈ X N (Ω):

Ω ε 0 ∂ 2 t (E α -E) • v dx + Ω 1 µ α curl (E α -E) • curl v dx = (6) m j=1 ( 1 µ 0 - 1 µ j ) z j +αB j curl E • curl v dx.
Let v α be defined by

v α ∈ X N (Ω), curl 1 µα curl v α = ∂ t (E α -E) in Ω. (7)
Then,

Ω 1 µ α curl (E α -E) • curl v α dx = - Ω ∂ t (E α -E) • (E α -E) dx = - 1 2 ∂ t Ω |E α -E| 2 dx
and by Green formula, relation [START_REF] Daveau | An inverse problem for the time dependent Maxwell's equations in the presence of imperfections of small volume[END_REF] gives:

Ω ∂ 2 t (E α -E) • v α dx = Ω curl 1 µ α curl ∂ t v α • v α dx = - Ω 1 µ α curl ∂ t v α • curl v α dx = - 1 2 ∂ t Ω 1 µ α | curl v α | 2 dx.
Thus, it follows from (6) that

ε 0 ∂ t Ω 1 µ α | curl v α | 2 dx + ∂ t Ω |E α -E| 2 dx = -2 m j=1 ( 1 µ 0 - 1 µ j ) z j +αB j curl E • curl v α dx. Next, | m j=1 ( 1 µ 0 - 1 µ j ) z j +αB j curl E • curl v α | ≤ C|| curl E|| L 2 (Bα) || curl v α || L 2 (Ω) . Since E ∈ C ∞ ([0, T ] × Ω) we have || curl E|| L 2 (Bα) ≤ || curl E|| L ∞ (Bα) α( m j=1 |B j |) 1 2 ≤ Cα, which gives | m j=1 ( 1 µ 0 - 1 µ j ) z j +αB j curl E • curl v α dx| ≤ Cα|| curl v α || L 2 (Ω)
and so,

ε 0 ∂ t Ω 1 µ α | curl v α | 2 dx+∂ t Ω |E α -E| 2 dx ≤ Cα( Ω 1 µ α | curl v α | 2 dx+ Ω |E α -E| 2 dx) 1/2 . (8) From the Gronwall Lemma it follows that ( Ω 1 µ α | curl v α | 2 dx) 1/2 + ( Ω |E α -E| 2 dx) 1/2 ≤ Cα. ( 9
)
Combining this last estimate [START_REF] Evans | Partial Differential Equations[END_REF] with the fact that

||∂ t (E α -E)|| L ∞ (0,T ;H -1 (Ω)) ≤ C|| curl v α || L ∞ (0,T ;L 2 (Ω))
the following estimate holds

||E α -E 0 || L ∞ (0,T ;L 2 (Ω)) + ||∂ t (E α -E 0 )|| L ∞ (0,T ;L 2 (Ω)) ≤ Cα. ( 10
)
Now, taking (formally 6) we arrive at

) v = ∂ t (E α -E) in (
ε 0 ∂ t Ω |∂ t (E α -E)| 2 + 1 µ α | curl (E α -E)| 2 dx = 2 m j=1 ( 1 µ 0 - 1 µ j ) z j +αB j curl E • curl ∂ t (E α -E) dx.
By using the regularity of E in Ω and estimate [START_REF] Friedman | Identification of small inhomogeneities of extreme conductivity by boundary measurements: a theorem on continuous dependence[END_REF] given above, we see that

| m j=1 ( 1 µ 0 - 1 µ j ) z j +αB j curl E • curl ∂ t (E α -E) dx| ≤ C|| curl E|| H 2 (Bα) ||∂ t (E α -E)|| H -1 (Ω) ≤ Cα 2 ,
where C is independent of t and α, and so, we obtain

∂ t Ω |∂ t (E α -E)| 2 + 1 µ α | curl (E α -E)| 2 dx ≤ Cα 2
which yields the following estimate

||∂ t (E α -E)|| L ∞ (0,T ;L 2 (Ω)) + ||E α -E|| L ∞ (0,T ;X N (Ω)) ≤ Cα,
where C is independent of α and the points {z j } m j=1 .

Now, we can estimate curl E α -curl E 0 as follows.

Proposition 3.1 Let E α and E be solutions to the problems ( 3) and ( 4) respectively. There exist constants 0 < α 0 , C such that for 0 < α < α 0 the following estimate holds:

|| curl (E α -E 0 )|| L ∞ (0,T ;L 2 (Ω)) ≤ Cα, ( 11
)
Proof. To prove estimate [START_REF] Kime | Boundary controllability of Maxwell's equations in a spherical region[END_REF] it is useful to introduce the following function

v(x) = T 0 v(x, t)z(t) dt ∈ L 2 (Ω), (12)
where

v ∈ L 1 (0, T ; L 2 (Ω)) and z(t) is a given function in C ∞ 0 (]0, T [). Then, Ê(x) = T 0 E(x, t)z(t) dt and Êα (x) = T 0 E α (x, t)z(t) dt ∈ X(Ω),
which by relation [START_REF] Bardos | Sharp sufficient conditions for the observation, control and stabilization of waves from the boundary[END_REF] give

             ( Êα -Ê) ∈ H 1 (Ω), curl curl ( Êα -Ê) = 0(α) in Ω, div ( Êα -Ê) = 0 in Ω, ( Êα -Ê) × n| ∂Ω = 0,
and so,

|| curl ( Êα -Ê)|| L 2 (Ω) = O(α). ( 13
)
The fact that curl (E α -E) belongs to L ∞ (0, T ; L 2 (Ω)) and by using estimate [START_REF] Lagnese | Exact boundary controllability of Maxwell's equations in a general region[END_REF] we deduce that

Ω | curl E α (x, t) -curl E(x, t)| 2 dx = O(α 2 ) a.e. in t ∈ (0, T ), which means that || curl (E α -E)|| L 2 (Ω) = O(α) a.e. in t ∈ (0, T ).
Thus, estimate [START_REF] Kime | Boundary controllability of Maxwell's equations in a spherical region[END_REF] follows immediately if we take the sup on t ∈ (0, T ) in the last relation.

Before formulating our main result in this section, let us denote Φ j , j = 1, . . . , m the unique vector-valued solution to

                   ∆Φ j = 0 in B j , and R d \ B j , Φ j is continuous across ∂B j , µ j µ 0 ∂Φ j ∂ν j | + - ∂Φ j ∂ν j | -= -ν j , lim |y|→+∞ |Φ j (y)| = 0, (14)
where ν j denotes the outward unit normal to ∂B j , and superscripts -and + indicate the limiting values as the point approaches ∂B j from outside B j , and from inside B j , respectively. The existence and uniqueness of this Φ j can be established using single layer potentials with suitably chosen densities, see [START_REF] Cedio-Fengya | Identification of conductivity imperfections of small diameter by boundary measurements. Continuous dependence and computional reconstruction[END_REF] for the case of conductivity problem. For each inhomogeneity z j + αB j we introduce the polarizability tensor M j which is a d × d, symmetric, positive definite matrix associated with the j-th inhomogeneity, given by

(M j ) k,l = e k • ( ∂B j (ν j + ( µ j µ 0 -1) ∂Φ j ∂ν j | + (y))y • e l dσ j (y)). (15)
Here (e 1 , . . . , e d ) is an orthonormal basis of R d . In terms of this function we are able to prove the following result about the asymptotic behavior of curl

E α • ν j | ∂(z j +αB j ) + .
Theorem 3.1 Suppose that ( 1) is satisfied and let Φ j , j = 1, . . . , m be given as in [START_REF] Lions | Contrôlabilité exacte, Perturbations et Stabilisation de Systèmes Distribués[END_REF].

Then, for the solutions E α , E of problems ( 3) and ( 4) respectively, and for y ∈ ∂B j we have

( curl E α (z j + αy) • ν j )| ∂(z j +αB j ) + = curl E(z j , t) • ν j (16) +(1 - µ j µ 0 ) ∂Φ j ∂ν j | + (y) • curl E(z j , t) + o(1).
The term o(1) uniform in y ∈ ∂B j and t ∈ (0, T ) and depends on the shape of

{B j } m j=1
and Ω, the constants c 0 , T , µ 0 , {µ j } m j=1 , the data ϕ, ψ, and f , but is otherwise independent of the points {z j } m j=1 .

Proof.

Let E α = curl E α (x, t) and E 0 = curl E(x, t). Then, according to (3)-( 4) we have

ε 0 ∂ 2 t E α -curl 1 µ α E α = 0 and curl E α = 0, for x ∈ Ω. (17)
We restrict, for simplicity, our attention to the case of a single inhomogeneity, i.e., the case m = 1. The proof for any fixed number m of well separated inhomogeneities follows by iteration of the argument that we will present for the case m = 1. In order to further simplify notation, we assume that the single inhomogeneity has the form αB, that is, we assume it is centered at the origin. We denote the electromagnetic permeability inside αB by µ * and define Φ * the same as Φ j , defined in ( 14), but with B j and µ j replaced by B and µ * , respectively. Define ν to be the outward unit normal to ∂B. Now, following a common practice in multiscale expansions we introduce the local variable y = x α , then the domain Ω = ( Ω α ) is well defined.

Next, let be given in

C ∞ 0 (]0, T [). For any function v ∈ L 1 (0, T ; L 2 (Ω)), we define v(x) = T 0 v(x, t) (t) dt ∈ L 2 (Ω).
We remark that

∂ t v(x) = - T 0 v(x, t) (t) dt. So that we deduce from (17) that Êα satisfies      curl 1 µ α Êα = T 0 E α (t) dt in Ω, curl Êα = 0 in Ω. Analogously, Ê satisfies      1 µ 0 curl Ê = T 0 E (t) dt in Ω, curl Ê = 0 in Ω.
Indeed, we have Êα × n = Ê × n = curl ∂Ω f × n on the boundary ∂Ω, where curl ∂Ω is the tangential curl. Following [START_REF] Ammari | Asymptotic formulas for perturbations in the electromagnetic fields due to the presence of inhomogeneities of small diameter II.The full Maxwell equations[END_REF] and [START_REF] Ammari | An inverse initial boundary value problem for the wave equation in the presence of imperfections of small volume[END_REF], we introduce q * α as the unique solution to the following problem

                   ∆q * α = 0 in Ω = ( Ω α ) \ B and in B, q * α is continuous across ∂B, µ 0 ∂q * α ∂ν | + -µ * ∂q * α ∂ν | -= -(µ 0 -µ * ) Ê(αy) • ν on ∂B, q * α = 0 on ∂ Ω.
The jump condition

µ 0 ∂q * α ∂ν | + -µ * ∂q * α ∂ν | -= -(µ 0 -µ * ) Ê(αy) • ν on ∂B
guarantees that Êα (x) -Ê(x) -grad y q * α ( x α ) belongs to the functional space X N (Ω), where grad ∂Ω is the tangential gradient. Since

             curl 1 µ α ( Êα -Ê -grad y q * α ( x α )) = T 0 E α -χ(Ω \ αB)E + µ * µ 0 χ(αB)E (t) dt in Ω, curl ( Êα -Ê -grad y q * α ( x α )) = 0 in Ω, ( Êα -Ê -grad y q * α ( x α )
) × n = 0 on ∂Ω, where χ(ω) is the characteristic function of the domain ω, we arrive, as a consequence of the energy estimate given by Lemma 3.1, at the following

                   ( Êα -Ê -grad y q * α ( x α )) ∈ X N (Ω), curl 1 µ α ( Êα -Ê -curl y q * α ( x α )) = 0(α) in Ω, curl ( Êα -Ê -grad y q * α ( x α )) = 0 in Ω, ( Êα -Ê -grad y q * α ( x α )) × n = 0 on ∂Ω.
From [START_REF] Ammari | Asymptotic formulas for perturbations in the electromagnetic fields due to the presence of inhomogeneities of small diameter II.The full Maxwell equations[END_REF] we know that this yields the following estimate

|| curl 1 µ α ( Êα -Ê -grad y q * α ( x α ))|| L 2 (Ω) + || Êα -Ê -grad y q * α ( x α )|| L 2 (Ω) ≤ Cα,
and so, ( Êα -Ê -grad y q * α (

x α )) • ν| + = 0(α) on ∂(αB).
Now, we denote by q * be the unique (scalar) solution to

                 ∆q * = 0 in R d \ B and in B, q * is continuous across ∂B, µ 0 ∂q * ∂ν | + -µ * ∂q * ∂ν | -= -(µ 0 -µ * ) Ê(0) • ν on ∂B, lim |y|→+∞ q * = 0.
In the spirit of Theorem 1 in [START_REF] Cedio-Fengya | Identification of conductivity imperfections of small diameter by boundary measurements. Continuous dependence and computional reconstruction[END_REF] it follows that ||( grad y q * -grad y q * α )(

x α )|| L 2 (Ω) ≤ Cα 1/2 , which yields ( Êα -Ê -grad y q * ( x α )) • ν = o(1) on ∂(αB).
Writing q * in terms of Φ * gives

T 0 ( curl E α (αy)•ν)| ∂(αB) + -ν• curl E(0, t)-( µ 0 µ * -1) ∂Φ * ∂ν | + (y)• curl E(0, t) (t) dt = o(1),
for any ∈ C ∞ 0 (]0, T [), and so, by iterating the same argument for the case of m (well separated) inhomogeneities z j + αB j , j = 1, . . . , m, we arrive at the promised asymptotic formula [START_REF] Mcdowall | An electromagnetic inverse problem in chiral media[END_REF].

The identification procedure

Before describing our identification procedure, let us introduce the following cutoff function β(x) ∈ C ∞ 0 (Ω) such that β ≡ 1 in a subdomain Ω of Ω that contains the inhomogeneities B α and let η ∈ R d . We will take in what follows E(x, t) = η ⊥ e iη•x-i √ µ 0 |η|t where η ⊥ is a unit vector that is orthogonal to η which corresponds to taking ϕ

(x) = η ⊥ e iη•x , ψ(x) = -i √ µ 0 |η|η ⊥ e iη•x , and f (x, t) = η ⊥ × ne iη•x-i √ µ 0 |η|t
and assume that we are in possession of the measurements of:

curl E α × n on Γ × (0, T ),
where Γ is an open part of ∂Ω. Suppose now that T and the part Γ of the boundary ∂Ω are such that they geometrically control Ω which roughly means that every geometrical optic ray, starting at any point x ∈ Ω at time t = 0 hits Γ before time T at a non diffractive point, see [START_REF] Bardos | Sharp sufficient conditions for the observation, control and stabilization of waves from the boundary[END_REF]. It follows from [START_REF] Nicaise | Exact boundary controllability of Maxwell's equations in heteregeneous media and an application to an inverse source problem[END_REF] (see also [START_REF] Lagnese | Exact boundary controllability of Maxwell's equations in a general region[END_REF], [START_REF] Kime | Boundary controllability of Maxwell's equations in a spherical region[END_REF] and [START_REF] Komornik | Boundary stabilization, observation and control of Maxwell's equations[END_REF]) that there exists (a unique) g η ∈ H 1 0 (0, T ; T L 2 (Γ)) (constructed by the Hilbert Uniqueness Method) such that the unique weak solution

w η to                    (∂ 2 t + curl curl )w η = 0 in Ω × (0, T ), div w η = 0 in Ω × (0, T ), w η | t=0 = β(x)η ⊥ e iη•x , ∂ t w η | t=0 = 0 in Ω, w η × n| ∂Ω\Γ×(0,T ) = 0, w η × n| Γ×(0,T ) = g η , (18) satisfies w η (T ) = ∂ t w η (T ) = 0 in Ω. Let θ η ∈ H 1 (0, T ; T L 2 (Γ)) denote the unique solution of the Volterra equation of second kind    ∂ t θ η (x, t) + T t e -i|η|(s-t) (θ η (x, s) -i|η|∂ t θ η (x, s)) ds = g η (x, t) for x ∈ Γ, t ∈ (0, T ), θ η (x, 0) = 0 for x ∈ Γ. ( 19
)
The existence and uniqueness of this θ η in H 1 (0, T ; T L 2 (Γ)) for any η ∈ R d can be established using the resolvent kernel. However, observing from differentiation of [START_REF] Sylvester | A global uniqueness theorem for an inverse boundary value problem[END_REF] with respect to t that θ η is the unique solution of the ODE:

∂ 2 t θ η -θ η = e i|η|t ∂ t (e -i|η|t g η ) for x ∈ Γ, t ∈ (0, T ), θ η (x, 0) = 0, ∂ t θ η (x, T ) = 0 for x ∈ Γ, ( 20
)
the function θ η may be find (in practice) explicitly with variation of parameters and it also immediately follows from this observation that θ η belongs to H 2 (0, T ; T L 2 (Γ)). We introduce v η as the unique weak solution (obtained by transposition as done in [START_REF] Lions | Nonhomogeneous Boundary Value Problems and Applications[END_REF] and in [START_REF] Lions | Contrôlabilité exacte, Perturbations et Stabilisation de Systèmes Distribués[END_REF] [Theorem 4.2, page 46] for the scalar function) in C 0 (0, T ;

X(Ω)) ∩ C 1 (0, T ; L 2 (Ω)) to the following problem                          (∂ 2 t + curl curl )v η = 0 in Ω × (0, T ), div v η = 0 in Ω × (0, T ), v η | t=0 = 0 in Ω, ∂ t v η | t=0 = m j=1 i(1 - µ 0 µ j )η × (ν j + ( µ 0 µ j -1) ∂Φ j ∂ν j | + )e iη•z j δ ∂(z j +αB j ) ∈ Y (Ω) in Ω, v η × n| ∂Ω×(0,T ) = 0.
Then, the following holds.

Proposition 4.1 Suppose that Γ and T geometrically control Ω. For any η ∈ R d we have

T 0 Γ g η • ( curl v η × n) dσ(x)dt = α 2 m j=1 µ 0 (1 - µ j µ 0 )e 2iη•z j η • ∂B j (ν j (21)
+(

µ j µ 0 -1) ∂Φ j ∂ν j | + (y))η • y ds j (y) + o(α 2 ).
Proof. Multiply the equation (∂ 2 t + curl curl )v η = 0 by w η and integrating by parts over (0, T ) × Ω, for any η ∈ R d we have

α m j=1 i(1 - µ j µ 0 )e 2iη•z j η • ∂B j (ν j + ( µ j µ 0 -1) ∂Φ j ∂ν j | + (y))e iαη•y ds(y) = -µ -1 0 T 0 Γ g η • ( curl v η × n) dσ(x)dt.
Now, we take the Taylor expansion of αe iαη•y in the left side of the last equation, we obtain the convenient asymptotic formula [START_REF] Vogelius | Asymptotic formulas for perturbations in the electromagnetic fields due to the presence of inhomogeneities[END_REF].

To identify the locations and certain properties of the small inhomogeneities B α let us view the averaging of the boundary measurements curl E α × n| Γ×(0,T ) , using the solution θ η to the Volterra equation [START_REF] Sylvester | A global uniqueness theorem for an inverse boundary value problem[END_REF] or equivalently the ODE [START_REF] Somersalo | A linearized inverse boundary value problem for Maxwell's equations[END_REF], as a function of η. The following holds.

Theorem 4.1 Let η ∈ R d . Let E α be the unique solution in C 0 (0, T ; X(Ω))∩C 1 (0, T ; L 2 (Ω)) to the Maxwell's equations (3) with ϕ(x) = η ⊥ e iη•x , ψ(x) = -i √ µ 0 |η|η ⊥ e iη•x , and f (x, t) = η ⊥ e iη•x-i √ µ 0 |η|t
. Suppose that Γ and T geometrically control Ω, then we have

T 0 Γ θ η • ( curl E α × n -curl E × n) + ∂ t θ η • ∂ t ( curl E α × n -curl E × n) dσ(x)dt = α 2 m j=1 (µ 0 -µ j )e 2iη•z j M j (η) • η + O(α 2 ), (22)
where θ η is the unique solution to the Volterra equation [START_REF] Somersalo | A linearized inverse boundary value problem for Maxwell's equations[END_REF] with g η defined as the boundary control in [START_REF] Romanov | Inverse Problems for Maxwell's Equations, Inverse and Ill-posed Problems Series[END_REF] and M j is the polarization tensor of B j , defined by

(M j ) k,l = e k • ( ∂B j (ν j + ( µ j µ 0 -1) ∂Φ j ∂ν j | + (y))y • e l ds j (y)). ( 23
)
Here (e 1 , e 2 ) is an orthonormal basis of R d . The term O(α 2 ) is independent of the points {z j , j = 1, • • • , m}.

Proof. From ∂ t θ η (T ) = 0 and ( curl

E α ×n-curl E ×n)| t=0 = 0 the term T 0 Γ ∂ t θ η • ∂ t ( curl E α × n -curl E × n) dσ(x)dt has to be interpreted as follows T 0 Γ ∂ t θ η • ∂ t ( curl E α × n -curl E × n) = - T 0 Γ ∂ 2 t θ η • ( curl E α × n -curl E × n). (24) Next, introduce Ẽα,η (x, t) = E(x, t) + t 0 e -i √ µ 0 |η|s v η (x, t -s) ds, x ∈ Ω, t ∈ (0, T ). ( 25
)
We have

T 0 Γ θ η • ( curl E α × n -curl E × n) + ∂ t θ η • ∂ t ( curl E α × n -curl E × n) = T 0 Γ θ η • ( curl E α × n -curl Ẽα,η × n) + ∂ t θ η • ∂ t ( curl E α × n -curl Ẽα,η × n) + T 0 Γ θ η • t 0 e -i √ µ 0 |η|s v η (x, t -s) × n ds + ∂ t θ η • ∂ t t 0 e -i √ µ 0 |η|s v η (x, t -s) × n ds .
Since θ η satisfies the Volterra equation [START_REF] Somersalo | A linearized inverse boundary value problem for Maxwell's equations[END_REF] and

∂ t ( t 0 e -i √ µ 0 |η|s v η (x, t -s) × n ds) = ∂ t (-e -i √ µ 0 |η|t t 0 e i √ µ 0 |η|s v η (x, s) × n ds) = i √ µ 0 |η|e -i √ µ 0 |η|t t 0 e i √ µ 0 |η|s v η (x, s) × n ds + v η (x, t) × n,
we obtain by integrating by parts over (0, T ) that

T 0 Γ θ η • t 0 e -i √ µ 0 |η|s v η (x, t -s) × n ds + ∂ t θ η • ∂ t t 0 e -i √ µ 0 |η|s v η (x, t -s) × n ds = T 0 Γ (v η (x, t) × n) • (∂ t θ η + T t θ η (s)e i √ µ 0 |η|(t-s) ds) -i √ µ 0 |η|(e -i √ µ 0 |η|t ∂ t θ η (t)) • t 0 e i √ µ 0 |η|s v η (x, s) × n ds dt = T 0 Γ v η (x, t) × n • (∂ t θ η + T t (θ η (s) -i √ µ 0 |η|∂ t θ η (s))e i √ µ 0 |η|(t-s) ds) dt = T 0 Γ g η (x, t) • ( curl v η (x, t) × n) dt
and so, from Proposition 4.1 we obtain

T 0 Γ θ η • ( curl E α × n -curl E × n) + ∂ t θ η • ∂ t ( curl E α × n -curl E × n) = α 2 m j=1 (1 - µ j µ 0 )e 2iη•z j η • ∂B j (ν j + ( µ j µ 0 -1) ∂Φ j ∂ν j | + (y))η • y ds j (y) + T 0 Γ θ η • ( curl E α × n -curl Ẽα,η × n) + ∂ t θ η • ∂ t ( curl E α × n -curl Ẽα,η × n) + o(α 2 ).
In order to prove Theorem 4.1 it suffices then to show that

T 0 Γ θ η •( curl E α ×n-curl Ẽα,η ×n)+∂ t θ η •∂ t ( curl E α ×n-curl Ẽα,η ×n) = o(α 2 ). (26) Since                              (∂ 2 t -curl 1 µ 0 curl )( t 0 e -i √ µ 0 |η|s v η (x, t -s) ds) = m j=1 i(1 - µ j µ 0 )η × (ν j + ( µ j µ 0 -1) ∂Φ j ∂ν j | + (y))e iη•z j δ ∂(z j +αB j ) e -i √ µ 0 |η|t in Ω × (0, T ), ( t 0 e -i √ µ 0 |η|s v η (x, t -s) ds)| t=0 = 0, ∂ t ( t 0 e -i √ µ 0 |η|s v η (x, t -s) ds)| t=0 = 0 in Ω, ( t 0 e -i √ µ 0 |η|s v η (x, t -s) ds) × n| ∂Ω×(0,T ) = 0, it follows from Theorem 3.1 that            (∂ 2 t -curl 1 µ 0 curl )(E α -Ẽα,η ) = o(α 2 ) in Ω × (0, T ), (E α -Ẽα,η )| t=0 = 0, ∂ t (E α -Ẽα,η )| t=0 = 0 in Ω, (E α -Ẽα,η ) × n| ∂Ω×(0,T ) = 0.
Following the proof of Proposition 3.1, we immediately obtain

||E α -Ẽα,η || L 2 (Ω) = o(α 2 ), t ∈ (0, T ), x ∈ Ω,
where o(α 2 ) is independent of the points {z j } m j=1 . To prove (26) it suffices then from (24) to show that the following estimate holds

|| curl E α × n -curl Ẽα,η × n|| L 2 (0,T ;T L 2 (Γ)) = o(α 2 ). Let θ be given in C ∞ 0 (]0, T [) and define Êα,η (x) = T 0 Ẽα,η (x, t)θ(t) dt and Êα (x) = T 0 E α (x, t)θ(t) dt.
From definition (25) we can write

                 ( Êα -Êα ) ∈ H 1 (Ω), curl curl ( Êα -Êα ) = 0(α) ∈ Y (Ω) in Ω, div ( Êα -Êα ) = 0 in Ω, ( Êα -Êα ) × n| ∂Ω = 0. (27)
In the spirit of the standard elliptic regularity [START_REF] Evans | Partial Differential Equations[END_REF] we deduce for the boundary value problem (27) that

|| curl ( Êα -Êα ) × n|| L 2 (Γ) = O(α 2 ), for all θ ∈ C ∞ 0 (]0, T [); whence || curl (E α -Êα ) × n|| L 2 (Γ) = o(α 2 ) a. e. in t ∈ (0, T ),
and so, the desired estimate [START_REF] Volkov | An inverse problem for the time harmonic Maxwell equations[END_REF] holds. The proof of Theorem 4.1 is then over.

Our identification procedure is deeply based on Theorem 4.1. Let us neglect the asymptotically small remainder in the asymptotic formula [START_REF] Volkov | An inverse problem for the time harmonic Maxwell equations[END_REF], and define ℵ α (η) by

ℵ α (η) = T 0 Γ θ η • ( curl (E α -E) × n) + ∂ t θ η • ∂ t ( curl (E α -E) × n) .
Recall that the function e 2iη•z j is exactly the Fourier Transform (up to a multiplicative constant) of the Dirac function δ -2z j (a point mass located at -2z j ). From Theorem 4.1 it follows that the function e 2iη•z j is (approximately) the Fourier Transform of a linear combination of derivatives of point masses, or

אα (η) ≈ α 2 m j=1 L j δ -2z j ,
where L j is a second order constant coefficient, differential operator whose coefficients depend on the polarization tensor M j defined by (23) (see [START_REF] Cedio-Fengya | Identification of conductivity imperfections of small diameter by boundary measurements. Continuous dependence and computional reconstruction[END_REF] for its properties) and אα (η) represents the inverse Fourier Transform of ℵ α (η). The reader is referred to [START_REF] Cedio-Fengya | Identification of conductivity imperfections of small diameter by boundary measurements. Continuous dependence and computional reconstruction[END_REF] for properties of the tensor polarization M j . The method of reconstruction consists in sampling values of אα (η) at some discrete set of points and then calculating the corresponding discrete inverse Fourier Transform. After a rescaling the support of this discrete inverse Fourier Transform yields the location of the small inhomogeneities B α . Once the locations are known we may calculate the polarization tensors (M j ) m j=1 by solving an appropriate linear system arising from [START_REF] Volkov | An inverse problem for the time harmonic Maxwell equations[END_REF]. This procedure generalizes the approach developed in [START_REF] Ammari | Boundary integral formulas for the reconstruction of electromagnetic imperfections of small diameter[END_REF] for the two-dimensional (time-independent) inverse conductivity problem and generalize the results in [START_REF] Ammari | An inverse initial boundary value problem for the wave equation in the presence of imperfections of small volume[END_REF] to the full time-dependent Maxwell's equations.

Conclusion

In this paper, we are convinced that the use of approximate formulae such as [START_REF] Volkov | An inverse problem for the time harmonic Maxwell equations[END_REF] represents a very promising approach to the dynamical identification of small inhomogeneities that are embedded in a homogeneous medium. We also believe that our method yields a good approximation to small amplitude perturbations in the electromagnetic parameters (for the example of electric permittivity ε α (x) = ε 0 + αε(x)) from the measurements: curl H α × n on Γ × (0, T ).

Our method may yield the Fourier transform of the amplitude perturbation ε(x). This issue will be considered in a forthcoming work [START_REF] Daveau | An inverse problem for the time dependent Maxwell's equations in the presence of imperfections of small volume[END_REF].