Generalized linearization techniques in electrical impedance tomography

TL;DR

This paper compares various linearization methods in electrical impedance tomography, introduces a new technique based on the logarithm of the Neumann-to-Dirichlet operator, and demonstrates its improved accuracy through numerical examples.

Contribution

It proposes a novel linearization approach using the logarithm of the Neumann-to-Dirichlet operator, enhancing accuracy in electrical impedance tomography reconstructions.

Findings

Conductivity linearization often reduces accuracy.

Logarithmic linearization improves forward and inverse problem solutions.

Numerical results confirm the effectiveness of the new method.

Abstract

Electrical impedance tomography aims at reconstructing the interior electrical conductivity from surface measurements of currents and voltages. As the current-voltage pairs depend nonlinearly on the conductivity, impedance tomography leads to a nonlinear inverse problem. Often, the forward problem is linearized with respect to the conductivity and the resulting linear inverse problem is regarded as a subproblem in an iterative algorithm or as a simple reconstruction method as such. In this paper, we compare this basic linearization approach to linearizations with respect to the resistivity or the logarithm of the conductivity. It is numerically demonstrated that the conductivity linearization often results in compromised accuracy in both forward and inverse computations. Inspired by these observations, we present and analyze a new linearization technique which is based on the logarithm of the Neumann-to-Dirichlet operator. The method is directly applicable to discrete settings, including the complete electrode model. We also consider Fr\'echet derivatives of the logarithmic operators. Numerical examples indicate that the proposed method is an accurate way of linearizing the problem of electrical impedance tomography.

Tables

Table 1. Table 1: Parameters defining the discrete random conductivity fields F1–F4 and the random contact conductances C1–C2 used in the numerical experiments.

	$N$	${\hat{\kappa }}_0$	${\varsigma }^2$	$\ell$	$𝐄\left[{\parallel \kappa \parallel }_{L^2\left(\Omega \right)}\right]$ \bigstrut[b]
F1	$1800$	$0$	$1/4$	$1/3$	$1.53$ \bigstrut[t]
F2	$1800$	$0$	$1/4$	$2/3$	$1.48$
F3	$1800$	$0$	$1$	$1/3$	$3.08$
F4	$1800$	$0$	$1$	$2/3$	$2.95$

Table 2. Table 2: Linearization errors e 𝑒 e for the continuum forward model (left) and for the CEM (right) with the different random distributions defined in Table  1 . The forward operators parametrized with respect to conductivity, resistivity and log-conductivity are denoted by Λ ^ ^ Λ \hat{\varLambda} , Λ ^ inv subscript ^ Λ inv \hat{\varLambda}_{\rm inv} and Λ ^ exp subscript ^ Λ \hat{\varLambda}_{\exp} , respectively, whereas the logarithmic forward operator is L ^ ^ 𝐿 \hat{L} .

	$e\left(\hat{\Lambda }\right)$	$e\left({\hat{\Lambda }}_{\mathrm{inv}}\right)$	$e\left({\hat{\Lambda }}_{\exp }\right)$	$e\left(\hat{L}\right)$ \bigstrut[b]
F1	$0.237$	$0.076$	$0.095$	$0.039$ \bigstrut[t]
F2	$0.294$	$0.043$	$0.123$	$0.027$
F3	$0.989$	$0.288$	$0.327$	$0.142$
F4	$1.481$	$0.146$	$0.434$	$0.094$

Table 3. Table 3: Conductivity reconstruction errors ι 𝜄 \iota for the linearized continuum model (left) and for the linearized CEM (right), together with the conductance errors d 𝑑 d for the CEM. The considered random fields and parameters are described in Table  1 .

	$\iota \left({\hat{\Lambda }}_{\exp }\right)$	$\iota \left(\hat{L}\right)$ \bigstrut[b]
F1	$0.583$	$0.362$ \bigstrut[t]
F2	$0.524$	$0.103$
F3	$2.134$	$1.028$
F4	$2.000$	$0.283$

Equations

$$\nabla \cdot (\sigma \nabla u)$$

$$\sigma \frac{\partial u}{\partial \nu }$$

$$L_{+}^{\infty }(\Omega ):=\left\{v\in L^{\infty }(\Omega ;\mathbf{R}):v\ge a\text{ a.e.leavevmodenobreak in }\Omega \text{ for some }a>0\right\},$$

$$H_{\diamond }^r(\partial \Omega ):=\left\{v\in H^r(\partial \Omega ):⟨1,v⟩=0\right\},r\in (-1,1),$$

$$\int_{\varOmega}\sigma\nabla u\cdot\nabla\overline{v}\,{\rm d}x=\big{\langle}f,v|_{\partial\varOmega}\big{\rangle}$$

$$U:=u{|}_{\partial \Omega }\in H_{\diamond }^{1/2}(\partial \Omega )$$

$$\Lambda (\sigma ):H_{\diamond }^{-1/2}(\partial \Omega )\to H_{\diamond }^{1/2}(\partial \Omega )$$

$$⟨f,\Lambda (\sigma )g⟩=\overline{⟨g,\Lambda (\sigma )f⟩}\mathrm{for}\mathrm{all}f,g\in H_{\diamond }^{-1/2}(\partial \Omega )$$

$$⟨f,\Lambda (\sigma )f⟩\ge c\Vert f{\Vert }_{H^{-1/2}(\partial \Omega )}^2\mathrm{for}\mathrm{all}f\in H_{\diamond }^{-1/2}(\partial \Omega )\mathrm{and}\mathrm{some}c>0,$$

$$\Lambda (\sigma ):f\mapsto U,L_{\diamond }^2(\partial \Omega )\to L_{\diamond }^2(\partial \Omega ),$$

$$\Lambda (\sigma )f=\sum _{k=1}^{\infty }{\lambda }_k(f,{\varphi }_k){\varphi }_k,$$

$$\Lambda :\sigma \mapsto \Lambda (\sigma ),L_{+}^{\infty }(\Omega )\to \mathcal{L}(L_{\diamond }^2(\partial \Omega ))$$

$${\Lambda }_{\mathrm{inv}}:=\Lambda \circ \mathrm{inv}:\rho \mapsto \Lambda ({\rho }^{-1}),L_{+}^{\infty }(\Omega )\to \mathcal{L}(L_{\diamond }^2(\partial \Omega )).$$

$${\Lambda }_{\exp }:=\Lambda \circ \exp :\kappa \mapsto \Lambda ({\mathrm{e}}^{\kappa }),L^{\infty }(\Omega )\to \mathcal{L}(L_{\diamond }^2(\partial \Omega )).$$

$$D\Lambda (\sigma ;\eta ):L_{\diamond }^2(\partial \Omega )\to L_{\diamond }^2(\partial \Omega )$$

$$\big{(}f,D\varLambda(\sigma;\eta)g\big{)}=-\int_{\varOmega}\eta\nabla u_{f}\cdot\nabla\overline{u}_{g}\,{\rm d}x,\qquad f,g\in L^{2}_{\diamond}(\partial\varOmega),$$

$$\big{(}f,D\varLambda_{\rm inv}(\rho;\eta)g\big{)}=\int_{\varOmega}\frac{\eta}{\rho^{2}}\nabla u_{f}\cdot\nabla\overline{u}_{g}\,{\rm d}x$$

$$\big{(}f,D\varLambda_{\exp}(\kappa;\eta)g\big{)}=-\int_{\varOmega}\eta\,\mathrm{e}^{\kappa}\nabla u_{f}\cdot\nabla\overline{u}_{g}\,{\rm d}x$$

$$\log \Lambda (\sigma ):f\mapsto \sum _{k=1}^{\infty }\log ({\lambda }_k)(f,{\varphi }_k){\varphi }_k,$$

$$\mathcal{D}\big{(}\log\!\varLambda(\sigma)\big{)}=\bigg{\{}g\in L^{2}_{\diamond}(\partial\varOmega)\ :\ \|\log\!\varLambda(\sigma)g\|_{L^{2}(\partial\varOmega)}^{2}=\sum_{k=1}^{\infty}\log^{2}(\lambda_{k})|(g,\phi_{k})|^{2}<\infty\bigg{\}}.$$

$$\log \Lambda (\sigma {)}^{\ast }:\mathcal{D}(\log \Lambda (\sigma {)}^{\ast })\to L_{\diamond }^2(\partial \Omega )$$

$$(\log \Lambda (\sigma )f,g)=(f,\log \Lambda (\sigma )g)$$

$$(\log \Lambda (\sigma )f,g)=(f,g^{\ast })\mathrm{for}\mathrm{all}f\in \mathcal{D}(\log \Lambda (\sigma )).$$

$$(f,g^{\ast })=(\log \Lambda (\sigma )f,g_0)+(\log \Lambda (\sigma )f,g_1)=(\log \Lambda (\sigma )f,g_1)$$

$$\Vert g{\Vert }_{\mathcal{G}(\log \Lambda (\sigma ))}={\left(\Vert g{\Vert }_{L^2(\partial \Omega )}^2+\Vert \log \Lambda (\sigma )g{\Vert }_{L^2(\partial \Omega )}^2\right)}^{1/2},$$

$$\log \Lambda (\sigma ):H_{\diamond }^{ϵ}(\partial \Omega )\to L_{\diamond }^2(\partial \Omega ),ϵ>0,$$

$$H^{s}_{\diamond}(\partial\varOmega)=\big{\{}g\in L^{2}_{\diamond}(\partial\varOmega)\ :\ \|g\|_{{s}}<\infty\big{\}},\qquad 0\leq s\leq\frac{1}{2},$$

$$\|g\|_{s}:=\bigg{(}\sum_{k=1}^{\infty}\frac{1}{\lambda_{k}^{2s}}|(g,\phi_{k})|^{2}\bigg{)}^{1/2}$$

$${\Lambda }^s(\sigma ):f\mapsto \sum _{k=1}^{\infty }{\lambda }_k^s(f,{\varphi }_k){\varphi }_k,L_{\diamond }^2(\partial \Omega )\to L_{\diamond }^2(\partial \Omega ),$$

$$\varLambda^{-s}(\sigma)\colon f\mapsto\sum_{k=1}^{\infty}\frac{1}{\lambda_{k}^{s}}(f,\phi_{k})\,\phi_{k},\quad\mathcal{D}\big{(}\varLambda^{-s}(\sigma)\big{)}\to L^{2}_{\diamond}(\partial\varOmega),\qquad s>0,$$

$$\mathcal{D}\big{(}\varLambda^{-s}(\sigma)\big{)}:=\mathcal{R}\big{(}\varLambda^{s}(\sigma)\big{)}=\big{\{}g\in L^{2}_{\diamond}(\partial\varOmega)\ :\ \|\varLambda^{-s}(\sigma)g\|_{L^{2}(\partial\varOmega)}=\|g\|_{s}<\infty\big{\}}.$$

Full text

Keywords: electrical impedance tomography, inverse elliptic boundary value problems, linearization, Neumann-to-Dirichlet map, operator logarithm, Fréchet differentiability

AMS subject classifications: 65N21, 47B15, 35R30

1 Introduction

The reconstruction task of electrical impedance tomography (EIT) is undoubtedly one of the most studied nonlinear inverse boundary value problems [22]. The mathematical foundations of the problem were introduced by Calderón in [3], which considers the conductivity reconstruction based on the Dirichlet-to-Neumann operator. As has become a more standard numerical approach, in this paper we consider the Neumann-to-Dirichlet operator, which maps the boundary currents to the corresponding boundary potentials for a given interior conductivity. The dependence of the Neumann-to-Dirichlet operator on the conductivity is nonlinear so that reconstructing the conductivity from the measurements is a nonlinear, and also illposed, inverse problem.

Nonlinear inverse problems can be straightforwardly approached with nonlinear least squares minimization algorithms, which rely on successive linearizations of the forward operator. In particular, the accuracy of the sequential linearizations certainly affects the performance of such an iterative reconstruction method. In EIT, one-step linearization can also be used to obtain an approximative reconstruction [1, 4, 11]. The Fréchet derivatives of the EIT forward operator with respect to the conductivity can be computed explicitly [4, 14, 18]. By applying the chain rule of Banach spaces, or via direct computation, it is possible to consider derivatives with respect to the resistivity or the logarithm of the conductivity as well. In this paper, we compare the linearization errors resulting from these different input parametrizations of the electrical properties. We also discuss the corresponding methods for the complete electrode model (CEM), which is an accurate model for practical EIT measurements [5, 21].

As the main novelty, we introduce the logarithm of the Neumann-to-Dirichlet operator and use its differentiability properties to construct a new linearization method for EIT. Loosely speaking, the traditional forward operator is replaced by an operator that maps the logarithm of the conductivity to the logarithm of the boundary potentials, the latter being understood in a sense of linear operators. The corresponding least squares inversion algorithm then fits the computed logarithmic potentials to the logarithm of the measurement matrix. Although this logarithmic forward operator is still nonlinear, numerical experiments show that the resulting linearization errors are in most cases smaller than with any other considered linearization method.

This paper is organized as follows. In Section 2, the continuum forward model of EIT and the related Neumann-to-Dirichlet operator are reviewed. We write the dependence of the Neumann-to-Dirichlet operator on the electrical properties in three different ways and recall also the corresponding Fréchet derivatives. Section 3 presents the formal definition for the logarithm of the Neumann-to-Dirichlet operator and studies its differentiability and other properties as an unbounded operator on the space of square-integrable functions. In Section 4, observations from Sections 2 and 3 are generalized for the complete electrode model. Numerical experiments concerning linearization errors for both forward and inverse computations are given in Section 5. Finally, Section 6 offers concluding remarks.

2 Input parametrization in EIT

The continuum forward model of EIT for the electric potential $u$ is written as

[TABLE]

where $\varOmega\subset\mathbf{R}^{d}$, $d\geq 2$, is a bounded domain with a Lipschitz boundary $\partial\varOmega$ and a connected complement, the electrical conductivity $\sigma\in L_{+}^{\infty}(\varOmega)$ is real-valued and isotropic, and $f\in H_{\diamond}^{-1/2}(\partial\varOmega)$ can be complex-valued. The conductivity is bounded from below by a positive constant, that is,

[TABLE]

and the boundary current density belongs to the mean-free Sobolev space $H_{\diamond}^{-1/2}(\partial\varOmega)$ defined via

[TABLE]

due to the conservation of electric charge. Here and in what follows, the bracket $\langle\cdot,\cdot\rangle\colon H^{-r}(\partial\varOmega)\times H^{r}(\partial\varOmega)\to\mathbf{C}$ denotes the sesquilinear dual pairing that has an interpretation as an extension of the $L^{2}(\partial\varOmega)$ inner product $(\cdot,\cdot)\colon L^{2}(\partial\varOmega)\times L^{2}(\partial\varOmega)\to\mathbf{C}$. In particular, apart from $L^{\infty}(\varOmega):=L^{\infty}(\varOmega;\mathbf{R})$, the multiplier field for all employed function spaces is $\mathbf{C}$.

The variational form of the Neumann problem (1) is to find $u\in H^{1}(\varOmega)$ such that

[TABLE]

holds for all $v\in H^{1}(\varOmega)$. The standard theory for elliptic partial differential equations states that there exists a unique solution for (5) in the quotient space $H^{1}(\varOmega)/\mathbf{C}$ for any given current density $f\in H^{-1/2}_{\diamond}(\varOmega)$. In particular, there is a unique mean-free boundary potential

[TABLE]

that depends linearly and boundedly on the corresponding $f\in H^{-1/2}_{\diamond}(\varOmega)$. [9]

The linear map $f\mapsto U$, which obviously depends on $\sigma$, is called the Neumann-to-Dirichlet operator and is denoted by $\varLambda(\sigma)$. For any given $\sigma\in L_{+}^{\infty}(\varOmega)$, the mapping

[TABLE]

is a linear isomorphism that is self-adjoint in the sense that

[TABLE]

and positive,

[TABLE]

as can be easily deduced from (5) and (Neumann) trace theorems for those elements of $H^{1}(\varOmega)/\mathbf{C}$ for which the range of $\nabla\cdot\sigma\nabla(\cdot)$ is a subspace of $L^{2}(\varOmega)$ (cf., e.g., [7, p. 381, Lemma 1]). Unless stated otherwise, we interpret $\varLambda(\sigma)$ as a self-adjoint operator from $L^{2}_{\diamond}(\partial\varOmega)$ to itself, that is, as a map

[TABLE]

which is compact as are the (dense) embeddings $H_{\diamond}^{1/2}(\partial\varOmega)\hookrightarrow L^{2}_{\diamond}(\partial\varOmega)\hookrightarrow H_{\diamond}^{-1/2}(\partial\varOmega)$. In particular, $\varLambda(\sigma)$ admits a spectral decomposition

[TABLE]

where the eigenvalues satisfy $\lambda_{k}\geq\lambda_{k+1}$ and $\mathbf{R}_{+}\ni\lambda_{k}\to 0$ as $k\to\infty$, and the corresponding eigenfunctions $\{\phi_{k}\}_{k=1}^{\infty}$ form an orthonormal basis for $L_{\diamond}^{2}(\partial\varOmega)$.

Let us denote by $\mathcal{L}(L_{\diamond}^{2}(\partial\varOmega))$ the Banach space of bounded linear operators from $L_{\diamond}^{2}(\partial\varOmega)$ to itself. The mapping

[TABLE]

is nonlinear and is called the (continuum model) forward operator of EIT. The inverse problem of EIT is to find $\sigma$ from the knowledge of $\varLambda(\sigma)$. In practice, this means that the determination of the conductivity is based on current-voltage pairs measured on the boundary. The following example shows that by composing the forward operator with an elementary function, one can obtain another forward operator that in some cases depends more linearly on the electrical input parameters.

Example 1** (Constant conductivity).**

If $\sigma\equiv\operatorname*{inv}(\rho)\coloneqq\rho^{-1}$ for $\rho\in\mathbf{R}_{+}$, then the potential $u$ from (1) solves the Laplace equation $\Delta u=0$ with the normal derivative $\rho f$ on the boundary. In particular, the mapping $\rho\mapsto\varLambda(\sigma)=(\varLambda\circ\operatorname*{inv})(\rho)$ is linear when restricted to spatially constant functions.

Of course, the example above is nothing but the Ohm’s law, stating that in a homogeneous medium the potential depends linearly on the resistivity $\rho$. Motivated by the example, we denote

[TABLE]

Similarly, the logarithm of the conductivity is denoted by $\kappa\coloneqq\log(\sigma)$ and the corresponding composite forward operator is defined as

[TABLE]

The advantage of considering the logarithm of the conductivity is that the domain of $\varLambda_{\exp}$ is the natural $L^{\infty}$-space without positivity constraints, simplifying many optimization schemes in numerical computations. On the other hand, the advantage of using the standard conductivity parametrization with the operator $\varLambda$ is the simplicity of the corresponding sesquilinear forms, for example in the context of stochastic Galerkin methods [19]. Note that considering the logarithm of the resistivity would not add anything new compared to (12), except for a sign change.

It is well known that the map $\sigma\mapsto\varLambda(\sigma)$ is Fréchet differentiable, i.e., that the bounded derivative

[TABLE]

exists for every $\sigma\in L^{\infty}_{+}(\varOmega)$ and depends linearly and boundedly on $\eta\in L^{\infty}(\varOmega)$ in the topology of $\mathcal{L}(L_{\diamond}^{2}(\partial\varOmega))$ (see, e.g., [18, Section 3]). In fact, the derivative can be obtained from the equation

[TABLE]

where $u_{f}\in H^{1}(\varOmega)$ denotes the solution to the forward problem (1) with the conductivity $\sigma$ and current density $f$. By using the chain rule of differentiation for Banach spaces and (13), it easily follows that the alternative parametrizations of the forward operator (11) and (12) are also Fréchet differentiable and the corresponding derivatives can be characterized by

[TABLE]

and

[TABLE]

for $f,g\in L^{2}_{\diamond}(\partial\varOmega)$.

3 Logarithmic forward operator

In this section, we first introduce the logarithm of the Neumann-to-Dirichlet map $\varLambda(\sigma)$ as an unbounded operator on $L^{2}_{\diamond}(\partial\varOmega)$ and subsequently consider its differentiability with respect to the conductivity.

3.1 Formal definition

The spectral representation (9) allows a simple way of defining a logarithm for the Neumann-to-Dirichlet operator:

[TABLE]

where $\log$ denotes the principal branch of the natural logarithm. As the eigenvalues $\{\lambda_{k}\}_{k=1}^{\infty}\subset\mathbf{R}_{+}$ accumulate at zero, $\log\!\varLambda(\sigma)$ is not bounded as an operator from $L^{2}_{\diamond}(\partial\varOmega)$ to itself. We define the domain of $\log\!\varLambda(\sigma)$ to simply be

[TABLE]

It is obvious that $\mathcal{D}(\log\!\varLambda(\sigma))$ is a dense linear subspace of $L^{2}_{\diamond}(\partial\varOmega)$ and $\log\!\varLambda(\sigma)f\in L^{2}_{\diamond}(\partial\varOmega)$ for any $f\in\mathcal{D}(\log\!\varLambda(\sigma))$.

Proposition 1**.**

The logarithmic Neumann-to-Dirichlet operator $\log\!\varLambda(\sigma)$ defined by (16) and (17) is a self-adjoint (unbounded) operator on $L^{2}_{\diamond}(\partial\varOmega)$.

Proof.

Since $\mathcal{D}(\log\!\varLambda(\sigma))$ is dense in $L^{2}_{\diamond}(\partial\varOmega)$, the adjoint operator

[TABLE]

is well defined [23, p. 196, Thm. 1]. With the help of the Cauchy–Schwarz inequality, it is easy to check that

[TABLE]

is finite for all $f,g\in\mathcal{D}(\log\!\varLambda(\sigma))$. In other words, $\log\!\varLambda(\sigma)$ is symmetric, that is, $\log\!\varLambda(\sigma)^{*}=\log\!\varLambda(\sigma)$ on $\mathcal{D}(\log\!\varLambda(\sigma))\subset\mathcal{D}(\log\!\varLambda(\sigma)^{*})$.

Let $V:=\mathcal{N}(\log\!\varLambda(\sigma))\subset\mathcal{D}(\log\!\varLambda(\sigma))$ be the nullspace of $\log\!\varLambda(\sigma)$. Note that $V$ is either finite-dimensional or the trivial subspace, depending on whether $\lambda=1$ is an eigenvalue of $\varLambda(\sigma)$. In consequence, $L^{2}_{\diamond}(\partial\varOmega)=V\oplus W$, where $W:=V^{\perp}$ is the orthogonal complement of $V$. It follows from (16) and (17) that the range of $\log\!\varLambda(\sigma)$ is exactly $W$, and thus $\log\!\varLambda(\sigma)$ is self-adjoint as an operator from $\mathcal{D}(\log\!\varLambda(\sigma))\cap W$ to $W$ [23, p. 199, Corollary].

Let $g\in\mathcal{D}(\log\!\varLambda(\sigma)^{*})$ be arbitrary, i.e., there exists $g^{*}\in L^{2}_{\diamond}(\partial\varOmega)$ such that

[TABLE]

Decompose $g=g_{0}+g_{1}$, with $g_{0}\in V\subset\mathcal{D}(\log\!\varLambda(\sigma))$ and $g_{1}\in W$, and observe that

[TABLE]

for all $f\in\mathcal{D}(\log\!\varLambda(\sigma))$ due to (18). Since $\log\!\varLambda(\sigma)$ is self-adjoint on $\mathcal{D}(\log\!\varLambda(\sigma))\cap W$, $g_{1}$ must belong to $\mathcal{D}(\log\!\varLambda(\sigma))\cap W$, and altogether we have $g=g_{0}+g_{1}\in\mathcal{D}(\log\!\varLambda(\sigma))\subset\mathcal{D}(\log\!\varLambda(\sigma)^{*})$. In other words, $\mathcal{D}(\log\!\varLambda(\sigma)^{*})=\mathcal{D}(\log\!\varLambda(\sigma))$, which completes the proof. ∎

Because any self-adjoint operator is also closed [23], $\mathcal{D}(\log\!\varLambda(\sigma))$ becomes a Hilbert space if equipped with the graph norm

[TABLE]

with respect to which $\log\!\varLambda(\sigma)$ is trivially a bounded operator.

Corollary 1**.**

The logarithmic Neumann-to-Dirichlet map $\log\!\varLambda(\sigma)$ defined by (16) and (17) can be interpreted as a compact operator

[TABLE]

that coincides with its $L^{2}_{\diamond}(\partial\varOmega)$-adjoint on $H^{\epsilon}_{\diamond}(\partial\varOmega)$.

The proof of Corollary 1 is based on the following lemma. Observe that the lemma could be extended (with obvious modifications) for $-1/2\leq s\leq 1/2$ by duality and for a wider scale of smoothness indices by utilizing some integer power of $\varLambda(\sigma)^{-1}$ in place of $\varLambda(\sigma)^{-1/2}$ in the following proof, assuming $\partial\varOmega$ is smooth enough.

Lemma 1**.**

For any fixed $\sigma\in L^{\infty}_{+}(\varOmega)$, it holds that

[TABLE]

where

[TABLE]

defines an equivalent norm for $H^{s}_{\diamond}(\partial\varOmega)$.

Proof.

Let us introduce the positive powers of $\varLambda(\sigma)\colon L^{2}_{\diamond}(\partial\varOmega)\to L^{2}_{\diamond}(\partial\varOmega)$ in the natural way, that is,

[TABLE]

which defines a compact, injective, self-adjoint operator with a dense range for any $s>0$. The negative powers are the corresponding inverse operators:

[TABLE]

where (cf., e.g., [8, Theorem 2.8])

[TABLE]

Using the same arguments as for $\log\!\varLambda(\sigma)$ above, it is trivial to deduce that $\varLambda^{-s}(\sigma)\colon\mathcal{D}(\varLambda^{-s}(\sigma))\to L^{2}_{\diamond}(\partial\varOmega)$ is self-adjoint and becomes an isomorphism between Hilbert spaces if its domain is equipped with the graph norm defined via

[TABLE]

for $g\in\mathcal{D}(\varLambda^{-s}(\sigma))$ and with $C>0$.

Since $\varLambda(\sigma)\colon H^{-1/2}_{\diamond}(\partial\varOmega)\to H^{1/2}_{\diamond}(\partial\varOmega)$ is positive and self-adjoint, the (positive) square root $\varLambda^{1/2}(\sigma)$ is, in fact, an isomorphism from $L^{2}_{\diamond}(\partial\varOmega)$ to $H^{1/2}_{\diamond}(\partial\varOmega)$; see, e.g., [2, Lemma 3.4] for a simple proof. In particular, $\mathcal{D}(\varLambda^{-1/2}(\sigma))=H^{1/2}_{\diamond}(\partial\varOmega)$. According to [20, p. 10, Definition 2.1 & Remark 2.3] and the definition of $H^{s}(\partial\varOmega)$ as a (complex) interpolation space (see, e.g., [20, p. 36, Theorem 7.7]), we thus have

[TABLE]

with the graph norm of $\mathcal{D}(\varLambda^{-s}(\sigma))$ being equivalent to that of $H^{s}_{\diamond}(\partial\varOmega)$. The claim now follows from (20). ∎

Proof of Corollary 1. By virtue of Lemma 1 and the asymptotic dominance of $\lambda_{k}^{-\epsilon}$ over $\log^{2}(\lambda_{k})$ as $k$ tends to infinity,

[TABLE]

for any $\epsilon>0$ and $f\in H^{\epsilon/2}_{\diamond}(\partial\varOmega)$, with $C=C(\epsilon)>0$ that may change between different occurrences. Since $\log\!\varLambda(\sigma)\colon L^{2}_{\diamond}(\partial\varOmega)\supset\mathcal{D}(\log\!\varLambda(\sigma))\to L^{2}_{\diamond}(\partial\varOmega)$ is self-adjoint and the embedding $H^{\epsilon}_{\diamond}(\partial\varOmega)\hookrightarrow H^{\epsilon/2}_{\diamond}(\partial\varOmega)$ is compact, the assertion follows. $\square$

So far we have discussed the logarithm of the linear Neumann-to-Dirichlet operator. A natural way to define the corresponding nonlinear forward operator is as follows (cf. (10) and (12)):

Definition 1**.**

The logarithmic forward operator of EIT is

[TABLE]

for some $\epsilon>0$.

Considering the conductivity logarithm $\kappa$ instead of $\sigma$ affords the same convenience as with the non-logarithmic forward operator $\varLambda_{\exp}$ in (12). Moreover, the linear nature of $\varLambda_{\rm inv}$ in (11) with constant conductivities (see Example 1) becomes now available as well, as shown next.

Example 2** (Constant conductivity revisited).**

If $\sigma\equiv\mathrm{e}^{\kappa}$ is constant in $\varOmega$, then the mapping $\kappa\mapsto L(\kappa)$ is affine.

Proof.

As in Example 1, the potential $u$ from (1) solves the Laplace equation with the normal derivative $\mathrm{e}^{-\kappa}f$. Thus, the Neumann-to-Dirichlet operator satisfies

[TABLE]

where the eigensystem corresponds to the unit conductivity. Taking the logarithm results in

[TABLE]

which proves the claim. ∎

The following example shows that the logarithmic forward operator depends “almost linearly” on the conductivity logarithm of an interior inclusion.

Example 3** (Nested concentric disk).**

Let $\varOmega\subset\mathbf{R}^{2}$ be the unit disk and let $\varOmega_{R}\subset\varOmega$ be an origin-centered disk with radius $0<R<1$. Define the conductivity via its logarithm $\kappa=\log(\sigma)$ as

[TABLE]

for some $\widetilde{\kappa}\in\mathbf{R}$. Now the logarithmic forward operator satisfies $L(\kappa)=L(0)+\widetilde{\kappa}\,L^{\prime}(0)+O(\widetilde{\kappa}^{3})$ for some $L^{\prime}(0)\in\mathcal{L}(L^{2}_{\diamond}(\partial\varOmega))$ and with $O(\widetilde{\kappa}^{3})$ referring to the topology of $\mathcal{L}(L^{2}_{\diamond}(\partial\varOmega))$.

Proof.

It is shown in [13] that in this rotationally symmetric case the eigenfunctions of $\varLambda(\sigma)$ do not depend on $\widetilde{\kappa}$, and the eigenvalues are

[TABLE]

where

[TABLE]

The logarithm of the eigenvalue $\lambda_{2k}$ satisfies

[TABLE]

which is an analytic, odd function with respect to $\tilde{\kappa}$. In particular, the second derivative of $\log(\lambda_{2k})$ with respect to $\widetilde{\kappa}$ vanishes at the origin, and it is also straightforward, yet tedious to check that the corresponding first and third derivatives are uniformly bounded over $k\in\mathbf{N}$ in a neighborhood of the origin. The claim thus follows by plugging (27) in (16) and using Taylor’s theorem. ∎

3.2 Fréchet differentiability

Let us next consider the differentiability of $\log\!\varLambda(\sigma)$ with respect to the conductivity. Since $\log\!\varLambda(\sigma)$ has not itself been introduced as a boundary operator corresponding to some (elliptic) partial differential equation, but its definition directly involves an eigensystem for $\varLambda(\sigma)$, the natural way to start would be to consider the differentiability of $\{\lambda_{k},\phi_{k}\}_{k=1}^{\infty}\subset\mathbf{R}_{+}\times L^{2}_{\diamond}(\partial\varOmega)$. Indeed, the Gateaux differentiability of (a suitable parametrization for) individual eigenvalues and eigenfunctions with respect to the conductivity could be proven by considering the corresponding eigensystem $\{\lambda_{k}^{-1},\phi_{k}\}_{k=1}^{\infty}$ of the (unbounded) Dirichlet-to-Neumann operator $\varLambda^{-1}(\sigma)$ and utilizing its analyticity with respect to perturbations in the conductivity; see, e.g., [3] for the analyticity of $\varLambda^{-1}(\,\cdot\,)$ and [16, p. 392, Theorem 3.9] or [17] for more information about the differentiability of an eigensystem for a self-adjoint, unbounded operator with respect to a real parameter that does not affect the domain of definition. However, due to the singularity of the logarithm at the origin and the infiniteness of the sum in (16), (the natural sense of) the differentiability of $\log\!\varLambda(\sigma)$ is not trivial to establish. Moreover, as our main motivation for differentiating $\log\!\varLambda(\sigma)$ with respect to the conductivity is numerical computing as in Section 5, we only consider derivatives for finite-dimensional approximations of $\sigma\mapsto\log\!\varLambda(\sigma)$ in what follows.

To be more precise, we investigate a finite-dimensional, positive and self-adjoint mapping

[TABLE]

where

[TABLE]

is an $L^{2}(\partial\varOmega)$-orthogonal projection composed of individual projections $P_{k}$ onto fixed orthonormal basis functions $\psi_{k}\in L^{2}_{\diamond}(\partial\varOmega)$, $k=1,\dots,n$, respectively. Furthermore, let

[TABLE]

be a spectral representation for $\varLambda^{(n)}(\sigma)$, with $\{\mu_{k}\}_{k=1}^{n}\subset\mathbf{R}_{+}$ and $\{\varphi_{k}\}_{k=1}^{n}\subset{\rm span}\{\psi_{1},\dots,\psi_{n}\}$ being orthonormal. Take note that $\{\mu_{k}\}_{k=1}^{n}$ and $\{\varphi_{k}\}_{k=1}^{n}$ depend on $\sigma$, but the basis for the discretization $\{\psi_{k}\}_{k=1}^{n}$ does not. The logarithm of $\varLambda^{(n)}(\sigma)$ is defined in the same way as for its infinite-dimensional counterpart:

[TABLE]

which is obviously a bounded, linear map from $L^{2}_{\diamond}(\partial\varOmega)$ to itself.

Theorem 1**.**

The mapping

[TABLE]

is Fréchet differentiable. The derivative at $\sigma\in L^{\infty}_{+}(\varOmega)$ in the direction $\eta\in L^{\infty}(\varOmega)$ is the element of $\mathcal{L}(L_{\diamond}^{2}(\partial\varOmega))$ defined via

[TABLE]

where $D\varLambda(\sigma;\eta):L^{2}_{\diamond}(\partial\varOmega)\to L^{2}_{\diamond}(\partial\varOmega)$ is the Fréchet derivative of $\varLambda(\sigma)$ and the divided difference is extended for recurrent eigenvalues in the standard limit sense, i.e.,

[TABLE]

if $\mu_{j}=\mu_{k}$.

Proof.

Obviously, the finite-dimensional approximation $\sigma\mapsto\varLambda^{(n)}(\sigma)$ inherits Fréchet differentiability from $\sigma\mapsto\varLambda(\sigma)$. In fact, the derivative of the former at $\sigma\in L^{\infty}_{+}(\varOmega)$ is simply

[TABLE]

The principal branch of a logarithm of a matrix (or of a finite-dimensional linear operator given in a fixed basis) is also Fréchet differentiable over the set of matrices that do not have eigenvalues on the closed negative real axis in the complex plain [12]. The Fréchet differentiability of the mapping (31) thus follows from the chain rule for Banach spaces, which also gives the representation (32) when combined with [12, Theorem 3.11 & Corollary 3.12] and (33). ∎

If one chooses the discretization frame to be $\psi_{k}=\phi_{k}\in H^{1/2}_{\diamond}(\partial\varOmega)$, $k=1,\dots,n$, i.e., the first $n$ orthonormal eigenfunctions of $\varLambda(\sigma)$ for a particular, fixed $\sigma\in L^{\infty}_{+}(\varOmega)$, then

[TABLE]

that is, $\log\!\varLambda^{(n)}(\sigma)$ inherits its eigenvalues and eigenfunctions from $\log\!\varLambda(\sigma)$. In this case,

[TABLE]

Recall that here the basis for the discretization is fixed. In other words, the projection $P^{(n)}$ in (34) is considered invariable when the Fréchet derivative is computed, that is, it corresponds all the time to the first $n$ orthonormal eigenfunctions of the unperturbed $\varLambda(\sigma)$, not to those of $\varLambda(\sigma+\eta)$.

Remark 1**.**

The derivative of the finite-dimensional logarithmic forward operator of EIT,

[TABLE]

is obtained by simply replacing $D\varLambda(\sigma;\eta)$ with $D\varLambda_{\exp}(\kappa;\eta)$ of (15) in (32).

4 Complete electrode model

Practical EIT measurement setups are characterized by $M\in\mathbf{N}\setminus\{1\}$ electrodes that are attached to the surface of the imaged object. Let $E_{m}\subset\partial\varOmega$, $m=1,\ldots,M$, denote these nonempty, open, connected, well-separated sets. In the CEM, the Neumann boundary condition in (1) is replaced by [5]

[TABLE]

where $\zeta=[\zeta_{m}]_{m=1}^{M}\in\mathbf{R}_{+}^{M}$ denote the (constant) contact conductances at the electrode-object interfaces. The current pattern $I=[I_{m}]_{m=1}^{M}$ is prescribed, whereas the electrode potentials $U=[U_{m}]_{m=1}^{M}$ are part of the solution. Both of these vectors belong to the mean-free subspace $\mathbf{C}_{\diamond}^{M}\subset\mathbf{C}^{M}$, the current pattern due to the conservation of electric charge and the electrode potentials as a result of our specific choice for the ground level of potential. The variational formulation corresponding to the conductivity equation from (1) with the boundary conditions (36) is to find $(u,U)\in\mathcal{H}^{1}(\varOmega):=H^{1}(\varOmega)\oplus\mathbf{C}^{M}_{\diamond}$ such that

[TABLE]

holds for all $(v,V)\in\mathcal{H}^{1}(\varOmega)$. The problem (39) is uniquely solvable [21].

For simplicity, we use the same notations for the forward operators of CEM as introduced previously for the continuum model in Section 2. However, as a disparity, the CEM forward operators depend on the contact conductances in addition to the conductivity. As an example, the CEM forward operator that takes the conductivity and contact conductances as its inputs is of the form

[TABLE]

where the output is called the resistance matrix (regardless of the used input parametrization), which is given with respect to a fixed orthonormal basis for $\mathbf{C}_{\diamond}^{M}$. In other words, if $\alpha\in\mathbf{C}^{M-1}$ carries the coordinates of an electrode current pattern $I\in\mathbf{C}_{\diamond}^{M}$ with respect to the employed basis, then $\varLambda(\sigma,\zeta)\alpha$ gives the coordinates of the resulting electrode potentials $U\in\mathbf{C}_{\diamond}^{M}$ in that same basis. Since $\sigma$ and $\zeta$ are real-valued and positive, it follows easily from (39) that $\varLambda(\sigma,\zeta)$ is Hermitian and, in particular, positive-definite. Similarly,

[TABLE]

and

[TABLE]

where we have abused the notation by defining $z^{-1}=[z_{m}^{-1}]_{m=1}^{M}$ and ${\rm e}^{\upsilon}=[{\rm e}^{\upsilon_{m}}]_{m=1}^{M}$ as vectors of $\mathbf{R}^{M}$. Finally, the logarithmic forward operator is defined in accordance with Definition 1, that is,

[TABLE]

where $\log\!\varLambda_{\exp}(\kappa,\upsilon)$ is the principal logarithm of the positive-definite matrix $\varLambda_{\exp}(\kappa,\upsilon)$.

The Fréchet derivative of the basic forward operator (41) with respect to $\sigma\in L^{\infty}_{+}(\varOmega)$ in the direction $\eta\in L^{\infty}(\varOmega)$ is characterized by [14]

[TABLE]

where $(u_{\beta},U_{\beta})$ and $(u_{\alpha},U_{\alpha})$ are the solutions of (39) corresponding to the current patterns with the coordinates $\alpha$ and $\beta$, respectively, in the employed orthonormal basis for $\mathbf{C}_{\diamond}^{M}$. Analogously, the derivative of $\zeta\mapsto\varLambda(\sigma,\zeta)$ in the direction $\xi\in\mathbf{R}^{M}$ can be assembled via

[TABLE]

Furthermore, the Fréchet derivatives of $\varLambda_{\rm inv}$ and $\varLambda_{\exp}$ can be easily obtained by utilizing the chain rule and the above two formulas for $\varLambda$ (cf. (14) and (15)), and those of the logarithmic forward operator $L$ can be deduced by writing an eigendecomposition for $\varLambda(\sigma,\zeta)$ and mimicking Theorem 1 and Remark 1.

Remark 2**.**

Applying the three basic operations ${\rm id}$, ${\rm inv}$ and $\exp$ to the conductivity and the contact conductances, one could come up with nine different parametrizations for the forward operator of the CEM. Here we only consider three of them, although it is probable that in some cases one of the other options is optimal from the standpoint of the linearization error.

Remark 3**.**

When the contact conductances are perfect, i.e., in the limit $\zeta_{m}\to\infty$ for all $m=1,\dots,M$, the resulting model is called the shunt model [5, 6]. On the other hand, when $\sigma/\zeta_{m}$ approaches infinity for all $m=1,\dots,M$, the model formally approaches a resistor network for which the resistance between electrodes $E_{l}$ and $E_{m}$ is $R_{lm}=(\zeta_{l}|E_{l}|)^{-1}+(\zeta_{m}|E_{m}|)^{-1}$. In particular, if $\zeta_{l}^{-1}=\zeta_{m}^{-1}=:\tilde{z}$ and $|E_{m}|=|E_{l}|=:A$, then $R_{lm}=2\tilde{z}A^{-1}$. Hence, for high contact resistances $z$ and low resistivity $\rho$, it is reasonable to expect that of the different forward operators $\varLambda_{\rm inv}$ is closest to being linear.

5 Numerical experiments

In this section, we study the accuracies of different linearization techniques by performing numerical experiments in the unit disk $\varOmega\subset\mathbf{R}^{2}$. The elliptic forward problems for both the continuum model and the CEM are solved with the finite element method (FEM) on meshes with approximately $30\,000$ nodes and piecewise linear basis functions. These meshes and corresponding FEM solutions are also employed when evaluating derivatives of forward maps (cf., e.g., (13)–(15)). The discretization errors of the FEM solutions can be regarded to be small. Throughout this section, the notation $\mathbf{E}[\cdot]$ is used to denote the sample mean. The size of the sample is $50\,000$ in each considered case and the random realizations of involved parameters are mutually independent.

5.1 Forward accuracy

We study the linearization errors of the EIT forward problem with discrete lognormal random conductivity fields. To this end, the domain $\varOmega$ is divided into $N=1800$ subdomains $\varOmega_{i}$, each having approximately the same area, and the piecewise constant conductivity is written as

[TABLE]

where $\chi_{i}$ is the indicator function of $\varOmega_{i}$ and the vector $\hat{\sigma}\in\mathbf{R}_{+}^{N}$ contains the corresponding conductivity values. The resistivity $\rho=\sigma^{-1}$ and the logarithm of the conductivity $\kappa=\log(\sigma)$ are defined similarly via vectors $\hat{\rho}\in\mathbf{R}_{+}^{N}$ and $\hat{\kappa}\in\mathbf{R}^{N}$, respectively. Furthermore, let $\hat{x}_{i}\in\varOmega_{i}$ denote the center of a subdomain. The random conductivities are drawn via their logarithms such that the vectors $\hat{\kappa}$ follow a Gaussian distribution with the probability density

[TABLE]

where $\hat{\kappa}_{0}\in\mathbf{R}^{N}$ specifies the discrete mean field and

[TABLE]

defines the covariance matrix for some variance parameter $\varsigma^{2}>0$ and correlation length $\ell>0$. Table 1 lists the used values for four different distributions.

First we consider the continuum model. For the purpose of numerical computations, the finite-dimensional approximation $\varLambda^{(n)}(\sigma)$ in (28) is considered with $n=16$ trigonometric basis functions

[TABLE]

where $\theta$ denotes the polar angle. In order to obtain a fully discrete forward operator, the conductivity $\sigma$ is replaced by $\hat{\sigma}$, resulting in

[TABLE]

where $\hat{\varLambda}(\hat{\sigma})$ is the matrix representation of $P^{(n)}\varLambda(\sigma)P^{(n)}$ for the conductivity (46); see (29). This nonlinear matrix-valued function is linearized around the unit conductivity $\sigma_{0}\equiv 1$ as

[TABLE]

where $\hat{\varLambda}^{\prime}(\hat{\sigma})\hat{\eta}$ is the Fréchet derivative of the mapping (50) in the direction $\hat{\eta}\in\mathbf{R}^{N}$. To be more precise,

[TABLE]

which can be evaluated with the help of (13); see also (33). The alternative fully discrete forward operators $\hat{\rho}\mapsto\hat{\varLambda}_{\rm inv}(\hat{\rho})$ and $\hat{\kappa}\mapsto\hat{\varLambda}_{\exp}(\hat{\kappa})$ are defined analogously, as are their linearizations around $\rho_{0}\equiv 1$ and $\kappa_{0}\equiv 0$, respectively, with the appropriate derivative from (14) or (15) replacing $D\varLambda(\sigma;\chi_{i})$ in (52). The same conclusions also apply to the fully discrete logarithmic forward operator $\hat{\kappa}\mapsto\hat{L}(\hat{\kappa})$; cf. Definition 1, Theorem 1 and Remark 1 with $\eta=\chi_{i}$. Note that for a rotationally symmetric conductivity, (49) are eigenfunctions of the Neumann-to-Dirichlet operator, meaning that $\varphi_{k}=\phi_{k}=\psi_{k}$ in the notation of Section 3 (cf. [13]). In particular, (34) and (35) hold.

The error indicator for the linearized $\hat{\varLambda}$ is

[TABLE]

where $\lVert\cdot\rVert_{F}$ is the Frobenius norm. The corresponding definitions of indicators for $\hat{\varLambda}_{\rm inv}$ and $\hat{\varLambda}_{\exp}$ should be obvious. For the logarithmic forward operator, the error is computed as

[TABLE]

where $\exp$ is just the ordinary matrix exponent, i.e., the inverse of (the principal branch of) the matrix logarithm.

The left-hand side of Table 2 lists the mean errors for the four linearization approaches and for the four different random fields defined in Table 1. It can be seen that linearizing with respect to the conductivity results in the worst accuracy in each case, whereas the novel logarithmic linearization method is clearly the best. Increasing the correlation length in the random field deteriorates the performance of the conductivity and log-conductivity linearizations, whereas $\hat{\varLambda}_{\rm inv}$ and $\hat{L}$ perform better with the smoother fields. An intuitive explanation for the latter phenomenon is that increasing the correlation length in our random model makes the corresponding realizations of the conductivity to be closer to constants, for which $\hat{\varLambda}_{\rm inv}$ and $\hat{L}$ are linear and affine, respectively; see Examples 1 and 2. It is not surprising that increasing the pointwise variance $\varsigma^{2}$ leads in every case to a higher sample mean for the linearization error.

For the CEM experiments we employ the same lognormal conductivity fields F3–F4 as with the continuum model, but now the contact conductances $\zeta_{m}$ are also random. The conductances are drawn independently and for each of the $M=16$ electrodes the conductance distribution is lognormal with the underlying Gaussian distribution having mean $\upsilon_{0}$ and variance $\gamma^{2}$ as listed on the right in Table 1. The identical electrodes are positioned uniformly and they cover approximately $46\%$ of the boundary $\partial\varOmega$. The orthonormal electrode current basis for $\mathbf{R}_{\diamond}^{M}$, employed when forming the resistance matrices, is the discrete counterpart of the Fourier basis (49), that is,

[TABLE]

and $I^{(M-1)}=[(-1)^{m-1}]_{m=1}^{M}/\sqrt{M}$. Here, $\theta_{m}=2\pi(m-1)/M$ is the central polar angle of the $m$th electrode. The definitions of the fully discrete CEM forward operators, their linearizations and the corresponding error indicators should be obvious. In particular, the forward operators are now mappings from $\mathbf{R}_{+}^{N+M}$ or $\mathbf{R}^{N+M}$ to $\mathbf{R}^{(M-1)\times(M-1)}$ and their derivatives can be assembled using the appropriate variants of (44), (45) and (35).

For the contact conductance distribution C1 with the lower expected value, the resistivity-based operator $\varLambda_{{\rm inv}}$ results in the best linearization accuracy, as can be seen on the right in Table 2. This is in line with Remark 3. On the other hand, when the contact conductance is high (i.e., the setting is close to the shunt model), the errors almost coincide with the corresponding values for the continuum model. That is, the rows corresponding to C2 are similar to the adjacent rows on the left-hand side of Table 2. In particular, the logarithmic forward operator $\hat{L}$ gives the best and the conductivity-based parametrization the worst linearization accuracy for the CEM as well.

5.2 Inverse accuracy

Regarding the accuracy of linearizations when solving the inverse problem of EIT, we only consider the forward operators $\hat{\varLambda}_{\exp}$ and $\hat{L}$ that are based on the logarithmic input parametrization. One reason for this is the fact that non-logarithmic parametrizations may lead to conductivity reconstructions that are not positive. In addition, it would be problematic to design equivalent regularizations for different input parametrizations, and thus the error indicators would not be directly comparable to each other. In the following, we give all definitions for the logarithmic forward operator $\hat{L}$, but the corresponding definitions for $\hat{\varLambda}_{\exp}$ should be obvious.

The simulated finite-dimensional Neumann-to-Dirichlet operators are denoted by $\widetilde{\varLambda}(\hat{\sigma})$. They are also computed with FEM, but in order to avoid inverse crimes, the matrices $\widetilde{\varLambda}(\hat{\sigma})$ are contaminated with Gaussian noise. More precisely, realizations of independent, zero-mean normal random variables with standard deviation

[TABLE]

are added to each element of a matrix $\hat{\varLambda}(\hat{\sigma})$, and the matrix is subsequently symmetricized. The simulated CEM measurements are obtained in the same way and they are denoted by $\widetilde{\varLambda}(\hat{\sigma},\zeta)$.

For a regularization parameter $t>0$, the reconstruction $\hat{\kappa}_{t}$ based on the logarithmic continuum model forward operator $\hat{L}$ is

[TABLE]

where $G\in\mathbf{R}^{N\times N}$ is a ‘Bayesian regularization matrix’ satisfying $\varGamma^{-1}=G^{\top}G$. Here, $\varGamma$ is the covariance matrix (48) of the distribution from which the logarithm of the vector $\hat{\sigma}$ is drawn, and the expected log-conductivity $\hat{\kappa}_{0}$ is as in Section 5.1. The reconstruction obtained from (56) approximates a maximum a posteriori (MAP) estimate under the assumption that the linearization error is negligible and the elements of the matrix $\log(\widetilde{\varLambda}(\hat{\sigma}))$ are contaminated by independent realizations of zero-mean Gaussian noise with standard deviation $t>0$ [15]. (Notice that the latter is not true for any $t$ due to the application of the matrix logarithm after adding the measurement noise.) Naturally, for $\hat{\varLambda}_{\exp}$ the minimization does not involve a logarithm, but the linearized operator is compared directly to the matrix $\widetilde{\varLambda}(\hat{\sigma})$ in the Frobenius norm. For the CEM forward operators, the reconstructions are computed according to

[TABLE]

or by replacing $\hat{L}_{\text{lin}}(\hat{y},w)$ with $(\hat{\varLambda}_{\exp})_{\rm lin}(\hat{y},w)$ and deleting the logarithm. In particular, no regularization is applied to the contact conductances since their estimation is a relatively stable task (cf, e.g., [10, (4.8)]). Observe that (56) and (57) correspond to minimizing quadratic Tikhonov functionals, which is a simple task.

The average reconstruction error for a fixed regularization parameter $t>0$ is

[TABLE]

where $\kappa=\log(\sigma)$ is the (random draw of a) true log-conductivity corresponding to the vector $\hat{\sigma}$ in (56); see (46). For a given conductivity/conductance distribution, the optimal regularization parameter is obtained by solving a one-dimensional optimization problem,

[TABLE]

and the corresponding error indicator is simply denoted by $\iota(\hat{L})\coloneqq\iota_{\tau(\hat{L})}(\hat{L})$. (Naturally, such an optimal regularization parameter cannot be found in practice when the true conductivity is not known, but the idea here is to compare upper bounds for the performance of the two considered one-step reconstruction algorithms.) In addition, for the CEM we denote the contact conductance reconstruction error by

[TABLE]

where $\upsilon=\log(\zeta)$ is defined componentwise, but it should be emphasized that the optimal regularization parameter $\tau(\hat{L})$ is still chosen so that it minimizes the $L^{2}(\varOmega)$ reconstruction error corresponding to the mere log-conductivity. Once again, the corresponding definitions for $\varLambda_{\exp}$ are analogous.

Table 3 reveals that the new logarithmic method is more accurate in each test case, and the difference is highlighted with the smoother fields F2 and F4, as could have been predicted based on the forward accuracy observations in Table 2. In practice, the choice of regularization parameter may well be suboptimal. However, Figure 1 demonstrates that the new method results in a smaller average error with almost any choices for the regularization parameters.

In order to illustrate the difference in the one-step inversion accuracies for $\hat{\varLambda}_{\exp}$ and $\hat{L}$, some example reconstructions, with regularization parameters that are optimized as described above, are shown in Figures 2 and 3. The first example, which is based on a realization of the random conductivity field F3 and the CEM with high contact conductances, demonstrates that the traditional linearization method is not capable of recovering the high contrast of the target conductivity. On the other hand, the new method reproduces the conductivity levels far more accurately and also recovers the resistive area close to the northwest boundary of the unit disk. The second example presented in Figure 3 considers a realization of the random conductivity field F4 with twice as long correlation length. The superiority of the new logarithmic method is even more clearly visible with this smoother conductivity field, as was expected based on the statistical evidence in Tables 2 and 3.

6 Discussion

We have reviewed existing linearization approaches for EIT and proposed a new logarithmic technique that seems to be the most accurate linearization method amongst those studied here. Although the EIT forward operator is somewhat frequently linearized with respect to the conductivity, it was numerically demonstrated that the linearization error in the measurement matrix becomes smaller if the linearization is based on the resistivity or the logarithm of the conductivity. In a sense, this is just a consequence of the Ohm’s law, and the conclusion applies to both the continuum model and to the CEM. Our novel method, which maps the logarithm of the conductivity to the logarithm of the Neumann-to-Dirichlet operator or the measurement matrix, produces the smallest linearization error in almost all tested cases.

Regarding the inverse problem of EIT, the proposed logarithmic method retains the positivity property of the log-conductivity parametrization, while preserving the accuracy of the resistivity linearization in the case of, e.g., constant conductivities. Numerical studies and example reconstructions demonstrate that when comparing to the traditional log-conductivity linearization approach, the proposed method is clearly more accurate, regardless of the parameters defining the lognormal random models for the conductivity and the contact conductances.

The above conclusions are valid if the measurements are modelled by the Neumann-to-Dirichlet map or its counterparts for electrode measurements. If the Dirichlet-to-Neumann map were employed, one would expect conductivity parametrizations to prevail over those based on the resistivity (cf. Ohm’s law). However, such a change would not affect the performance of the novel method since the logarithms of the Neumann-to-Dirichlet and Dirichlet-to-Neumann operators only differ by their sign.

Acknowledgments

This work was supported by the Academy of Finland (decision 267789) and the Finnish Foundation for Technology Promotion TES.