INSTITUT FÜR WISSENSCHAFTLICHES RECHNEN UND MATHEMATISCHE MODELLBILDUNG Resolution-Controlled Conductivity Discretization in Electrical Impedance Tomography
Abstract
This work contributes to the numerical solution of the inverse problem of determining an isotropic conductivity from boundary measurements, known as Electrical Impedance Tomography. To this end, we first investigate the imaging resolution of the Complete Electrode Model in a circular geometry using analytic solutions of the forward problem and conformal maps. Based on this information we propose a discretization of the conductivity space. Roughly speaking, the resulting conductivity meshes comply with the maximal resolution provided by discrete data with a known noise level. We heuristically extend this approach to domains of arbitrary shape and present its performance under a Newton-type inversion algorithm. 1. Background and motivation Electrical impedance tomography (EIT) is an imaging method for determining the electrical conductivity of an object from measurements on its surface. The underlying mathematical model is an elliptic boundary value problem where the Dirichlet data represents the potential on the object surface, the Neumann data represents the normal current flow through the surface and an elliptic PDE models the flow of current inside the object. The coefficient of the PDE is the searched-for conductivity. The evaluation of the Neumann-to-Dirichlet (ND) operator for a given conductivity is called the forward problem of EIT. The task in EIT is usually to solve the inverse conductivity problem (ICP) of determining the conductivity from the knowledge of the ND operator. Electrode models and inversion. The pioneering work of Calderón [Cal80] started the investigation of the continuum boundary model, where a complete knowledge of the continuum ND operator is assumed. However in practice, one can only apply currents and measure potentials through finitely many electrodes. As a consequence, the available boundary data is finitedimensional – the ND map of this discrete model is a matrix. Moreover, due to electrochemical effects and the conducting nature of electrodes, it is impossible to prescribe Neumann data or access Dirichlet data explicitly anywhere at the boundary. Instead, the accurate Complete Electrode Model (CEM) [SCI92] described in section 2.2 uses Robin-type boundary conditions to account for these effects. Analytic properties of electrode models are rare. [LR08] show the injectivity of the Fréchet derivative of the ND map for piecewise polynomial conductivities on a triangulation of the domain for the CEM. The number of electrodes necessary for injectivity is finite, but unknown for any fixed triangulation. [BDGVM12] use a model reduction approach to transform the ICP into the problem of determining the resistors in a network, which is uniquely solvable under certain conditions. However, the reduction approach is theoretically justified for radially symmetric near-constant conductivities only and the boundary model is slightly different from commonly used electrode models. Stability of the reconstruction scheme is a critical matter for solving the ICP with measured data. Due to Alessandrini’s logarithmic stability estimate [Ale88, Theorem 1] in the operator norm of the Dirichlet-to-Neumann (DN) map, the ICP for the continuum model is usually 2010 Mathematics Subject Classification. 65N20, 65F22, 92C55. The work of the authors was supported by the German Research Foundation (DFG) under grant RI 975/7-1. 1 2 ROBERT WINKLER AND ANDREAS RIEDER said to be exponentially unstable. Later results (e.g. [Dob92], [Pal02], [MMM04]) show that the stability decays rapidly away from the boundary in typical norms. The consequence of all results is that regularization is necessary to achieve a stable inversion scheme for the ICP. In the absence of a general analytic result, stability studies for electrode models are usually done by investigating what can actually be detected in a fixed domain geometry and electrode setting. Here the notion of sensitivity (or distinguishability [Isa86], detectability), i.e. the effect of perturbations in conductivity to boundary measurements, is important. Most algorithms for the computation of the ICP are either motivated by the continuum model (for example complex geometrical optics solutions [Nac96], factorization method [BH00]) or regularized iterative methods for an electrode model. An overview can be found in [MS12]. In this work we consider the iterative Newton-type method REGINN for the CEM [LR06]. All numerical inversion methods require a discretization of the problem, in particular of the underlying conductivity space. Aim and structure of this work. We have pointed out that it is necessary to discretize the conductivity space and to find a suitable regularization scheme when solving the ICP for an electrode model and noisy measurements. While the former is traditionally done by a simple triangulation of the domain, the choice of regularization is delicate as it imposes assumptions on the searched-for conductivity – usually on its smoothness – to control noise amplification during inversion. As a consequence, EIT images are typically heavily blurred. Our tailored discretization scheme for the conductivity space is piecewise constant on a partition of the domain with a locally adaptive, noise-dependent mesh size which we call an optimal resolution mesh. It inherits the smallest possible details which are recoverable from boundary data with a given noise level. Solving the ICP on these meshes with a Newton-type algorithm can be done efficiently with the estimated noise level as the only free design parameter for regularization. In section 2, we introduce the relevant EIT models and notations. In section 3, we derive an analytic solution for the CEM forward problem on a disk in the presence of a centered circular perturbation in conductivity. Next, we use conformal maps to extend this result to perturbations at arbitrary locations in section 4. We use these analytic solutions to get local resolution information of the CEM. With this information, we design optimal resolution meshes to discretize the conductivity space in section 5. We outline the connections between the CEM and the continuum model and also compare our results to the results of [MMM04] and the meshes arising from resistor networks in [BDGVM12]. We then heuristically extend our results to domains of arbitrary shape. Finally, we demonstrate the performance of these discretizations using the regularized inexact Newton scheme REGINN [Rie99] in section 6. 2. Preliminaries A potential u ∈ H1(Ω) on a source-free, simply connected domain Ω ⊂ R2 with piecewise Lipschitz boundary is a solution of the elliptic equation −∇ · (σ∇u) = 0 on Ω, (1) where σ∈L+ (Ω)={φ∈L∞(Ω)|φ > 0 a.e.} is the isotropic conductivity coefficient. 2.1. Continuum boundary model. For given Dirichlet boundary data u|∂Ω = f ∈ H(∂Ω), (2) the problem (1),(2) has a unique solution. Denoting by ν the outer normal on ∂Ω, we call σ ∂u ∂ν =: iν ∈ H (∂Ω) = { φ ∈ H−1/2(∂Ω) |〈φ, 1〉∂Ω = 0 } (3) RESOLUTION-CONTROLLED DISCRETIZATION IN EIT 3 the corresponding Neumann data. Conversely, the problem (1),(3) has a solution which is unique if we require the trace of u to have vanishing mean, that is u|∂Ω ∈ H 1/2 (∂Ω) = { φ ∈ H(∂Ω) |〈φ, 1〉∂Ω = 0 } . Hence, the ND operator Λσ : H −1/2 (∂Ω)→ H (∂Ω), iν 7→ f, and its inverse are well-defined and one-to-one. The ICP in this setting is the inversion of the operator σ 7→ Λσ. Uniqueness is shown in [AP06]. 2.2. Electrode boundary models. Electrode boundary models describe the injection of currents and the measurement of potentials through a finite number of electrodes. Assume we have electrodes E1, . . . , EL ⊂ ∂Ω, L ∈ N≥2, at which we can inject current patterns I ∈ RL = { V ∈ RL |〈V, 1〉 = 0 } and measure the resulting potential vectors U ∈ RL . We assume that normal current iν at the domain boundary only occurs at electrodes and that the electrodes are perfect conductors, i.e. the potential is constant on each electrode. Furthermore, we model contact impedances zl ≥ 0 that occur at each electrode-domain interface and cause a potential drop zliν , cf. [SCI92]. This leads to Robintype boundary conditions iν = 0 on ∂Ω \ E, E = E1 ∪ . . . ∪ EL, (4a) u+ zliν = Ul on El, l = 1, . . . , L, (4b) ∫ El iν dS = Il, l = 1, . . . , L. (4c) For zl ≡ 0, (1),(4) is called the Shunt Model, for zl > 0 it is called Complete Electrode Model (CEM). The CEM is usually defined weakly as the unique solution (u, U)∈H1 (Ω)× RL of a( (u, U), (w,W ) ) = L ∑ l=1 IlWl for all (w,W ) ∈ H(Ω)× R , (5) where a : ( H1(Ω)× RL ) × ( H1(Ω)× RL ) → R, a ( (v,W ), (w,W ) ) = ∫