Smoothened complete electrode model

TL;DR

This paper introduces a smoothened boundary conductance model for electrical impedance tomography that enhances numerical solution efficiency and convergence by improving the regularity of the electromagnetic potential.

Contribution

It proposes a novel smoothened boundary conductance model, proving its unique solvability and better regularity, leading to faster finite element method convergence.

Findings

The smoothened model is computationally feasible.

It improves convergence speed of numerical methods.

Compatible with real-world measurement data.

Abstract

This work reformulates the complete electrode model of electrical impedance tomography in order to enable more efficient numerical solution. The model traditionally assumes constant contact conductances on all electrodes, which leads to a discontinuous Robin boundary condition since the gaps between the electrodes can be described by vanishing conductance. As a consequence, the regularity of the electromagnetic potential is limited to less than two square-integrable weak derivatives, which negatively affects the convergence of, e.g., the finite element method. In this paper, a smoothened model for the boundary conductance is proposed, and the unique solvability and improved regularity of the ensuing boundary value problem are proven. Numerical experiments demonstrate that the proposed model is both computationally feasible and also compatible with real-world measurements. In particular, the new model allows faster convergence of the finite element method.

Equations

$$\mathbb{C}^{M}_{\diamond}\,:=\,\Big{\{}J\in\mathbb{C}^{M}\,\Big{|}\,\sum_{m=1}^{M}J_{m}=0\Big{\}}.$$

$$U=\sum _{m=1}^M{U}_m{\chi }_m$$

$$\mathrm{Re}(\sigma \xi \cdot \overline{\xi })\ge {\varsigma }_{-}\Vert \xi {\Vert }_2^2,$$

$$\begin{array}[]{ll}\displaystyle{\nabla\cdot(\sigma\nabla u)=0}&{\rm in}\;\;\Omega,\\[6.0pt] {\displaystyle{\nu\cdot\sigma\nabla u}=\zeta(U-u)}&{\rm on}\;\;\partial\Omega,\\[2.0pt] {\displaystyle\int_{E_{m}}\nu\cdot\sigma\nabla u\,{\rm d}S}=I_{m},&m=1,\ldots,M,\\[4.0pt] \end{array}$$

$${\rm Re}(\zeta)\geq 0,\qquad\zeta_{\partial\Omega\setminus{\overline{E}}}\equiv 0,\qquad{\rm Re}\big{(}\zeta|_{E_{m}}\big{)}\not\equiv 0$$

$$\begin{array}[]{ll}{\displaystyle{\nu\cdot\sigma\nabla u}=0}&{\rm on}\;\;\partial\Omega\setminus\mkern 1.0mu\overline{\mkern-1.0muE\mkern-1.0mu}\mkern 1.0mu,\\[6.0pt] {\displaystyle{u+z_{m}\nu\cdot\sigma\nabla u}=U_{m}}&{\rm on}\;\;E_{m},\quad m=1,\dots,M,\end{array}$$

$${\mathcal{H}}^s$$

$$\displaystyle\simeq\Big{\{}\big{\{}\big{(}v+c,\sum_{m}(V_{m}+c)\chi_{m}\big{)}\,|\,c\in\mathbb{C}\big{\}}\,\Big{|}\,(v,V)\in H^{s}(\Omega)\oplus\mathbb{C}^{M}\Big{\}}$$

$$\|(v,V)\|_{\mathcal{H}^{1}}=\inf_{c\in\mathbb{C}}\Big{(}\|v-c\|_{H^{1}(\Omega)}^{2}+\|V-c{\bf 1}\|_{2}^{2}\Big{)}^{1/2}.$$

$$B\big{(}(u,U),(v,V)\big{)}\,=\,I\cdot\mkern 1.0mu\overline{\mkern-1.0muV\mkern-1.0mu}\mkern 1.0mu\qquad{\rm for}\ {\rm all}\ (v,V)\in\mathcal{H}^{1},$$

$$B\big{(}(w,W),(v,V)\big{)}=\int_{\Omega}\sigma\nabla w\cdot\nabla\mkern 1.0mu\overline{\mkern-1.0muv\mkern-1.0mu}\mkern 1.0mu\,{\rm d}x+\int_{\partial\Omega}\zeta(W-w)(\mkern 1.0mu\overline{\mkern-1.0muV\mkern-1.0mu}\mkern 1.0mu-\mkern 1.0mu\overline{\mkern-1.0muv\mkern-1.0mu}\mkern 1.0mu)\,{\rm d}S,$$

$$\big{|}B\big{(}(w,W),(v,V)\big{)}\big{|}\,\leq\,C\|(w,W)\|_{\mathcal{H}^{1}}\|(v,V)\|_{\mathcal{H}^{1}}$$

$${\rm Re}\Big{(}B\big{(}(v,V),(v,V)\big{)}\Big{)}\,\geq\,c\|(v,V)\|_{\mathcal{H}^{1}}^{2},$$

$$\mathrm{Re}(\zeta )\ge {\zeta }_{-}a.e.\mathrm{on}e:=\bigcup _{m=1}^M{e}_m.$$

$$\|V-c{\bf 1}\|_{2}^{2}\leq C\|V-c\|_{L^{2}(e)}^{2}\leq C\Big{(}\|V-v\|_{L^{2}(e)}^{2}+\|v-c\|_{L^{2}(\partial\Omega)}^{2}\Big{)}$$

$$\Vert (v,V){\Vert }_{{\mathcal{H}}^1}^2$$

$$\displaystyle\leq C\Big{(}\|\nabla v\|_{L^{2}(\Omega)}^{2}+\|V-v\|_{L^{2}(e)}^{2}\Big{)},$$

$$\displaystyle{\rm Re}\Big{(}B\big{(}(v,V),(v,V)\big{)}\Big{)}$$

$$\Vert (u,U){\Vert }_{{\mathcal{H}}^1}\le C\Vert I{\Vert }_2,$$

$${\mathcal{H}}^1\ni (v,V)\mapsto I\cdot \overline{V}\in \mathbb{C}$$

$$\big{|}I\cdot\overline{V}\big{|}\leq\|I\|_{2}\inf_{c\in\mathbb{C}}\|V-c{\bf 1}\|_{2}\leq\|I\|_{2}\,\|(v,V)\|_{\mathcal{H}^{1}},$$

$${\mathcal{M}}_{\eta }:v\mapsto \eta v,H^r(\Gamma )\to H^r(\Gamma ),$$

$$\Vert \eta v{\Vert }_{H^r(\Gamma )}\le C\Vert \eta {\Vert }_{H^s(\Gamma )}\Vert v{\Vert }_{H^r(\Gamma )},-s<r\le s,$$

$$\Vert W-w{\Vert }_{H^r(E)}\le C\Vert (w,W){\Vert }_{{\mathcal{H}}^{r+1/2}},$$

$$\Vert W-w{\Vert }_{H^r(E)}$$

$$\displaystyle\leq 2\inf_{c\in\mathbb{C}}\Big{(}\|w-c\|_{H^{r}(\partial\Omega)}^{2}+\|W-c\|_{H^{r}(E)}^{2}\Big{)}^{1/2}$$

$$\le C(r,E,\Omega )\Vert (w,W){\Vert }_{{\mathcal{H}}^{r+1/2}},$$

$$\Vert (u,U){\Vert }_{{\mathcal{H}}^{s+3/2}}\le C\Vert I{\Vert }_2.$$

$$\nabla \cdot (\sigma \nabla v)=0\mathrm{in}\Omega ,\nu \cdot \sigma \nabla v=f\mathrm{on}\partial \Omega$$

$$H_{\diamond }^r(\partial \Omega )=\{g\in H^r(\partial \Omega )|⟨g,1{⟩}_{\partial \Omega }=0\},r\in \mathbb{R},$$

$$\Vert v{\Vert }_{H^{r+3/2}(\Omega )/\mathbb{C}}\le C\Vert f{\Vert }_{H^r(\partial \Omega )}.$$

$${\Lambda }_{\sigma }:f\mapsto v{|}_{\partial \Omega },H_{\diamond }^r(\partial \Omega )\to H^{r+1}(\partial \Omega )/\mathbb{C},$$

Full text

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

1. Introduction

Electrical impedance tomography (EIT) is a noninvasive imaging technique based on controlling and measuring electric currents and voltages on the surface of the imaged object. The aim is to reconstruct the electrical conductivity (or admittivity, resistivity, or impedivity) inside the object. Applications of EIT include biomedical imaging [1], nondestructive testing [17], and process tomography [27]. The reconstruction task is an extensively studied inverse problem for which both direct and iterative methods have been proposed [16, 22, 30].

Several mathematical models have been applied to incorporating boundary conditions in the forward problem of EIT including those discussed in [2, 28]. The simplest one is the continuum model, which assumes a (typically smooth) Neumann boundary condition. This is useful for theoretical considerations and convenient in numerical computations, but does not typically result in accurate reconstructions because practical measurement setups employ a finite number of electrodes, which is not taken into account in the continuum model. The point electrode model is mainly useful when the electrodes are small and the reconstruction is based on difference measurements [10]. The so called shunt model correctly models the geometry of the electrodes but neglects the thin resistive layer that may appear at the contact between the electrodes and the object. The presence of this layer is included in the complete electrode model (CEM) that has become the standard for computing reconstructions in practical applications [2]. Loosely speaking, the shunt model can be regarded as the limit of the CEM when the contact resistances tend to zero [8].

When the CEM is employed, the contact conductances (or admittances, resistances, or impedances) are usually not known, but they are estimated along with the interior conductivity [11]. One may also simultaneously reconstruct the electrode locations and the shape of the imaged object [6, 7, 24]. However, even if all these parameters were known, an inherent property of the traditional CEM is that the employed “discontinuous” Robin-type boundary condition causes the regularity of the electromagnetic potential to be limited, namely, it is of the Sobolev class $H^{2-\epsilon}$, $\epsilon>0$. The same conclusion applies to an even greater extent to the shunt model for which the potential only exhibits $H^{3/2-\epsilon}$-regularity. Regarding numerical computations, the lack of smoothness makes it difficult to construct forward solvers that converge fast; in particular, the efficient use of finite element methods (FEM) of higher order is prevented. This is an issue for iterative reconstruction algorithms that require repetitive and accurate solutions of the forward problem. In addition, certain quantities that are derived from forward solutions of the CEM, such as those needed when computing shape derivatives of electrode measurements (cf. [6, 7]), are expected to suffer in accuracy even when first-order FEM is used [7, Remark 2.4].

To overcome these problems arising from the discontinuity of a conductance coefficient in the Robin boundary condition, we propose a smoothened version of the CEM that exploits nonhomogeneous conductances and is more suitable for numerical computations. It is shown that arbitrarily high regularity for the (interior) electromagnetic potential in EIT can be achieved while only slightly deviating from the standard CEM, assuming the conductivity and the object boundary are smooth enough. Moreover, it is numerically demonstrated that forward computations exhibit faster convergence (both in practice and asymptotically). According to our preliminary tests, a version of the smoothened model is in approximately as good agreement with experimental data as the standard CEM and reconstructions based on the two models are almost indistinguishable.

It should be emphasized that we do not claim that the smoothened CEM is a more appropriate model for EIT from the standpoint of the physical phenomena occurring at the electrode-object interface. Our assertion is merely that the new model is computationally more efficient and predicts the electrode measurements of EIT with accuracy comparable to the traditional CEM. How well the new model — or the new family of models — predicts the behavior of the electromagnetic potential in the interior of the imaged object remains an open question; the same is actually true also for the traditional CEM as its performance has only been validated in regard to electrode measurements [2]. (Observe that there is no actual reason to expect that the traditional assumption of having a constant contact conductance on each electrode is completely accurate either. In fact, the standard CEM can be viewed as a special case of the smoothened CEM.)

Although we only discuss EIT, it is worth noting that same kinds of electrode models can be used in, e.g., electrical capacitance tomography (ECT) [9] and electroencephalography (EEG) [25]. In particular, since ECT is mathematically equivalent to EIT, our theoretical results apply directly to ECT as well.

This text is organized as follows. In the next section, we introduce the smoothened CEM for EIT and prove the unique solvability of the corresponding elliptic boundary value problem. Our main theoretical result addressing the regularity of the electromagnetic potential is also formulated and proven in that section. Section 3 reviews the Fréchet differentiability of electrode potentials with respect to the shape of the imaged object, which provides an example of a setting where the new smoothened model clearly prevails. Numerical examples demonstrating the improved convergence are presented in Section 4, where the smoothened model is also compared to the traditional CEM and to experimental data. In addition, example reconstructions based on water tank measurements are presented. Finally, conclusions are drawn in Section 5.

2. Smoothened complete electrode model

A physical body imaged by EIT is modeled as a bounded Lipschitz domain $\Omega\subset\mathbb{R}^{n}$, $n=2$ or $3$. The boundary $\partial\Omega$ is partially covered by $M\in\mathbb{N}\setminus\{1\}$ electrodes $\{E_{m}\}_{m=1}^{M}$ that are identified with nonempty, connected, open surface patches and assumed to be well-separated, i.e., $\overline{E}_{m}\cap\overline{E}_{l}=\emptyset$ if $m\not=l$. We denote $E=\cup E_{m}$. The net currents $I_{m}\in\mathbb{C}$, $m=1,\dots,M$, are driven through the corresponding electrodes and the resulting constant electrode potentials $U_{m}\in\mathbb{C}$, $m=1,\dots,M$, are measured. Due to the conservation of electric charge and under the reasonable assumption that there are no sinks or sources inside $\Omega$, any realizable current pattern $I=[I_{1},\dots,I_{M}]^{\rm T}$ belongs to the subspace

[TABLE]

The electrode potential vector $U=[U_{1},\dots,U_{M}]^{\rm T}$ is identified with

[TABLE]

where $\chi_{m}$ is the characteristic function of $E_{m}\subset\partial\Omega$. Whether $U$ refers to such a piecewise constant function supported on $\overline{E}$ or to a vector of $\mathbb{C}^{M}$ should be clear from the context; in particular, under an integral a capital letter always refers to a piecewise constant function vanishing in between the electrodes. The real-symmetric admittivity distribution $\sigma\in L^{\infty}(\Omega,\mathbb{C}^{n\times n})$ inside $\Omega$ is assumed to satisfy

[TABLE]

for all $\xi\in\mathbb{C}^{n}$ almost everywhere in $\Omega$ with $\varsigma_{-}>0$ being some positive constant. In other words, the imaged body is allowed to be characterized by anisotropic conductivity and permittivity but these coefficients are required to be symmetric and the conductivity, in addition, strictly positive definite (cf., e.g., [31]). Here and in what follows, $\|\cdot\|_{2}$ denotes the Euclidean norm of a finite-dimensional vector.

The CEM is a mathematical model that accurately predicts real-life EIT measurements, i.e., its validity has been confirmed in regard to data collected at electrodes [2]. We consider a nonstandard formulation of the CEM: The electromagnetic potential $u$ inside $\Omega$ and the piecewise constant electrode potential $U$ satisfy

[TABLE]

interpreted in the weak sense. Here, $\nu\in L^{\infty}(\partial\Omega,\mathbb{R}^{n})$ denotes the exterior unit normal of $\partial\Omega$ and $\zeta\in L^{\infty}(\partial\Omega)$ describes the contact admittance over $\partial\Omega$. The gaps between the electrodes can be characterized by vanishing admittance. Moreover, the conductance, i.e., the real part of the admittance, cannot be negative, and to be able to drive currents through the electrodes, the conductance must not vanish everywhere on any of the electrodes. To summarize, it is physically reasonable to assume

[TABLE]

for all $m=1,\dots,M$ in the topology of $L^{\infty}(\partial\Omega)$. Take note that the second assumption on $\zeta$ reduces the second equation of (3) into a homogeneous Neumann condition on $\partial\Omega\setminus{\overline{E}}$, i.e., no current flows through the object boundary in between the electrodes.

A physical justification of (3) can be found in [2], where the second condition is divided into two parts as

[TABLE]

and the contact impedances $z_{m}:=(1/\zeta)|_{E_{m}}$, $m=1,\dots,M$, are assumed to be constants. Although the CEM has previously been analyzed also for nonconstant contact impedances (cf. [12, 32]), to the authors’ knowledge all previous mathematical works on the CEM assume the impedances are bounded away from infinity. As the assumptions (4) allow the contact admittances to vanish on some subsets of the electrodes, the unique solvability of (3) does not directly follow from previous analyses but a bit of extra work is required.

We look for the solution of (3) in the quotient space $\mathcal{H}^{1}$, with the definition

[TABLE]

for $s\in\mathbb{R}$ with ${\bf 1}:=[1,\dots,1]^{\rm T}\in\mathbb{R}^{M}$. Here, “$\simeq$” is to be understood via an isomorphic identification of vectors and piecewise constant functions on the electrodes. The use of a quotient structure reflects the freedom in the choice of the ground level of potential: All elements of $H^{s}(\Omega)\oplus\mathbb{C}^{M}$ that differ by an additive constant are identified as an equivalence class. In particular, when the second component of an element of $\mathcal{H}^{s}$ is interpreted as a piecewise constant function on the electrodes, the additive constant is also supported on $\overline{E}$.

To prove the unique solvability of (3), notice first that the standard quotient norm for $\mathcal{H}^{1}$ is defined by

[TABLE]

Moreover, by following the line of reasoning in [28], one sees that the variational formulation of (3) is to find $(u,U)\in\mathcal{H}^{1}$ that satisfies

[TABLE]

where “$\,\cdot\,$” denotes the real inner product and the sesquilinear form $B:\mathcal{H}^{1}\times\mathcal{H}^{1}\to\mathbb{C}$ is defined by

[TABLE]

with $W,V\in\mathbb{C}^{M}$ identified with the corresponding piecewise constant functions.

Lemma 2.1**.**

Under the assumptions (2) and (4), the sesquilinear form $B:\mathcal{H}^{1}\times\mathcal{H}^{1}\to\mathbb{C}$ is well-defined, bounded and coercive, that is,

[TABLE]

and

[TABLE]

where $c,C>0$ do not depend on $(w,W),(v,V)\in\mathcal{H}^{1}$.

Proof.

First of all, due to the second condition of (4), there is no ambiguity in the definition of $B$ on $\mathcal{H}^{1}\times\mathcal{H}^{1}$, i.e., the value $B((w,W),(v,V))$ does not depend on the particular representatives of the equivalence classes $(w,W),(v,V)\in\mathcal{H}^{1}$. Moreover, the continuity of $B$ can be proved by following the argumentation in [12, proof of Lemma 2.5].

Since $\zeta$ does not vanish identically almost everywhere on any electrode (cf. (4)), there exist open subsets $e_{m}\subset E_{m}$, $m=1,\dots,M$, of nonzero measure and a constant $\zeta_{-}>0$ such that

[TABLE]

To deduce the coercivity, note first that

[TABLE]

due to the triangle inequality. Hence, by the trace theorem,

[TABLE]

where the second step is a consequence of the Poincaré inequality. The assertion now follows by combining (2) with the estimate

[TABLE]

that is induced by (2), (4) and (8). ∎

Theorem 2.2**.**

Under the assumptions (2) and (4), the problem (3) has a unique solution $(u,U)\in\mathcal{H}^{1}$ for any current pattern $I\in\mathbb{C}_{\diamond}^{M}$. Moreover,

[TABLE]

where $C=C(\Omega,E,\sigma,\zeta)>0$ is independent of $I$.

Proof.

Since $I\in\mathbb{C}_{\diamond}^{M}$ has vanishing mean, the map

[TABLE]

is well-defined (and antilinear). Moreover, since

[TABLE]

the claim follows by applying Lemma 2.1 and the Lax–Milgram theorem to the variational formulation (6). ∎

If one resorts to the standard formulation of the CEM and replaces the second condition of (3) by (5) with contact impedances that are bounded away from infinity, the highest Sobolev regularity that the interior potential can, in general, exhibit is $u\in H^{2-\epsilon}(\Omega)$ due to the abrupt change of the boundary condition from Robin to homogeneous Neumann at $\partial E$ (see, e.g., [4]). With this in mind, an intriguing property of (3) is that the smoothness of the interior electromagnetic potential is directly controlled by the smoothness of the contact admittance $\zeta:\partial\Omega\to\mathbb{C}$, assuming that the boundary $\partial\Omega$ and the admittivity $\sigma$ are smooth enough. We start with two simple lemmas:

Lemma 2.3**.**

Assume that $\partial\Omega$ is of class $C^{\infty}$, $\Gamma\subset\partial\Omega$ is open with a Lipschitz boundary or $\Gamma=\partial\Omega$, and $\eta\in H^{s}(\Gamma)$ for some $s>(n-1)/2$. Then the multiplication operator

[TABLE]

is bounded for any $-s<r\leq s$. More precisely,

[TABLE]

where $C=C(s,r,\Gamma)>0$ is independent of $\eta$ and $v$.

Proof.

The claim follows from, e.g., [26, p. 190, Theorem 1, (i)] by choosing $s_{1}=r$, $s_{2}=s$, $p=q=q_{1}=q_{2}=2$; see also [26, Propositions on p. 14 and p. 150]. ∎

Notice that the assumption $s>(n-1)/2$, which recurs many times in the following, ensures that $H^{s}(\partial\Omega)$ is continuously embedded in the Banach space of continuous functions $\mathcal{C}(\partial\Omega)$ by virtue of the Sobolev embedding theorem.

Lemma 2.4**.**

Assume that $\partial\Omega$ is of class $\mathcal{C}^{\infty}$. Then, for any $(w,W)\in\mathcal{H}^{r+1/2}$ with $r>0$,

[TABLE]

where $C=C(r,E,\Omega)>0$.

Proof.

We may estimate as follows:

[TABLE]

where the last step follows from the trace theorem and the equivalence of norms on a finite-dimensional space. ∎

The following theorem is the main result of this section. Take note that one could also prove a version where the required regularity of $\sigma$ and $\partial\Omega$ depends on the smoothness of $\zeta$, i.e. on $s$, but such a generalization would demand extra work without adding anything to the intuitive contents of the result. To be more precise, the proof would be analogous, but one would have to keep carefully track of how much regularity is required of $\partial\Omega$ and $\sigma$ to guarantee a certain smoothness for $u$.

Theorem 2.5**.**

Assume that (2) and (4) hold, $\sigma\in\mathcal{C}^{\infty}(\overline{\Omega},\mathbb{C}^{n\times n})$ and $\partial\Omega$ is of class $\mathcal{C}^{\infty}$. If furthermore $\zeta\in H^{s}(\partial\Omega)$ for some $s>(n-1)/2$, then the solution to (6) satisfies

[TABLE]

where $C=C(\Omega,E,\sigma,\zeta,s)>0$ is independent of $I\in\mathbb{C}_{\diamond}^{M}$.

Proof.

To begin with, consider the Neumann boundary value problem

[TABLE]

for an arbitrary $f$ in the mean-free Sobolev space

[TABLE]

where $\langle\,\cdot\,,\,\cdot\,\rangle_{\partial\Omega}:H^{r}(\partial\Omega)\times H^{-r}(\partial\Omega)\to\mathbb{C}$ denotes the bilinear dual evaluation between Sobolev spaces on $\partial\Omega$. Due to the infinite smoothness of $\sigma$ and $\partial\Omega$, for any $r\in\mathbb{R}$ the problem (11) has a unique solution $v\in H^{r+3/2}(\Omega)/\mathbb{C}$ satisfying [23]

[TABLE]

In particular, the Neumann-to-Dirichlet map

[TABLE]

is bounded for any $r\in\mathbb{R}$ by virtue of trace theorems for those elements of $H^{r+3/2}(\Omega)$ for which the range of the differential operator $\nabla\cdot\sigma\nabla(\,\cdot\,)$ is a subspace of $L^{2}(\Omega)$ [23]. To complete this introductory part of the proof, let

[TABLE]

be a linear and bounded extension operator (with respect to any given norm for the finite-dimensional space ${\rm span}\{\chi_{1},\dots,\chi_{M}\}$), i.e., such that $(\mathcal{E}\eta)|_{E}=\eta|_{E}$ for all $\eta\in{\rm span}\{\chi_{1},\dots,\chi_{M}\}$.

Let $(u,U)\in\mathcal{H}^{1}$ be the solution to (3) for some $I\in\mathbb{C}_{\diamond}^{M}$. The second equation of (3) immediately leads to the preliminary estimate

[TABLE]

by virtue of Lemmas 2.3 and 2.4, and Theorem 2.2.

Set $r=\min\{1,s\}$ and denote by $(\tilde{u},\tilde{U})$ the particular representative of the equivalence class $(u,U)\in\mathcal{H}^{1}$ for which $\tilde{u}|_{\partial\Omega}\in H^{1/2}_{\diamond}(\partial\Omega)$. Since $\zeta$ is supported on $\overline{E}$ by assumption (4), the second condition of (3) can be rewritten in an alternative form

[TABLE]

Because of Lemma 2.3, the vanishing mean of $\tilde{u}$, the boundedness of the operators $\Lambda_{\sigma}$ and $\mathcal{E}$, and (14), we thus have

[TABLE]

Moreover, since $\zeta\in H^{s}(\partial\Omega)\subset\mathcal{C}(\partial\Omega)$ is nonzero on some set of nonzero measure on each $E_{m}$, $m=1,\dots,M$, due to (4), the second equation of (3) leads to

[TABLE]

where $C>0$ does not depend on $I$ (or $\tilde{U}$) and the last step is a consequence of (14) and the boundedness of $\Lambda_{\sigma}:L^{2}_{\diamond}(\partial\Omega)\to L^{2}(\partial\Omega)/\mathbb{C}$. Combining the previous estimate with (12) and (2), results in

[TABLE]

In particular, it is straightforward to deduce that

[TABLE]

which proves the claim if $s\leq 1$.

If $s>1$, one can repeat the above argument with $r=\min\{2,s\}$, using the previous (16) in place of (14) in the new (2), to end up once again with (17), which this time around corresponds to (10) if $s\leq 2$. The complete assertion follows by inductively reiterating this argument. ∎

Remark 2.6**.**

The assumption on the smoothness of the boundary $\partial\Omega$ can be relaxed in many ways. For example, if $n=2$, the domain is a square and the corners are well separated from the electrodes, one can locally straighten the corners by introducing a suitable conformal map ($z\mapsto z^{2}$). The crucial property is the vanishing Neumann boundary condition that is preserved when the domain is deformed.

3. Shape derivatives

As an example of a setting where the higher regularity of the smoothened CEM may be useful, we consider the Fréchet derivative of electrode measurements with respect to the shape of $\Omega$. Such a derivative is needed in absolute EIT if the body shape is not accurately known (cf., e.g., medical imaging) and must be reconstructed at the same time as the admittivity to avoid severe artifacts [6, 7, 24]. For simplicity and to be able to use the results in [5, 6] without further generalizations, we assume throughout this section that $\partial\Omega$ and $\partial E$ are smooth.

The measurement, or current-to-voltage map of the CEM is usually defined via

[TABLE]

where $U$ is the second component of the solution to (3) and $\mathbb{C}^{M}/\mathbb{C}=\{\{V+c\mathbf{1}\,|\,c\in\mathbb{C}\}\,|\,V\in\mathbb{C}^{M}\}$. Due to a certain symmetry of (6), $R$ can be represented as a real-symmetric, complex $(M-1)\times(M-1)$ matrix in terms of any orthonormal basis for $\mathbb{C}^{M}_{\diamond}\sim\mathbb{C}^{M}/\mathbb{C}$. However, assuming that the admittivity $\sigma$ is defined in some neighborhood of $\Omega$, $R$ can be interpreted as a function of two variables,

[TABLE]

where $\mathcal{B}_{d}\subset[\mathcal{C}^{1}(\partial\Omega)]^{n}$ is an origin-centered open ball of radius $d>0$. Moreover, $(u(I,h),U(I,h))$ is the solution of (3) when $\partial\Omega$ and $E_{m}$, $m=1,\dots,M$, are replaced by the perturbed versions

[TABLE]

respectively. The surface admittance $\zeta$ is assumed to stretch accordingly, i.e., the electrode contacts on $\partial\Omega^{h}$ are characterized by the function

[TABLE]

If $d>0$ is chosen small enough, the above definitions are unambiguous in the sense that $\partial\Omega^{h}$ defines a $\mathcal{C}^{1}$ boundary of a bounded domain $\Omega^{h}$, $\{E^{h}_{m}\}_{m=1}^{M}$ is a proper set of connected, well-separated electrodes and $\zeta^{h}$ defines an element of $L^{\infty}(\partial\Omega^{h})$ satisfying (4) with respect to the perturbed boundary and electrodes [6]. In the following, we will implicitly assume that $d>0$ is chosen small enough in this sense.

For the traditional CEM, the elliptic boundary value problem defining the Fréchet derivative of $R:\mathbb{C}^{M}_{\diamond}\times\mathcal{B}_{d}\to\mathbb{C}^{M}/\mathbb{C}$ with respect to its second variable falls outside the $H^{1}$-based variational theory; to be more precise, the solution of the shape derivative problem belongs to $\mathcal{H}^{1-\epsilon}$ for any $\epsilon>0$ [6, Theorem 3.3]. However, for the smoothened version (3), the problem defining the shape derivative admits a $\mathcal{H}^{1}$-based variational formulation with the same left-hand side as in (6), if $\zeta$ is regular enough.

Assume that $\zeta\in H^{s}(\partial\Omega)$ with some $s>(n-1)/2$ and let $(u,U)\in\mathcal{H}^{s+3/2}$ be the solution of (3) for some $I\in\mathbb{C}_{\diamond}^{M}$. Consider the problem of finding $(u^{\prime}[h],U^{\prime}[h])\in\mathcal{H}^{1}$ such that

[TABLE]

where the antilinear functional $F[h]:\mathcal{H}^{1}\to\mathbb{C}$ is defined by

[TABLE]

Here, ${\rm Grad}$ denotes the surface gradient [3], $\kappa:\partial\Omega\to\mathbb{R}$ is the sum of principal curvatures on $\partial\Omega$, and $h_{\nu}=h\cdot\nu$ and $h_{\tau}=h-h_{\nu}\nu$ are the normal (scalar) and tangential (vector) components of $h:\partial\Omega\to\mathbb{R}^{n}$, respectively. Moreover, $\langle\,\cdot\,,\,\cdot\,\rangle_{E}:H^{-r}(E)\times H^{r}(E)\to\mathbb{C}$, $|r|<1/2$, denotes the bilinear dual evaluation between Sobolev spaces on the electrodes (see, e.g., [23]). Finally, note that $F[h]$ depends linearly on $h$, and so the same conclusion also holds for the solution $(u^{\prime}[h],U^{\prime}[h])\in\mathcal{H}^{1}$, if it uniquely exists.

Lemma 3.1**.**

Assume that (2) and (4) hold, $\sigma\in\mathcal{C}^{\infty}(\overline{\Omega},\mathbb{C}^{n\times n})$, and $\partial\Omega$ and $\partial E$ are of class $\mathcal{C}^{\infty}$. If furthermore $\zeta\in H^{s}(\partial\Omega)$ for some $s>(n-1)/2$, then the variational problem (20) has a unique solution $(u^{\prime}[h],U^{\prime}[h])\in\mathcal{H}^{1}$ for any $h\in[\mathcal{C}^{1}(\partial\Omega)]^{n}$ and $I\in\mathbb{C}_{\diamond}^{M}$, satisfying

[TABLE]

where $C>0$ is independent of $h$ and $I$.

Proof.

By virtue of Lemma 2.1 and the Lax–Milgram theorem, it is enough to prove that

[TABLE]

for any $h\in[\mathcal{C}^{1}(\partial\Omega)]^{n}$ and $I\in\mathbb{C}_{\diamond}^{M}$. By resorting to the triangle inequality, we may handle the terms on the right-hand side of (3) one by one.

According to the Cauchy–Schwarz inequality and Lemma 2.4,

[TABLE]

Applying the triangle inequality and Lemma 2.3 (with $\eta=\zeta|_{E}\in H^{s}(E)$, $\eta=h_{\nu}|_{E}\in\mathcal{C}^{1}(\overline{E})\hookrightarrow H^{1}(E)$ and $\eta=\kappa|_{E}\in\mathcal{C}^{\infty}(\overline{E})$) to the first factor on the right-hand side gives

[TABLE]

where the last two steps follow from the trace theorem, Lemma 2.4 and Theorem 2.5. This takes care of the first term on the right-hand side of (3).

Set $\delta=\min\{1/2,(s-1/2)/2\}$, so that $0<\delta\leq 1/2$ and $1/2+\delta<s$. Lemmas 2.3 and 2.4 yield

[TABLE]

where the last step corresponds to a weaker version of Theorem 2.5. Moreover, by the continuity of ${\rm Grad}:H^{s}(E)\supset H^{1/2+\delta}(E)\to[H^{-1/2+\delta}(E)]^{n}$ (cf. [23, p. 85, Prop. 12.1]) and Lemma 2.3 with $\eta=[h_{\tau}|_{E}]_{j}\in\mathcal{C}^{1}(\overline{E})\hookrightarrow H^{1}(E)$, $j=1,\dots,n-1$,

[TABLE]

In consequence,

[TABLE]

which handles the second term on the right-hand side of (3).

Finally, by using the continuity of ${\rm Grad}:H^{1/2}(\partial\Omega)/\mathbb{C}\to[H^{-1/2}(\partial\Omega)]^{n}$, Lemma 2.3 with $\eta=h_{\nu}\in\mathcal{C}^{1}(\partial\Omega)\hookrightarrow H^{1}(\partial\Omega)$ and the trace theorem,

[TABLE]

where the final step follows from Theorem 2.5. This completes the proof. ∎

As hinted above, the problem (20) defines the Fréchet derivative of the measurement map $R$ with respect to its second variable, that is, with respect to the shape of $\Omega$.

Theorem 3.2**.**

Assume that (2) and (4) hold, $\sigma\in\mathcal{C}^{\infty}(\overline{\Omega},\mathbb{C}^{n\times n})$, and $\partial\Omega$ and $\partial E$ are of class $\mathcal{C}^{\infty}$. If furthermore $\zeta\in H^{s}(\partial\Omega)$ for some $s>(n-1)/2$, then $R:\mathbb{C}^{M}_{\diamond}\times\mathcal{B}_{d}\to\mathbb{C}^{M}/\mathbb{C}$ is Fréchet differentiable with respect to its second variable at the origin, with the Fréchet derivative given by the linear and bounded map

[TABLE]

where $U^{\prime}[h]$ is the second component of the solution to (20).

Proof.

The assertion follows by slightly modifying (or simplifying) the proofs in [5, 6] taking into account that now $(u,U)\in\mathcal{H}^{2}$ unlike in [5, 6]. ∎

If $I,\tilde{I}\in\mathbb{C}^{M}_{\diamond}$ are electrode current patterns and $(u,U),(\tilde{u},\tilde{U})$ are the respective solutions of (3), then $U^{\prime}[h]$ can be assembled using the relation

[TABLE]

where it is assumed that $\zeta\in H^{\max\{1,s\}}(\partial\Omega)$ for some $s>(n-1)/2$. This formula can be straightforwardly deduced from the variational problems (6) and (20): Use the complex conjugate of $(u^{\prime}[h],U^{\prime}[h])$ as the test element in (6) with $I=\tilde{I}$, then interpret the complex conjugate of $(\tilde{u},\tilde{U})$ as the test element in (20), and finally employ the higher smoothness of the pairs $(u,U)$ and $(\tilde{u},\tilde{U})$ guaranteed by Theorem 2.5 to interpret the dual evaluations in (3) as regular integrals (cf. [5, 6]). If $\sigma$ is scalar-valued (i.e., isotropic), the integral (22) can be rewritten as

[TABLE]

which reveals a symmetry between $(u,U)$ and $(\tilde{u},\tilde{U})$, that is, the roles of $(u,U)$ and $(\tilde{u},\tilde{U})$ can be reversed on the right-hand side of (22)–(24) without altering its value.

To complete this section, take note that the sampling formula (22)–(24) also holds if $\zeta$ is less regular, in which case (23) just needs to be interpreted as an appropriate dual evaluation. In particular, if the second equation of (3) is replaced by the conditions (5) of the traditional CEM, the term $h_{\tau}\cdot{\rm Grad}\,\zeta$ becomes a weighted delta distribution supported on $\partial E$ and (23) turns into an integral over $\partial E$ [6, Corollary 3.4].

4. Numerical experiments

The numerical implementation of the proposed smoothened CEM for EIT is a direct generalization of the traditional CEM implementation presented in, e.g., [16]. Merely, the nonhomogeneous admittance $\zeta$ has to be incorporated into the associated integrals. As discussed above, the same also holds true for the numerical shape derivatives, with the exception of (23) that reduces to an integral over the electrode boundaries in the traditional CEM.

In this section, we restrict our attention to isotropic, real-valued material parameters and currents. In particular, $\sigma$ denotes the electrical conductivity and $\zeta$ is the contact conductance. This is the most commonly considered setting in practical EIT: For anisotropic conductivities, the inverse problem of EIT suffers from an inherent nonuniqueness [29] and, on the other hand, the temporal frequencies employed by EIT devices are often so low that one can ignore the capacitive effects, which leads to a real-valued $\sigma$ (see [31] and references therein). Be that as it may, some of the conclusions drawn in what follows may not apply as such to anisotropic or/and complex conductivities. We only work in two dimensions for the sake of computational simplicity.

Our numerical experiments are divided into four parts: First, we investigate the discrepancy between electrode measurements predicted by the traditional CEM and by a particular version of the smoothened CEM (see Figure 1). Next, we study the convergence of FEM for the two forward models as well as for the shape derivative integrals (22)–(24). The third test compares predictions of the two models to experimental EIT data from a homogeneous water tank. Finally, using experimental data corresponding to an insulating inclusion in the aforementioned water tank, we demonstrate that the quality of a reconstruction produced by a Bayesian algorithm, which aims at finding a maximum a posteriori (MAP) estimate for the conductivity and the contact conductance parameters, is independent of the choice between the two models.

4.1. Model differences

In order to study the properties of the smoothened CEM, we set up a test geometry on the unit square $\Omega=(0,1)^{2}$ with $M=8$ electrodes as shown on the left in Figure 1. A square is our choice of test domain because it allows an exact representation as a union of triangles, and thus the discretization of the domain boundary does not induce any extra error in the considered FEM solutions; cf. Remark 2.6. In what follows, the numerical solutions for the elliptic problem (6) are computed by FEM on uniform meshes similar to the one shown in Figure 1. The number of nodes ranges from $9^{2}$ to $2049^{2}$ and in each case the electrodes cover an even number of element edges. Unless otherwise stated, the discretization is piecewise linear so that the number of nodes equals the number of degrees of freedom. The width of a boundary edge element is denoted by $h>0$.

Recall that the standard CEM is recovered from our setting by defining $\zeta$ to be constant on each electrode and zero elsewhere, cf. (5). Let us replace these box-shaped functions by hat functions with an equal area under their graphs as depicted on the right in Figure 1. Such a conductance function satisfies $\zeta\in H^{3/2-\epsilon}(\Gamma)$ for any $\epsilon>0$ and $\Gamma\subset\partial\Omega$ that does not contain any of the four corners. According to Theorem 2.5 and Remark 2.6, the corresponding solution $(u,U)$ thus belongs to $\mathcal{H}^{3-\epsilon}$ provided that the conductivity $\sigma$ is smooth enough.

Before considering convergence of FEM or real measurements, a natural question to ask is how the electrode potentials produced by the new model differ from the standard CEM when both forward problems are solved numerically. To this end, we choose a finite element mesh with $1025^{2}$ nodes and define the relative difference

[TABLE]

where $U^{(m)}$ and $U^{(m)}_{\text{CEM}}$ denote the electrode potential vectors for the hat model and the standard model, respectively, corresponding to a common current pattern $I^{(m)}=\mathrm{e}_{M}-\mathrm{e}_{m}\in\mathbb{R}^{M}_{\diamond}$. Here, $\mathrm{e}_{m}$ denotes the $m$th Cartesian basis vector of $\mathbb{R}^{M}$. The unit or the magnitude of the current are not relevant since they cancel out due to linearity. It is assumed that the conductivity $\sigma>0$ is constant and that the contact conductance, characterized by the height of the box function (or the half-height of the hat function), is the same for each electrode. We denote this conductance value by $\zeta_{\text{el}}>0$. Now the relative difference (25) actually depends only on the ratio $\sigma/\zeta_{\text{el}}\in\mathbb{R}$, as can be easily deduced from (3). The quantity $\sigma/\zeta_{\text{el}}$ has the unit of length (e.g., meters).

It turns out that the relative difference $d_{U}$ is around $9\%$ for some values of $\sigma/\zeta_{\text{el}}$, but it clearly diminishes if the conductance value is either very high or very low. This behavior is illustrated in the top left image of Figure 2. The next step is to study whether the discrepancy between the models can be kept smaller if the hat functions are scaled differently than in Figure 1. Given a ratio $\sigma/\zeta_{\text{el}}$ for the standard CEM, we thus try to find $\zeta_{\text{el}}^{\prime}>0$, the half-height of $\zeta$ in our smoothened model, that minimizes the corresponding difference $d_{U}^{\prime}(\sigma,\zeta_{\text{el}},\zeta_{\text{el}}^{\prime})$, which is computed as (25) but now with this new scaling for the hat function. For given values of $\sigma$ and $\zeta_{\text{el}}$, finding $\zeta_{\text{el}}^{\prime}$ is a well-defined one-dimensional minimization problem. The plot at the bottom in Figure 2 presents the dependence between $\zeta_{\text{el}}$ and $\zeta_{\text{el}}^{\prime}$, whereas the top right image shows the smallest relative difference $d_{U}^{\prime}$ as a function of $\sigma/\zeta_{\text{el}}$. The difference attains its maximum value of about $5.8\cdot 10^{-3}$ at $\sigma/\zeta_{\text{el}}\approx 50\cdot 10^{-3}\,\mathrm{m}$ and clearly decreases when the ratio is changed to either direction. Some numerical instabilities occur when the ratio gets very low, that is, when the setting approaches the shunt model [8].

When dealing with real-world measurements, the scaling of hat functions is most likely not an issue since the values of the contact conductances are usually not interesting per se, but they are merely estimated in order to obtain a more accurate reconstruction for the conductivity. Moreover, there is no reason to expect that the standard CEM with piecewise continuous conductances would be a totally accurate measurement model either, which is why the quantity (25) is not called an “error” but a “difference”.

In [2], the ratio $\sigma/\zeta_{\text{el}}\approx 2.4\cdot 10^{-3}\,\mathrm{m}$ was reported in several saline experiments. Regarding EIT with tap water environment, the prior means used in [7] correspond to $1.3\cdot 10^{-3}\,\mathrm{m}$, whereas in [13] the initial guess for the conductivity and the contact conductance results in $\sigma/\zeta_{\text{el}}\approx 0.7\cdot 10^{-3}\,\mathrm{m}$. In [14], the permissible ratio varies between $1.1\cdot 10^{-5}\,\mathrm{m}$ and $1.1\cdot 10^{-3}\,\mathrm{m}$. All these experiments were performed with cylindrically symmetric water tanks with circumferences of about $1\,\mathrm{m}$, i.e., roughly one fourth of the circumference of our unit square. Thus, the corresponding points in Figure 2 can be found by multiplying these ratios by 4, assuming that our observations can be generalized to cases where the geometry and the number of electrodes differ from our test setup. In any case, it seems that the lowest reported values for $\sigma/\zeta_{\text{el}}$ are well below the peak in the top right plot of Figure 2, actually approaching the shunt model. However, one cannot exclude the possibility of encountering conductivity-conductance ratios that correspond to the largest mismatch between the two models because the contact conductances can vary significantly, e.g., in EIT imaging of concrete [17]. On the other hand, typical measurement noise levels are clearly nonnegligible compared to any value on the modeling difference graph in the top right image of Figure 2; cf., e.g., [19] where the performance of the EIT unit used in Sections 4.3 and 4.4 is analyzed.

4.2. Convergence of FEM

Let us continue using the test setup of Figure 1 for both the standard CEM and for the hat function conductance model. We study the convergence of finite element solutions toward “exact” reference solutions that are computed for both electrode models by using the finest mesh with $2049^{2}$ nodes. In addition to the piecewise linear FEM basis functions, the quadratic Lagrange elements are used for comparison; for general information on properties and advantages of FEMs of different types and order, we refer to the textbook [20] and the references therein. Two pairs of conductivities and contact conductances are considered: the first one corresponds to $\sigma/\zeta_{\text{el}}\approx 50\cdot 10^{-3}\,\mathrm{m}$, i.e., the peak in the top right plot in Figure 2, whereas the second pair results in $\sigma/\zeta_{\text{el}}\approx 4\cdot 10^{-3}\,\mathrm{m}$, which roughly corresponds to the values mentioned in [7]. The conductance half-heights $\zeta_{\text{el}}^{\prime}$ for the smoothened model are computed as in the bottom plot of Figure 2 and they are $\sigma/\zeta_{\text{el}}^{\prime}\approx 30\cdot 10^{-3}\,\mathrm{m}$ and $\sigma/\zeta_{\text{el}}^{\prime}\approx 0.5\cdot 10^{-3}\,\mathrm{m}$, respectively. The relative errors are computed in the same way as in (25), but now separately for each conductance model and comparing against the respective reference solution.

The left-hand plot in Figure 3 demonstrates that when using piecewise linear finite elements, the smoothened model converges faster toward its reference solution, but the difference cannot be considered significant. Both models seem to exhibit an asymptotic decay rate of $h^{2}$ or a bit less. For the standard CEM, this is in line with the material in [8], where the rate $h^{2-\epsilon}$, $\epsilon>0$, is predicted; for the smoothened hat-CEM, one would expect to gain the aforementioned $\epsilon>0$ that is due to the solution of the standard CEM only lying in $\mathcal{H}^{2-\epsilon}$. On the other hand, there is a big difference between the two conductance values: the further the values are from the shunt model (i.e., from $\zeta=\infty$), the more accurate the numerical solutions are for a given mesh. This observation is also in tune with [8]. The dashed curves, which correspond to the values in [7], clearly show the numerical difficulties that appear when approaching the shunt model.

There is a greater difference between the models when quadratic elements are used. The right-hand plot in Figure 3 indicates that the asymptotic convergence rate of the smoothened model increases to approximately $h^{3}$, whereas the rate for the standard model stays at about $h^{2}$. In fact, an extrapolation of the argumentation in [8] suggests that the asymptotic rates should be $h^{3-\epsilon}$ and $h^{2-\epsilon}$, respectively, with the smoothness of the forward solutions ($\mathcal{H}^{3-\epsilon}$ and $\mathcal{H}^{2-\epsilon}$, respectively) imposing in both cases the upper bound on the speed of convergence. In particular, it is to be expected that the use of, say, third order elements would improve the asymptotic convergence rate only if the model for the contact conductance were also further smoothened; see Theorem 2.5 and [20]. Once again, more accurate solutions are obtained for the higher ratio $\sigma/\zeta_{\text{el}}$, that is, when the contacts are worse and the setting is further away from the shunt model (cf. [8]).

To confirm the above findings about the convergence of FEM, we repeat the numerical experiment with the constant conductivity replaced by the phantom shown on left in Figure 4. This time there are twelve electrodes attached to the boundary of the unit square, and the contact conductances shown on right in Figure 4 are chosen randomly so that their ratios with the mean of the conductivity phantom cover approximately the whole scale on the horizontal axes of the images in Figure 2. Figure 5 presents the corresponding results. It is organized in the same way as Figure 3, that is, the left-hand image illustrates the convergence for the two models with piecewise linear FEM, whereas the right-hand image considers quadratic FEM. The conclusions are the same as for the first experiment: For piecewise linear basis functions both models exhibit convergence rates of approximately $h^{2}$, with the smoothened CEM converging slightly faster asymptotically. For the quadratic basis functions, the smoothened model clearly prevails with a decay rate $h^{3}$ compared to $h^{2}$ for the traditional model.

The difference between the standard CEM and the new smoothened version becomes apparent even for the piecewise linear FEM when studying the convergence of the numerical Fréchet derivatives. To demonstrate this, we consider the integrals (22)–(24) with the choices $\sigma=1$ and $h_{\nu}=\|h_{\tau}\|_{2}=1$ as well as $\kappa|_{E}=0$ dictated by the geometry. More precisely, we numerically evaluate the integrals

[TABLE]

based on approximations of the solutions to (3) on different finite element meshes. Here $\dot{\zeta}$ denotes the derivative of the contact conductance with respect to the arclength parameter. Notice that $\dot{\zeta}$ becomes a linear combination of delta distributions and $\mathcal{I}_{2}^{(m,n)}$ a linear combination of pointwise evaluations of the integrand $(U^{(m)}-u^{(m)})(U^{(n)}-u^{(n)})$ at the end points of the electrodes when the traditional CEM is considered [6, 7].

We only consider the setup with eight electrodes depicted in Figure 1 as well as a constant conductivity and two sets of identical contact conductances so that either $\sigma/\zeta_{\text{el}}\approx 50\cdot 10^{-3}\,\mathrm{m}$ or $\sigma/\zeta_{\text{el}}\approx 4\cdot 10^{-3}\,\mathrm{m}$, which correspond to $\sigma/\zeta_{\text{el}}^{\prime}\approx 30\cdot 10^{-3}\,\mathrm{m}$ or $\sigma/\zeta_{\text{el}}^{\prime}\approx 0.5\cdot 10^{-3}\,\mathrm{m}$, respectively (cf. the bottom image of Figure 2). In other words, the parameter choices are as in Figure 3. The corresponding results for the inhomogeneous phantom and the contact conductances of Figure 4 would be analogous. The employed current patterns and their correspondence to the forward solutions of (3) in (26)–(28) are as before, and the relative errors are computed via

[TABLE]

where tilde denotes an integral computed using the appropriate reference solutions on the densest mesh with $2049^{2}$ nodes. The results for $\delta_{i},i=1,2,3$ are plotted in Figure 6, which demonstrates that the smoothened model is considerably more accurate than the traditional CEM, except for very coarse meshes. In addition, the difference between the two conductance values is again clearly visible.

4.3. Comparison to experimental EIT data

The compatibility of the proposed smoothened model with real measurements is studied with data from a water tank shown on left in Figure 7. There are 16 electrodes of width $2\,\mathrm{cm}$ attached to the interior lateral surface of the tank, extending from the bottom all the way up to the water surface. The circumference of the tank is $106\,\mathrm{cm}$. The actual measurements were performed with low-frequency (1 kHz) alternating current using the Kuopio impedance tomography (KIT4) device [19]. The phase information is ignored and the data interpreted as if resulting from the use of direct current, which is reasonable due to the insignificance of capacitive effects at low temporal frequencies [31], as explained at the beginning of Section 4. The geometry, which is essentially two-dimensional because the electrodes are of the same height as the water layer (see, e.g., [21]), is known up to unavoidable mismodeling due to small imperfections in the construction of the tank. The domain is discretized for the potential with triangular elements having $12\,218$ nodes. Each electrode is divided into approximately $10$ element edges.

We make the reasonable assumption that the conductivity of the water layer within the tank is homogeneous and aim at reconstructing its value and the contact conductance parameters by solving the minimization problem

[TABLE]

where $y\in\mathbb{R}_{+}^{1+16}$ represents the unknowns. The vector $\tilde{\mathcal{U}}\in\mathbb{R}^{240}$ contains all measured electrode potentials and the values in $\mathcal{U}(y)$ are the computed potentials for the employed fifteen linearly independent current patterns and a given parameter vector $y$. For the traditional CEM, the latter sixteen components of $y$ are the contact conductances, whereas for the smoothened CEM they represent the half-heights of the hat-shaped conductance functions.

The employed minimization algorithm is based on Levenberg–Marquardt method. The initial guesses for the constant conductivity and contact conductances were chosen as $\sigma=0.25\,\mathrm{mS}/\mathrm{cm}$ and $\zeta_{\text{el}}=10\,\mathrm{mS}/\mathrm{cm}^{2}$, respectively, for the traditional CEM. For the smoothened CEM, the values $\sigma=0.25\,\mathrm{mS}/\mathrm{cm}$ and $\zeta_{\text{el}}^{\prime}=70\,\mathrm{mS}/\mathrm{cm}^{2}$ were used. The chosen initial values for $\sigma$ and $\zeta_{\text{el}}$ are close to the ones used in [13] and the ratio $\zeta_{\text{el}}/\zeta_{\text{el}}^{\prime}$ agrees with the appropriate value on the graph in the bottom image of Figure 2. In practice, one need not have such a graph available since the minimization process does not seem to be sensitive to the initial guess for $\zeta_{\text{el}}^{\prime}$; this conclusion also applies to the considerations in Section 4.4. The three-dimensional quantities for the conductivity and the conductance can be converted to their two-dimensional counterparts by multiplying with the height of the tank which is $5\;\mathrm{cm}$. For both models, the minimization algorithm converged without any complications.

The optimal values for the constant conductivity level are $\sigma=0.22722\,\mathrm{mS}/\mathrm{cm}$ for the traditional model and $\sigma=0.22720\,\mathrm{mS}/\mathrm{cm}$ for the smoothened one. The corresponding contact conductance parameters for the two models are shown on right in Figure 7. The results are well in line with the bottom image of Figure 2: As the contacts are good in water tank experiments, i.e., the ratio $\sigma/\zeta$ is low, one expects the half-heights of the hat-like conductances to be substantially higher than the constant contact conductances for the traditional CEM. The minimal relative discrepancies between the measurements and the two models, i.e.,

[TABLE]

where $y^{*}\in\mathbb{R}^{17}$ is the respective solution of (29), was $1.20\%$ for the traditional CEM and $1.21\%$ for the smoothened CEM with hat-like conductances. If one considers the discrepancy in comparison to the maximal variation in the measurements, i.e., replaces $\lVert\tilde{\mathcal{U}}\rVert_{2}$ in the denominator by $\sqrt{240}\,\max_{1\leq i,j\leq 240}\lvert\tilde{\mathcal{U}}_{i}-\tilde{\mathcal{U}}_{j}\rvert$, the numbers are $0.220\%$ for the traditional model and $0.221\%$ for the smoothened one. In particular, according to this single experiment, the two models seem to be in approximately as good accordance with real-world EIT data.

The relative discrepancies listed above are somewhat higher than the expected noise level in the data (cf. [19]). Part of this extra mismatch probably originates from small errors in the model for the measurement configuration; absolute EIT is known to be extremely sensitive to geometric mismodeling [18]. In addition, treating the data as if it originated from direct current measurements may have a small effect on the results. In particular, the evidence about the validity of the smoothened model presented here is not conclusive, but it needs to be confirmed by more carefully designed experimental studies in the future.

4.4. Reconstructions from experimental EIT data

Let us next consider the same water tank as in the previous section, but this time with one embedded insulating cylindrical inclusion of radius $3.5\,{\rm cm}$ made out of plastic; see the top image in Figure 8. This time the (two-dimensional) domain is discretized for the potential and the conductivity with triangular elements having $8004$ nodes. Each electrode is again divided into approximately $10$ element edges. Our goal is to compute reconstructions of the conductivity phantom by using the traditional and the smoothened CEM and to demonstrate that the choice between the two models does not have a significant effect on the outcome. In particular, our objective is not to exploit the faster convergence of FEM approximations for the smoothened model by employing sparser meshes or/and higher order elements, but such considerations are left for future studies.

The reconstruction algorithm aims at computing a MAP estimate

[TABLE]

for the unknown parameter vector $y\in\mathbb{R}_{+}^{8004+16}$ that represents the discretized conductivity field and the contact conductance values and whose expected value is $y_{0}$. As in the previous section, the vector $\tilde{\mathcal{U}}\in\mathbb{R}^{240}$ contains all measured electrode potentials and the values in $\mathcal{U}(y)$ are the computed potentials for the employed fifteen linearly independent current patterns and a given parameter vector $y$. The positive semidefinite matrix $G$ originates from combining a prior distribution for the conductivity with the assumption that each measured electrode potential is corrupted by an independent realization of a normally distributed random variable with zero mean and standard deviation $2\cdot 10^{-3}\max_{1\leq i,j\leq 240}\lvert\tilde{\mathcal{U}}_{i}-\tilde{\mathcal{U}}_{j}\rvert$, the choice of which is motivated by the relative discrepancies listed in the previous section. More precisely, it is formally assumed that the conductivity is a priori a Gaussian random field with correlation length $4\,\mathrm{cm}$, pointwise standard deviation $0.25\,\mathrm{mS}/\mathrm{cm}$, and a constant expectation function $0.25\,\mathrm{mS}/\mathrm{cm}$ that is close to the conductivity of Finnish tap water. The correlation length is intentionally chosen to be of the same order as the radius of the inclusion, i.e., it corresponds to the size of inhomogeneities we expect to find inside the tank. Choosing a longer correlation length would blur the reconstruction of the inclusion, whereas employing a significantly shorter correlation length would lead to unwanted oscillations in the reconstruction of the background conductivity level. The contact conductances are estimated without further prior knowledge and thus the corresponding parts of $G$ and $y_{0}$ are empty. The employed minimization algorithm is again based on Levenberg–Marquardt method. We refer to, e.g., [7, 15] for more details on Bayesian inversion.

The initial values for the conductivity and the contact conductance parameters are the same as in the previous section. In particular, the initial guess/expected value for the conductivity is homogeneous. The reconstruction corresponding to the traditional CEM is presented on the bottom left in Figure 8 and the one corresponding to the smoothened model with hat-like conductances on the bottom right. The results are remarkably similar, which is what we have observed with other water tank experiments as well. In fact, the reconstruction obtained by using the smoothened model is practically indistinguishable from the one corresponding to the traditional CEM: their relative discrepancy in the $L^{2}(\Omega)$ norm is just $0.3\%$. Both reconstructions also clearly indicate the location of the insulating inclusion as a region of almost vanishing conductivity and the background conductivity level is close to the values obtained in the previous section. Thus, the proposed model seems to be as compatible with experimental measurements as the traditional CEM, at least in the examined water tank setting where the contact conductances are relatively high. The reason may simply be that the difference between the two models is quite small when the conductance functions are scaled properly, as demonstrated by the top right plot in Figure 2. On the other hand, it cannot be ruled out that the real physical phenomena at the electrode contacts may even be more accurately described by some variant of the proposed smoothened CEM.

5. Concluding remarks

We have introduced a smoothened version of the CEM for EIT. The new model retains the essential solvability and differentiability properties, while the regularity of the solution can be arbitrarily improved by choosing an appropriate smoothened conductance for the electrode contacts. It was numerically demonstrated that at least the simplest, piecewise linear smoothening is almost equivalent to the standard CEM, if the hat-shaped conductance functions are scaled properly. What is more, the presented EIT reconstruction with the proposed model is almost indistinguishable from the one obtained with the traditional CEM. The computational feasibility of the new model was also demonstrated by the superior convergence of the FEM.

To summarize, we recommend using some variant of the smoothened CEM especially if one wants to exploit a higher-order FEM solver or needs to compute numerical shape derivatives. However, even the proposed smoothened model cannot overcome the numerical problems related to very high contact conductances that appear when the measurement setting approaches the so called shunt model.