Reconstructing conductivities in three dimensions using a - - PowerPoint PPT Presentation

Reconstructing conductivities in three dimensions using a non-physical scattering transform

Kim Knudsen Department of Mathematics Technical University of Denmark Joint with Jutta Bikowski and Jennifer Mueller, Colorado State University AIP2009 Vienna July 24

Direct reconstruction in 3D 1/26

Outline

• 1. Introduction to the problem
• 2. Calderon’s approximate reconstruction algorithm
• 3. Exact high complex frequency reconstruction algorithm
• 4. Implementation details and numerical results

Direct reconstruction in 3D 2/26

• 1. Introduction to the problem

Direct reconstruction in 3D 3/26

The conductivity equation

Smooth domain Ω ⊂ R3; conductivity γ ∈ L∞(Ω), C−1 ≤ γ ≤ C, for C > 0. A voltage potential u in Ω : ∇ · γ∇u = 0 in Ω. Measurements:

ν f, g ∂Ω Ω γ(x)

Apply voltage potential: f = u|∂Ω Measure current flux : g = γ∂νu|∂Ω. Dirichlet to Neumann (Voltage to Current) map Λγ : f → g, weakly defined by Λγf, h =

• ∂Ω

(Λγf)hdσ(x) =

• Ω

γ∇u · ∇vdx, v|∂Ω = h ∈ H1/2(∂Ω).

Direct reconstruction in 3D 4/26

Inverse problem

Consider the (non-linear) mapping Λ: γ → Λγ. This mapping encodes the direct problem. The Calder´

• n problem (the inverse conductivity problem):
• Uniqueness: Is Λ injective?
• Reconstruction: How can γ be computed from Λγ?

Direct reconstruction in 3D 5/26

Transformation to Schr¨

• dinger equation

Suppose u solves ∇ · γ∇u = 0 in Ω, u|∂Ω = f. Then v = γ−1/2u solves (∆ + q)v = 0 in Ω v|∂Ω = γ−1/2f, with q = −∆γ1/2/γ1/2 ⇔ (∆ + q)γ1/2 = 0. Dirichlet to Neumann map Λqf = ∂νv. If γ = 1 near ∂Ω then Λq = Λγ. Equivalent inverse problem:

• Does Λq determine q?
• How can q be computed from Λq?

Direct reconstruction in 3D 6/26

Short and incomplete history

1980 Calder´

• n: Problem posed, uniqueness for linearized

problem, and linear, approximate reconstruction algorithm 1987 Sylvester and Uhlmann: Uniqueness for smooth

• conductivities. Implicit reconstruction algorithm

1987-88 Novikov, Nachman-Sylvester-Uhlmann, Nachman: Uniqueness for conductivities with 2 derivatives and explicit high frequency reconstruction algorithm. Multidimensional D-bar equation. 2003 Brown-Torres, P¨ aiv¨ arinta-Panchenko-Uhlmann: Uniqueness for conductivities with 3/2 derivatives. 2006 Cornean-Knudsen-Siltanen: Low frequency reconstruction algorithm 2009 Bikowski-Knudsen-Mueller: Numerical implementation

• f reconstruction algorithms

Direct reconstruction in 3D 7/26

• 2. Calderon’s approximate reconstruction algorithm

Direct reconstruction in 3D 8/26

Calder´

• n’s reconstruction method

Integration by parts (Λγ − Λ1)f, h =

• Ω

(γ − 1)∇u · ∇vdx, where ∇ · γ∇u = 0, u|∂Ω = f ∆v = 0 v|∂Ω = h. Idea: take h, f as restrictions of harmonic functions f = eix·ζ h = e−ix·(ξ+ζ) with ξ ∈ R3 and ζ ∈ C3 s.t. (ξ + ζ)2 = ζ2 = 0. The near-field scattering transform texp(ξ, ζ) =

• (Λγ − Λ1)eix·ζ, e−ix·(ξ+ζ)

=

• Ω

(γ − 1)∇uexp(x, ζ) · ∇e−ix·(ξ+ζ)dx, with

Direct reconstruction in 3D 9/26

Writing uexp = eix·ζ + δu yields texp(ξ, ζ) =

• Ω

(γ − 1)∇uexp(x, ζ) · ∇e−ix·(ζ+ξ)dx = −|ξ|2 2

• Ω

(γ − 1)e−ix·ξdx + R(ξ, ζ). Estimate: |R(ξ, ζ)| ≤ Cγ − 12

L∞(Ω)(1 + |ζ|)2e2R|ζ|.

Ω ⊂ BR Apply low frequency filter χK and invert Fourier transform to get Calder´

• n approximation formula

γapp(x) = 1 − 1 2(2π)n texp(ξ, ζ) |ξ|2 eix·ξχK (ξ)dξ.

Direct reconstruction in 3D 10/26

Writing uexp = eix·ζ + δu yields texp(ξ, ζ) =

• Ω

(γ − 1)∇uexp(x, ζ) · ∇e−ix·(ζ+ξ)dx = −|ξ|2 2

• Ω

(γ − 1)e−ix·ξdx + R(ξ, ζ). Estimate: |R(ξ, ζ)| ≤ Cγ − 12

L∞(Ω)(1 + |ζ|)2e2R|ζ|.

Ω ⊂ BR Apply low frequency filter χK and invert Fourier transform to get Calder´

• n approximation formula

γapp(x) = 1 − 1 2(2π)n texp(ξ, ζ) |ξ|2 eix·ξχK (ξ)dξ. Remark: lim

|ζ|→0 texp(ξ, ζ) = −|ξ|2

γ1/2 − 1(ξ).

Direct reconstruction in 3D 10/26

Calder´

• n reconstruction of potential

With ξ ∈ R3 and ζ ∈ C3 with (ξ + ζ)2 = ζ2 = 0, compute from Λq texp(ξ, ζ) =

• (Λq − Λ0)eix·ζ, e−ix·(ζ+ξ)

=

• Ω

qv exp(x, ζ)e−ix·(ζ+ξ)dx = ˆ q(ξ) + R(ξ, ζ), where (∆ + q)v exp = 0 in Ω and v exp|∂Ω = eix·ξ. Estimate of R as before... Inversion formula q ≈ qapp(x) =

• R3 texp(ξ, ζ)eix·ξχR(ξ)dξ.

Good approximation when q is small.

Direct reconstruction in 3D 11/26

• 3. Exact high complex frequency reconstruction

algorithm

Direct reconstruction in 3D 12/26

Complex geometrical optics

Complex geometrical optics: Let ζ ∈ C3 such that ζ · ζ = 0. For sufficiently large ζ there is a unique solution to the problem (∆ + q)ψ(x, ζ) = 0 in R3, ψ(x, ζ) ∼ eix·ζ for large |x| or |ζ|. ψ can be found by solving the Lippmann-Schwinger-Faddeev (LSF) equation ψ(x, ζ) = eix·ζ +

• Ω

Gζ(x − y)q(y)ψ(y, ζ)dx, ∆Gζ = δ, Gζ ∼ eix·ζ. Moreover, ψ|∂Ω satisfies the solvable Fredholm equation ψ(x, ζ) +

• Ω

Gζ(x − y)(Λq − Λ0)ψ(y, ζ)dσ(y) = eix·ζ, x ∈ ∂Ω.

Direct reconstruction in 3D 13/26

Exceptional points

Exceptional points: ζ ∈ C3 for which there is no unique complex geometrical optics.

• In 2D: q = −∆γ1/2/γ1/2 if and only if there are no exceptional
• points. In 3D?
• Theorem (Cornean - K - Siltanen 2003): If q is small there are no

exceptional points, and γ can be reconstructed by a low complex frequency method.

Direct reconstruction in 3D 14/26

The scattering transform

The key intermediate object, the non-physical scattering transform, t(ξ, ζ) =

• Ω

e−ix·(ξ+ζ)q(x)ψ(x, ζ)dx =

• ∂Ω

e−ix·(ξ+ζ)(Λq − Λ0)ψ(x, ζ)|∂Ωdσ(x), (ξ + ζ)2 = 0. Facts

• 1. ψ|∂Ω can be computed from boundary measurements by solving

ψ + Sζ(Λq − Λ0)ψ = eix·ζ, x ∈ ∂Ω

• 2. t can be computed from boundary data
• 3. q can be computed by using the estimate

|ˆ q(ξ) − t(ξ, ζ)| = O(1/|ζ|) Nonlinear, direct reconstruction algorithm: Λq → t(ξ, ζ) → q(x)(→ γ(x))

Direct reconstruction in 3D 15/26

Connection to Calder´

• n reconstruction

Near-field scattering transform: texp(ξ, ζ) =

• (Λq − Λ0)eix·ζ, e−ix·(ζ+ξ)

=

• Ω

e−ix·(ξ+ζ)q(x)v exp(x, ζ)dx, with (∆ + q)v exp = 0 in Ω and v exp|∂Ω = eix·ζ. Scattering transform: t(ξ, ζ) =

• (Λq − Λ0)ψ, e−ix·(ζ+ξ)

=

• Ω

e−ix·(ξ+ζ)q(x)ψ(x, ζ)dx, where (∆ + q)ψ = 0 in Rn and ψ ∼ eix·ζ near ∞. In 2D (Siltanen-Isaacson-Mueller, 2001) t was replaced by texp.

Direct reconstruction in 3D 16/26

• 4. Implementation details and numerical results

Direct reconstruction in 3D 17/26

Implementation details t

exp

• 1. Solve numerically (∆ + q)v exp = 0 with v exp|∂Ω = eix·ζ
• 2. Integrate numerically
• Ω

qv expeix·(ξ+ζ)dx. Spherically symmetric q(x) = q(|x|) in Ω = B(0, 1) : eix·ζ =

• k,l

al,kY k

l ,

e−ix·(ξ+ζ) =

• k,l

bl,kY k

l

Implies seperation of variables v exp(x) =

• k,l

al,kRl(r)Y k

l (x/r),

r = |x|. Finally texp(ξ, ζ) =

• k,l

al,kbl,k 1 Rlq(r)r2dr.

Direct reconstruction in 3D 18/26

Implementation details t

Computation of Green’s function Gζ(x) = eix·ζgζ(x) : ge1+ie2(x) = e−r+x2−ix1 4πr − 1 4π 1

s

e−ru+x2−ix1 √ 1 − u2 J1(r

• 1 − u2)du,

|x| < 2R from [Newton, 1989] + symmetry.

Direct reconstruction in 3D 19/26

Implementation details t

Computation of Green’s function Gζ(x) = eix·ζgζ(x) : ge1+ie2(x) = e−r+x2−ix1 4πr − 1 4π 1

s

e−ru+x2−ix1 √ 1 − u2 J1(r

• 1 − u2)du,

|x| < 2R from [Newton, 1989] + symmetry. Computation of ψ : technique of Vainikko for solving Lippman-Schwinger eq. µ(x, ζ) = ψ(x, ζ)e−ix·ζ gζ(x) = Gζ(x)e−ix·ζ. Then µ(x, ζ) = 1 +

• Ω

gζ(x − y)q(y)µ(y, ζ)dy.

Direct reconstruction in 3D 19/26

Implementation details t

µ(x, ζ) −

• Ω

gζ(x − y)q(y)µ(y, ζ)dy = 1. Note

• RHS is periodic
• Integral is on bounded domain (compact support of q)

Periodic equation for µp : µp(x, ζ) −

• R3 g

p

ζ(x − y)qp(y)µp(y, ζ)dy = 1.

• Periodic equation is uniquely solvable and on Ω µp(x, ζ) = µ(x, ζ)
• Solved using FFT and an GMRES

Direct reconstruction in 3D 20/26

Implementation details t

µ(x, ζ) −

• Ω

gζ(x − y)q(y)µ(y, ζ)dy = 1. Note

• RHS is periodic
• Integral is on bounded domain (compact support of q)

Periodic equation for µp : µp(x, ζ) −

• R3 g

p

ζ(x − y)qp(y)µp(y, ζ)dy = 1.

• Periodic equation is uniquely solvable and on Ω µp(x, ζ) = µ(x, ζ)
• Solved using FFT and an GMRES

Numerical integration t(ξ, ζ) =

• Ω

e−ix·(ξ+ζ)q(x)ψ(x, ζ)dx.

Direct reconstruction in 3D 20/26

Examples

Take Ω = B(0, 1) and define for α ∈ R+ and 0 < d < 1 γ(x) = ⎧ ⎪ ⎨ ⎪ ⎩

• 1 + αe

−

|x|2 (|x|2−d2)2

2 , |x| ≤ d 0, d < |x| ≤ 1. Convenient with radially symmetric conductivity, since ΛγY k

l = λl,kY k l .

Explicit calculation of eigenvalues Moreover, ˆ q, t and texp have symmetries.

Direct reconstruction in 3D 21/26

Convergence of scattering transform

Particular example with α = 2, d = .9 :

2 4 6 8 10 12 14 16 18 20 −8 −6 −4 −2 2 4 6 8 10 12

|ξ|

Scattering transform t(ξ,ζ) for different values of ζ. Compared to FT q(ξ)

|ζ|=2 |ζ|=4 |ζ| = 8 |ζ|=16 FT(q) Error: |ζ|=2 : 54% |ζ|=4 : 22% |ζ|=8 : 8% |ζ| = 16: 4%

Direct reconstruction in 3D 22/26

Numerical results 1: low contrast

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −2 −1.5 −1 −0.5 0.5 1 1.5 2 2.5

x

Conductivity and potential for α = .1 and d = .9

2 4 6 8 10 12 14 16 18 20 −0.6 −0.5 −0.4 −0.3 −0.2 −0.1 0.1 0.2 0.3

Plot of FT(q) and T(ξ,ζ) and Texp( ξ,ζ)

|ξ| texp FT(q) t, |ζ|=1

Direct reconstruction in 3D 23/26

Numerical results 2: high contrast

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −20 −10 10 20 30 40 50

|x|

Conductivity and potential for α = 2 and d = .9

q(x) γ(x)

2 4 6 8 10 12 14 16 18 20 −10 −5 5 10 15

|ξ|

Plot of FT(q) and t(ξ,ζ) and texp(ξ,ζ)

FT(q) |ζ|=10 texp

Direct reconstruction in 3D 24/26

Calderon reconstruction of conductivity

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 |x| Calderon reconstruction (blue) campared to exact conductivty (dashed) and low frequency approximation of exact conductivty (red). Spectral cutoff at |ξ| = K = 20

Direct reconstruction in 3D 25/26

Conclusion

• Implementation of linear approximate reconstruction algorithm

and non-linear exact method for computing potential/conductivity.

• Linear works well for small potentials
• texp contains medium frequency information of potential
• Suffices to compute t(ξ, ζ) for medium sized |ξ| and |ζ|
• Close connection between Calder ´
• n’s reconstruction and direct,

non-linear method

Direct reconstruction in 3D 26/26

Conclusion

• Implementation of linear approximate reconstruction algorithm

and non-linear exact method for computing potential/conductivity.

• Linear works well for small potentials
• texp contains medium frequency information of potential
• Suffices to compute t(ξ, ζ) for medium sized |ξ| and |ζ|
• Close connection between Calder ´
• n’s reconstruction and direct,

non-linear method

Outlook

• Study numerically 3D examples
• Solution of boundary integral equation for solving inverse

problem

• Understanding spectral cut-off as regularization strategy
• Investigate exceptional points in 3D
• Real data
• The reconstruction problem with data measured only on part of

the boundary

Direct reconstruction in 3D 26/26

Conclusion

• Implementation of linear approximate reconstruction algorithm

and non-linear exact method for computing potential/conductivity.

• Linear works well for small potentials
• texp contains medium frequency information of potential
• Suffices to compute t(ξ, ζ) for medium sized |ξ| and |ζ|
• Close connection between Calder ´
• n’s reconstruction and direct,

non-linear method

Outlook

• Study numerically 3D examples
• Solution of boundary integral equation for solving inverse

problem

• Understanding spectral cut-off as regularization strategy
• Investigate exceptional points in 3D
• Real data
• The reconstruction problem with data measured only on part of

the boundary Thank you

Direct reconstruction in 3D 26/26