Applications of microlocal analysis to inverse problems

TL;DR

This paper discusses how microlocal analysis techniques are applied to solve inverse problems, providing insights and methods for improving problem-solving strategies in this area.

Contribution

It introduces new applications of microlocal analysis specifically tailored for inverse problems, expanding the toolkit for researchers.

Findings

Enhanced understanding of wavefront set analysis in inverse problems

Development of new microlocal techniques for imaging reconstruction

Improved stability estimates for inverse problem solutions

Abstract

These are lecture notes for a minicourse on applications of microlocal analysis in inverse problems, given in Helsinki and Shanghai in June 2019.

Equations

$${\partial }_j=\frac{\partial }{\partial x_j},D_j=\frac{1}{i}{\partial }_j,D=\frac{1}{i}\nabla ,D^{\alpha }=D_1^{{\alpha }_1}\cdots D_n^{{\alpha }_n}.$$

$$\hat{u}(\xi ):={\int }_{{\mathbb{R}}^n}{e}^{-ix\cdot \xi }u(x)dx,\xi \in {\mathbb{R}}^n.$$

$$(D_{j}u)\hat{}(\xi )={\xi }_j\hat{u}(\xi ).$$

$${\mathcal{F}}^{-1}v(x):=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }v(\xi )d\xi$$

$$-\Delta u=f\text{ in Rn}.$$

$$-\Delta u=f$$

$$⟺\hat{u}(\xi )=\frac{1}{|\xi {|}^2}\hat{f}(\xi )$$

$$⟺u(x)={\mathcal{F}}^{-1}\left\{\frac{1}{|\xi {|}^2}\hat{f}(\xi )\right\}=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }\frac{1}{|\xi {|}^2}\hat{f}(\xi )d\xi .$$

$$a(D)u=f\text{ in Rn},a(D)=\sum _{|\alpha |\le m}{a}_{\alpha }{D}^{\alpha },$$

$$a(D)u(x)={\mathcal{F}}^{-1}\left\{a(\xi )\hat{u}(\xi )\right\}=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }a(\xi )\hat{f}(\xi )d\xi .$$

$$u(x)={\mathcal{F}}^{-1}\left\{\frac{1}{a(\xi )}\hat{u}(\xi )\right\}=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }\frac{1}{a(\xi )}\hat{f}(\xi )d\xi .$$

$$Au=f\text{ in Rn},A=a(x,D)=\sum _{|\alpha |\le m}{a}_{\alpha }(x)D^{\alpha },$$

$$Au(x)$$

$$=\sum _{|\alpha |\le m}{a}_{\alpha }(x)D^{\alpha }\left[(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }\hat{u}(\xi )d\xi \right]$$

$$=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }\left[\sum _{|\alpha |\le m}{a}_{\alpha }(x){\xi }^{\alpha }\right]\hat{u}(\xi )d\xi$$

$$=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }a(x,\xi )\hat{u}(\xi )d\xi$$

$$a(x,\xi ):=\sum _{|\alpha |\le m}{a}_{\alpha }(x){\xi }^{\alpha }$$

$$u(x)=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }\frac{1}{a(x,\xi )}\hat{f}(\xi )d\xi .$$

$${\sigma }_{\mathrm{pr}}(A):=\sum _{|\alpha |=m}{a}_{\alpha }(x){\xi }^{\alpha }.$$

$$u_1(x):=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }b(x,\xi )\hat{f}(\xi )d\xi$$

$$b(x,\xi ):=\frac{1-\psi (\xi )}{a(x,\xi )},$$

$$A{u}_1=f+f_1$$

$$A{u}_{\mathrm{app}}=f+r,r\in C^{\infty }({\mathbb{R}}^n).$$

$$Au(x)=(2\pi {)}^{-n}{\int }_{{\mathbb{R}}^n}{e}^{ix\cdot \xi }a(x,\xi )\hat{u}(\xi )d\xi$$

$$|{\partial }_x^{\alpha }{\partial }_{\xi }^{\beta }a(x,\xi )|\le C_{\alpha ,\beta }(1+|\xi |{)}^{m-|\beta |},\xi \in {\mathbb{R}}^n.$$

$$G:f\mapsto {\mathcal{F}}^{-1}\left\{\frac{1}{|\xi {|}^2}\hat{f}(\xi )\right\}.$$

$$G=G_1+G_2,G_1:={\mathcal{F}}^{-1}\left\{\frac{1-\psi (\xi )}{|\xi {|}^2}\hat{f}(\xi )\right\},G_2:={\mathcal{F}}^{-1}\left\{\frac{\psi (\xi )}{|\xi {|}^2}\hat{f}(\xi )\right\},$$

$$b(x,\xi )\sim \sum _{j=0}^{\infty }{b}_{m-j}(x,\xi )$$

$${\sigma }_{\mathrm{pr}}(A)(x,\xi )=\underset{\lambda \to \infty }{\lim }{\lambda }^{-m}{e}^{-i\lambda x\cdot \xi }A(e^{i\lambda x\cdot \xi });$$

$$\mathrm{sing}\mathrm{supp}(u)$$

$${\mathrm{sing}\mathrm{supp}}^{\alpha }(u)$$

$$u=\sum _{j=1}^N{c}_j{\chi }_{D_j}$$

$${\mathrm{sing}\mathrm{supp}}^{\alpha }(u)=\varnothing \text{ for α<1/2}$$

Full text

Preface

Microlocal analysis originated in the 1950s, and by now it is a substantial mathematical theory with many different facets and applications. One might view microlocal analysis as

• •

a kind of ”variable coefficient Fourier analysis” for solving variable coefficient PDEs; or

• •

as a theory of pseudodifferential operators ($\Psi$DOs) and Fourier integral operators (FIOs); or

• •

as a phase space (or time-frequency) approach to studying functions, operators and their singularities (wave front sets).

$\Psi$DOs were introduced by Kohn and Nirenberg [KN65], and FIOs and wave front sets were studied systematically by Hörmander [Hö71]. Much of the theory up to the early 1980s is summarized in the four volume treatise of Hörmander [Hö85]. There are remarkable applications of microlocal analysis and related ideas in many fields of mathematics. Classical examples include spectral theory and the Atiyah-Singer index theorem, and more recent examples include scattering theory, behavior of chaotic systems, inverse problems, and general relativity.

In this minicourse we will try to describe some classical applications of microlocal analysis to inverse problems, together with a very rough non-technical overview of relevant parts of microlocal analysis. In a nutshell, here are a few typical applications:

Computed tomography / X-ray transform: the X-ray transform is an FIO, and under certain conditions its normal operator is an elliptic $\Psi$DO. Microlocal analysis can be used to predict which sharp features (singularities) of the image can be reconstructed in a stable way from limited data measurements. Microlocal analysis is also a powerful tool in the study of geodesic X-ray transforms related to seismic imaging applications. 2. 2.

Calderón problem / Electrical Impedance Tomography: the boundary measurement map (Dirichlet-to-Neumann map) is a $\Psi$DO, and the boundary values of the conductivity as well as its derivatives can be computed from the symbol of this $\Psi$DO. 3. 3.

Gel’fand problem / seismic imaging: the boundary measurement operator (hyperbolic Dirichlet-to-Neumann map) is an FIO, and the scattering relation of the sound speed as well as certain X-ray transforms of the coefficients can be computed from the canonical relation and the symbol of this FIO.

These notes are organized as follows. In Section 1, we will motivate the theory of $\Psi$DOs and discuss some of its properties without giving proofs. Section 2 will continue with a brief introduction to wave front sets and FIOs (again with no proofs). The rest of the notes is concerned with applications to inverse problems. Section 3 considers the Radon transform in $\mathbb{R}^{2}$ and its normal operator, and describes what kind of information about the singularities of $f$ can be stably recovered from the Radon transform. Sections 4 and 5 discuss the Gel’fand and Calderón problems, and prove results related to recovering X-ray transforms or boundary determination. The treatment is motivated by $\Psi$DO and FIO theory, but we give direct and (in principle) elementary proofs based on a quasimode constructions. The results discussed in these notes are classical. For more recent results, we refer to the surveys [IM19, KQ15, La18, Uh14].

Notation

We will use multi-index notation. Let $\mathbb{N}_{0}=\{0,1,2,\ldots\}$ be the set natural numbers. Then $\mathbb{N}_{0}^{n}$ consists of all $n$-tuples $\alpha=(\alpha_{1},\ldots,\alpha_{n})$ where the $\alpha_{j}$ are nonnegative integers. Such an $n$-tuple $\alpha$ is called a multi-index. We write $|\alpha|=\alpha_{1}+\ldots+\alpha_{n}$ and $\xi^{\alpha}=\xi_{1}^{\alpha_{1}}\cdots\xi_{n}^{\alpha_{n}}$ for $\xi\in\mathbb{R}^{n}$. For partial derivatives, we will write

[TABLE]

If $\Omega\subset\mathbb{R}^{n}$ is a bounded domain with $C^{\infty}$ boundary, we denote by $C^{\infty}(\overline{\Omega})$ the set of infinitely differentiable functions in $\Omega$ whose all derivatives extend continuously to $\overline{\Omega}$. The space $C^{\infty}_{c}(\Omega)$ consist of $C^{\infty}$ functions having compact support in $\Omega$. The standard $L^{2}$ based Sobolev spaces are denoted by $H^{s}(\mathbb{R}^{n})$ with norm $\lVert f\rVert_{H^{s}(\mathbb{R}^{n})}=\lVert(1+\lvert\xi\rvert^{2})^{s/2}\hat{f}\rVert_{L^{2}(\mathbb{R}^{n})}$, with $\hat{f}$ denoting the Fourier transform. We also write $\lVert f\rVert_{W^{k,\infty}}=\sum_{\lvert\alpha\rvert\leq k}\lVert D^{\alpha}f\rVert_{L^{\infty}}$. The notation $A\lesssim B$ means that $A\leq CB$ for some uniform (with respect to the relevant parameters) constant $C$. In general, all coefficients, boundaries etc are assumed to be $C^{\infty}$ for ease of presentation.

1. Pseudodifferential operators

In this minicourse we will try to give a very brief idea of the different points of view to microlocal analysis mentioned in the introduction (and repeated below), as

1. (1)

a kind of ”variable coefficient Fourier analysis” for solving variable coefficient PDEs; or 2. (2)

a theory of $\Psi$DOs and FIOs; or 3. (3)

a phase space (or time-frequency) approach to studying functions, operators and their singularities (wave front sets).

In this section we will discuss (1) and (2) in the context of $\Psi$DOs (we will continue with (2) and (3) in the context of FIOs in Section 2). The treatment is mostly formal and we will give no proofs whatsoever. A complete reference for the results in this section is [Hö85, Section 18.1].

1.1. Constant coefficient PDEs

We recall the following facts about the Fourier transform (valid for sufficiently nice functions):

If $u$ is a function in $\mathbb{R}^{n}$, its Fourier transform $\hat{u}=\mathscr{F}u$ is the function

[TABLE] 2. 2.

The Fourier transform converts derivatives to polynomials (this is why it is useful for solving PDEs):

[TABLE] 3. 3.

A function $u$ can be recovered from $\hat{u}$ by the Fourier inversion formula $u=\mathscr{F}^{-1}\{\hat{u}\}$, where

[TABLE]

is the inverse Fourier transform.

As a motivating example, let us solve formally (i.e. without worrying about how to precisely justify each step) the equation

[TABLE]

This is a constant coefficient PDE, and such equations can be studied with the help of the Fourier transform. We formally compute

[TABLE]

The same formal argument applies to a general constant coefficient PDE

[TABLE]

where $a_{\alpha}\in\mathbb{C}$. Then $(a(D)u)\,\hat{\rule{0.0pt}{6.0pt}}\,(\xi)=a(\xi)\hat{u}(\xi)$ where $a(\xi)=\sum_{\lvert\alpha\rvert\leq m}a_{\alpha}\xi^{\alpha}$ is the symbol of $a(D)$. Moreover, one has

[TABLE]

The argument leading to (1.1) gives a formal solution of $a(D)u=f$:

[TABLE]

Thus formally $a(D)u=f$ can be solved by dividing by the symbol $a(\xi)$ on the Fourier side. Of course, to make this precise one would need to show that the division by $a(\xi)$ (which may have zeros) is somehow justified.

1.2. Variable coefficient PDEs

We now try to use a similar idea to solve the variable coefficient PDE

[TABLE]

where $a_{\alpha}(x)\in C^{\infty}(\mathbb{R}^{n})$ and $D^{\beta}a_{\alpha}\in L^{\infty}(\mathbb{R}^{n})$ for all multi-indices $\alpha,\beta$. Since the coefficients $a_{\alpha}$ depend on $x$, Fourier transforming the equation $Au=f$ is not immediately helpful. However, we can compute an analogue of (1.2):

[TABLE]

where

[TABLE]

is the (full) symbol of $A=a(x,D)$.

Now, we could try to obtain a solution to $a(x,D)u=f$ in $\mathbb{R}^{n}$ by dividing by the symbol $a(x,\xi)$ as in (1.3):

[TABLE]

Again, this is only formal since the division by $a(x,\xi)$ needs to be justified. However, this can be done in a certain sense if $A$ is elliptic:

Definition**.**

The principal symbol (i.e. the part containing the highest order derivatives) of the differential operator $A=a(x,D)$ is

[TABLE]

We say that $A$ is elliptic if its principal symbol is nonvanishing for $\xi\neq 0$.

A basic result of microlocal analysis states that the function

[TABLE]

with

[TABLE]

where $\psi\in C^{\infty}_{c}(\mathbb{R}^{n})$ is a cutoff with $\psi(\xi)=1$ in a sufficiently large neighborhood of $\xi=0$ (so that $a(x,\xi)$ does not vanish outside this neighborhood), is an approximate solution of $Au=f$ in the sense that

[TABLE]

where $f_{1}$ is one derivative smoother than $f$. Moreover, it is possible to construct an approximate solution $u_{\mathrm{app}}$ so that

[TABLE]

1.3. Pseudodifferential operators

In analogy with the formula (1.4), a pseudodifferential operator ($\Psi$DO) is an operator $A$ of the form

[TABLE]

where $a(x,\xi)$ is a symbol with certain properties. The most standard symbol class $S^{m}=S^{m}_{1,0}(\mathbb{R}^{n})$ is defined as follows:

Definition**.**

The symbol class $S^{m}$ consists of functions $a\in C^{\infty}(\mathbb{R}^{n}\times\mathbb{R}^{n})$ such that for any $\alpha,\beta\in\mathbb{N}_{0}^{n}$ there is $C_{\alpha,\beta}>0$ with

[TABLE]

If $a\in S^{m}$, the corresponding $\Psi$DO $A=\mathrm{Op}(a)$ is defined by (1.7). We denote by $\Psi^{m}$ the set of $\Psi$DOs corresponding to $S^{m}$.

Note that symbols in $S^{m}$ behave roughly like polynomials of order $m$ in the $\xi$-variable. In particular, the symbols $a(x,\xi)$ in (1.5) belong to $S^{m}$ and the corresponding differential operators $a(x,D)$ belong to $\Psi^{m}$. Moreover, if $a(x,D)$ is elliptic, then the symbol $b(x,\xi)=\frac{1-\psi(\xi)}{a(x,\xi)}$ as in (1.6) belongs to $S^{-m}$. Thus the class of $\Psi$DOs is large enough to include differential operators as well as approximate inverses of elliptic operators. Also normal operators of the X-ray transform or Radon transform in $\mathbb{R}^{n}$ are $\Psi$DOs.

Remark 1.1** (Homogeneous symbols).**

We saw in Section 1.1 that the elliptic operator $-\Delta$ has the inverse

[TABLE]

The symbol $\frac{1}{\lvert\xi\rvert^{2}}$ is not in $S^{-2}$, since it is not smooth near [math]. However, one often thinks of $G$ as a $\Psi$DO by writing

[TABLE]

where $\psi\in C^{\infty}_{c}(\mathbb{R}^{n})$ satisfies $\psi=1$ near [math]. Now $G_{1}$ is a $\Psi$DO in $\Psi^{-2}$, since $\frac{1-\psi(\xi)}{\lvert\xi\rvert^{2}}\in S^{-2}$, and $G_{2}$ is smoothing in the sense that it maps any $L^{1}$ function into a $C^{\infty}$ function (at least if $n\geq 3$).

In general, in $\Psi$DO theory smoothing operators are considered to be negligible (since at least they do not introduce new singularities), and many computations in $\Psi$DO calculus are made only modulo smoothing error terms. In this sense one often views $G$ as a $\Psi$DO by identifying it with $G_{1}$. The same kind of identification is done for operators whose symbol $a(x,\xi)$ is homogenous of some order $m$ in $\xi$. More generally one can consider polyhomogeneous symbols $b\in S^{m}$ having the form

[TABLE]

where each $b_{m-j}$ is homogeneous of order $m-j$ in $\xi$, and $\sim$ is a certain asympotic summation. Corresponding $\Psi$DOs are called classical $\Psi$DOs.

It is very important that one can compute with $\Psi$DOs in much the same way as with differential operators. One often says that $\Psi$DOs have a calculus. The following theorem lists typical rules of computation (it is instructive to think first why such rules are valid for differential operators):

Theorem 1.2** ($\Psi$DO calculus).**

**

1. (a)

(Principal symbol) There is a one-to-one correspondence between operators in $\Psi^{m}$ and (full) symbols in $S^{m}$, and each operator $A\in\Psi^{m}$ has a well defined principal symbol $\sigma_{\mathrm{pr}}(A)$. The principal symbol may be computed by testing $A$ against highly oscillatory functions111This is valid if $A$ is a classical $\Psi$DO.:

[TABLE] 2. (b)

(Composition) If $A\in\Psi^{m}$ and $B\in\Psi^{m^{\prime}}$, then $AB\in\Psi^{m+m^{\prime}}$ and $\sigma_{\mathrm{pr}}(AB)=\sigma_{\mathrm{pr}}(A)\sigma_{\mathrm{pr}}(B)$; 3. (c)

(Sobolev mapping properties) Each $A\in\Psi^{m}$ is a bounded operator $H^{s}(\mathbb{R}^{n})\to H^{s-m}(\mathbb{R}^{n})$ for any $s\in\mathbb{R}$; 4. (d)

(Elliptic operators have approximate inverses) If $A\in\Psi^{m}$ is elliptic, there is $B\in\Psi^{-m}$ so that $AB=\mathrm{Id}+K$ and $BA=\mathrm{Id}+L$ where $K,L\in\Psi^{-\infty}$, i.e. $K,L$ are smoothing (they map any $H^{-s}$ function to $H^{t}$ for any $t$, hence also to $C^{\infty}$ by Sobolev embedding).

The above properties are valid in the standard $\Psi$DO calculus in $\mathbb{R}^{n}$. However, motivated by different applications, $\Psi$DOs have been considered in various other settings. Each of these settings comes with an associated calculus whose rules of computation are similar but adapted to the situation at hand. Examples of different settings for $\Psi$DO calculus include

1. (1)

open sets in $\mathbb{R}^{n}$ (local setting); 2. (2)

compact manifolds without boundary, possibly acting on sections of vector bundles; 3. (3)

compact manifolds with boundary (transmission condition / Boutet de Monvel calculus); 4. (4)

non-compact manifolds (e.g. Melrose scattering calculus); and 5. (5)

operators with a small or large parameter (semiclassical calculus).

2. Wave front sets and Fourier integral operators

For a reference to wave front sets, see [Hö85, Chapter 8]. Sobolev wave front sets are considered in [Hö85, Section 18.1]. FIOs are discussed in [Hö85, Chapter 25].

2.1. The role of singularities

We first discuss the singular support of $u$, which consists of those points $x_{0}$ such that $u$ is not a smooth function in any neighborhood of $x_{0}$. We also consider the Sobolev singular support, which also measures the ”strength” of the singularity (in the $L^{2}$ Sobolev scale).

Definition** (Singular support).**

We say that a function or distribution $u$ is $C^{\infty}$ (resp. $H^{\alpha}$) near $x_{0}$ if there is $\varphi\in C^{\infty}_{c}(\mathbb{R}^{n})$ with $\varphi=1$ near $x_{0}$ such that $\varphi u$ is in $C^{\infty}(\mathbb{R}^{n})$ (resp. in $H^{\alpha}(\mathbb{R}^{n})$). We define

[TABLE]

Example 2.1**.**

Let $D_{1},\ldots,D_{N}$ be bounded domains with $C^{\infty}$ boundary in $\mathbb{R}^{n}$ so that $\overline{D}_{j}\cap\overline{D}_{k}=\emptyset$ for $j\neq k$, and define

[TABLE]

where $c_{j}\neq 0$ are constants, and $\chi_{D_{j}}$ is the characteristic function of $D_{j}$. Then

[TABLE]

since $u\in H^{\alpha}$ for $\alpha<1/2$, but

[TABLE]

since $u$ is not $H^{1/2}$ near any boundary point. Thus in this case the singularities of $u$ are exactly at the points where $u$ has a jump discontinuity, and their strength is precisely $H^{1/2}$. Knowing the singularities of $u$ can already be useful in applications. For instance, if $u$ represents some internal medium properties in medical imaging, the singularities of $u$ could determine the location of interfaces between different tissues. On the other hand, if $u$ represents an image, then the singularities in some sense determine the ”sharp features” of the image.

Next we discuss the wave front set which is a more refined notion of a singularity. For example, if $f=\chi_{D}$ is the characteristic function of a bounded strictly convex $C^{\infty}$ domain $D$ and if $x_{0}\in\partial D$, one could think that $f$ is in some sense smooth in tangential directions at $x_{0}$ (since $f$ restricted to a tangent hyperplane is identically zero, except possibly at $x_{0}$), but that $f$ is not smooth in normal directions at $x_{0}$ since in these directions there is a jump. The wave front set is a subset of $T^{*}\mathbb{R}^{n}\setminus 0$, the cotangent space with the zero section removed:

[TABLE]

Definition** (Wave front set).**

Let $u$ be a distribution in $\mathbb{R}^{n}$. We say that $u$ is (microlocally) $C^{\infty}$ (resp. $H^{\alpha}$) near $(x_{0},\xi_{0})$ if there exist $\varphi\in C^{\infty}_{c}(\mathbb{R}^{n})$ with $\varphi=1$ near $x_{0}$ and $\psi\in C^{\infty}(\mathbb{R}^{n}\setminus\{0\})$ so that $\psi=1$ near $\xi_{0}$ and $\psi$ is homogeneous of degree [math], such that

[TABLE]

(resp. $\mathscr{F}^{-1}\{\psi(\xi)(\varphi u)\,\hat{\rule{0.0pt}{6.0pt}}\,(\xi)\}\in H^{\alpha}(\mathbb{R}^{n})$). The wave front set $WF(u)$ (resp. $H^{\alpha}$ wave front set $WF^{\alpha}(u)$) consists of those points $(x_{0},\xi_{0})$ where $u$ is not microlocally $C^{\infty}$ (resp. $H^{\alpha}$).

Example 2.2**.**

The wave front set of the function $u$ in Example 2.1 is

[TABLE]

where $N^{*}(D_{j})$ is the conormal bundle of $D_{j}$,

[TABLE]

The wave front set describes singularities more precisely than the singular support, since one always has

[TABLE]

where $\pi:(x,\xi)\mapsto x$ is the projection to $x$-space.

It is an important fact that applying a $\Psi$DO to a function or distribution never creates new singularities:

Theorem 2.3** (Pseudolocal/microlocal property of $\Psi$DOs).**

Any $A\in\Psi^{m}$ has the pseudolocal property

[TABLE]

and the microlocal property

[TABLE]

Elliptic operators are those that completely preserve singularities:

Theorem 2.4**.**

(Elliptic regularity) Let $A\in\Psi^{m}$ be elliptic. Then, for any $u$,

[TABLE]

Thus any solution $u$ of $Au=f$ is singular precisely at those points where $f$ is singular. There are corresponding statements for Sobolev singularities.

Proof.

First note that by Theorem 2.3,

[TABLE]

Conversely, since $A\in\Psi^{m}$ is elliptic, by Theorem 1.2(d) there is $B\in\Psi^{-m}$ so that

[TABLE]

Thus for any $u$ one has

[TABLE]

Since $L$ is smoothing, $Lu\in C^{\infty}$, which implies that $u=BAu$ modulo $C^{\infty}$. Thus it follows that

[TABLE]

Thus $WF(Au)=WF(u)$. The claim for singular supports follows by (2.1). ∎

2.2. Fourier integral operators

We have seen in Section 1.3 that the class of pseudodifferential operators includes approximate inverses of elliptic operators. In order to handle approximate inverses for hyperbolic and transport equations, it is required to work with a larger class of operators.

Motivation**.**

Consider the initial value problem for the wave equation,

[TABLE]

This is again a constant coefficient PDE, and we will solve this formally by taking the Fourier transform in space,

[TABLE]

After taking Fourier transforms in space, the above equation becomes

[TABLE]

For each fixed $\xi$ this is an ODE in $t$, and the solution is

[TABLE]

Taking inverse Fourier transforms in space, we obtain

[TABLE]

Generalizing (2.2), we can consider operators of the form

[TABLE]

where $a(x,\xi)$ is a symbol (for instance in $S^{m}$), and $\varphi(x,\xi)$ is a real valued phase function. Such operators are examples of Fourier integral operators (more precisely, FIOs whose canonical relation is locally the graph of a canonical transformation, see [Hö85, Section 25.3]). For $\Psi$DOs the phase function is always $\varphi(x,\xi)=x\cdot\xi$, but for FIOs the phase function can be quite general, though it is usually required to be homogeneous of degree $1$ in $\xi$, and to satisfy the non-degeneracy condition $\det(\partial_{x_{j}\xi_{k}}\varphi)\neq 0$.

We will not go into precise definitions, but only remark that the class of FIOs includes pseudodifferential operators as well as approximate inverses of hyperbolic and transport operators (or more generally real principal type operators). There is a calculus for FIOs, analogous to the pseudodifferential calculus, under certain conditions in various settings. An important property of FIOs is that they, unlike pseudodifferential operators, can move singularities. This aspect will be discussed next.

2.3. Propagation of singularities

Example 2.5**.**

Let $t>0$ be fixed, and consider the operators from (2.2),

[TABLE]

Then

[TABLE]

Using FIO theory, since the phase functions are $\varphi(x,\xi)=x\cdot\xi\mp t\lvert\xi\rvert$, it follows that

[TABLE]

where $\chi_{\pm t}$ is the canonical transformation (i.e. diffeomorphism of $T^{*}\mathbb{R}^{n}\setminus 0$ that preserves the symplectic structure) given by

[TABLE]

This means that the FIO $A_{\pm}$ takes a singularity $(x,\xi)$ of the initial data $f$ and moves it along the line through $x$ in direction $\pm\xi/\lvert\xi\rvert$ to $(x\pm t\xi/\lvert\xi\rvert,\xi)$. Thus singularities of solutions of the wave equation $(\partial_{t}^{2}-\Delta)u=0$ propagate along straight lines with constant speed one.

Remark 2.6**.**

In general, any FIO has an associated canonical relation that describes what the FIO does to singularities. The canonical relation of the FIO $A$ defined in (2.3) is (see [Hö85, Section 25.3])

[TABLE]

and $A$ moves singularities according to the rule

[TABLE]

where

[TABLE]

Using these formulas, it is easy to check that the canonical relation $C_{\pm}$ of $A_{\pm t}$ in Example 2.5 is the graph of $\chi_{\pm t}$ in the sense that

[TABLE]

and one indeed has $WF(A_{\pm t}u)\subset C_{\pm}(WF(u))=\chi_{\pm t}(WF(u))$.

There is a far reaching extension of Example 2.5, which shows that the singularities of a solution of $Pu=0$ propagate along certain curves in phase space (so called null bicharacteristic curves) as long as $P$ has real valued principal symbol.

Theorem 2.7** (Propagation of singularities).**

Let $P\in\Psi^{m}$ have real principal symbol $p_{m}$ that is homogeneous of degree $m$ in $\xi$. If

[TABLE]

then $WF(u)\setminus WF(f)$ is contained in the characteristic set $p_{m}^{-1}(0)$. Moreover, if $(x_{0},\xi_{0})\in WF(u)\setminus WF(f)$, then the whole null bicharacteristic curve $(x(t),\xi(t))$ through $(x_{0},\xi_{0})$ is in $WF(u)\setminus WF(f)$, where

[TABLE]

Example 2.8**.**

We compute the null bicharacteristic curves for the wave operator $P=\frac{1}{2}(\Delta-\partial_{t}^{2})$. The principal symbol of $P$ is

[TABLE]

The characteristic set is

[TABLE]

which consists of light-like cotangent vectors on $\mathbb{R}^{n+1}_{x,t}$. The equations for the null bicharacteristic curves are

[TABLE]

Thus, if $\lvert\xi_{0}\rvert=1$, then the null bicharacteristic curve through $(x_{0},t_{0},\xi_{0},\pm 1)$ is

[TABLE]

The result of Example 2.5 may thus be interpreted so that singularities of solutions of the wave equation propagate along null bicharacteristic curves for the wave operator.

3. The Radon transform in the plane

In this section we outline some applications of microlocal analysis to the study of the Radon transform in the plane. Similar ideas apply to X-ray and Radon transforms in higher dimensions and Riemannian manifolds as well. The microlocal approach to Radon transforms was introduced by Guillemin [Gu75]. We refer to [Qu06], [KQ15] and references therein for a more detailed treatment of the material in this section.

3.1. Basic properties of the Radon transform

The X-ray transform $If$ of a function $f$ in $\mathbb{R}^{n}$ encodes the integrals of $f$ over all straight lines, whereas the Radon transform $Rf$ encodes the integrals of $f$ over $(n-1)$-dimensional planes. We will focus on the case $n=2$, where the two transforms coincide. There are many ways to parametrize the set of lines in $\mathbb{R}^{2}$. We will parametrize lines by their direction vector $\omega$ and distance $s$ from the origin.

Definition**.**

If $f\in C^{\infty}_{c}(\mathbb{R}^{2})$, the Radon transform of $f$ is the function

[TABLE]

Here $\omega^{\perp}$ is the vector in $S^{1}$ obtained by rotating $\omega$ counterclockwise by $90^{\circ}$.

There is a well-known relation between $Rf$ and the Fourier transform $\hat{f}$. We denote by $(Rf)\,\tilde{\rule{0.0pt}{6.0pt}}\,(\,\cdot\,,\omega)$ the Fourier transform of $Rf$ with respect to $s$.

Theorem 3.1**.**

(Fourier slice theorem)

[TABLE]

Proof.

Parametrizing $\mathbb{R}^{2}$ by $y=s\omega^{\perp}+t\omega$, we have

[TABLE]

This result gives the first proof of injectivity of the Radon transform:

Corollary 3.2**.**

If $f\in C^{\infty}_{c}(\mathbb{R}^{2})$ is such that $Rf\equiv 0$, then $f\equiv 0$.

Proof.

If $Rf\equiv 0$, then $\hat{f}\equiv 0$ by Theorem 3.1 and consequently $f\equiv 0$. ∎

To obtain a different inversion method, and for later purposes, we will consider the adjoint of $R$. The formal adjoint of $R$ is the backprojection operator222The formula for $R^{*}$ is obtained as follows: if $f\in C^{\infty}_{c}(\mathbb{R}^{2})$, $h\in C^{\infty}(\mathbb{R}\times S^{1})$ one has

$\displaystyle(Rf,h)_{L^{2}(\mathbb{R}\times S^{1})}$ $\displaystyle=\int_{-\infty}^{\infty}\int_{S^{1}}Rf(s,\omega)\overline{h(s,\omega)}\,d\omega\,ds$

$\displaystyle=\int_{-\infty}^{\infty}\int_{S^{1}}\int_{-\infty}^{\infty}f(s\omega^{\perp}+t\omega)\overline{h(s,\omega)}\,dt\,d\omega\,ds$

$\displaystyle=\int_{\mathbb{R}^{2}}f(y)\left(\int_{S^{1}}\overline{h(y\cdot\omega^{\perp},\omega)}\,d\omega\right)\,dy.$

[TABLE]

The following result shows that the normal operator $R^{*}R$ is a classical $\Psi$DO of order $-1$ in $\mathbb{R}^{2}$, and also gives an inversion formula.

Theorem 3.3**.**

(Normal operator) One has

[TABLE]

and $f$ can be recovered from $Rf$ by the formula

[TABLE]

Remark 3.4**.**

Above we have written, for $\alpha\in\mathbb{R}$,

[TABLE]

The notation $(-\Delta)^{\alpha/2}=\lvert D\rvert^{\alpha}$ is also used.

Proof.

The proof is based on computing $(Rf,Rg)_{L^{2}(\mathbb{R}\times S^{1})}$ using the Parseval identity, Fourier slice theorem, symmetry and polar coordinates:

[TABLE]

The same argument, based on computing $(\lvert D_{s}\rvert^{1/2}Rf,\lvert D_{s}\rvert^{1/2}Rg)_{L^{2}(\mathbb{R}\times S^{1})}$ instead of $(Rf,Rg)_{L^{2}(\mathbb{R}\times S^{1})}$, leads to the famous filtered backprojection (FBP) inversion formula:

[TABLE]

where $\lvert D_{s}\rvert Rf=\mathscr{F}^{-1}\{\lvert\sigma\rvert(Rf)\,\tilde{\rule{0.0pt}{6.0pt}}\,\}$. This formula is efficient to implement and gives good reconstructions when one has complete X-ray data and relatively small noise, and hence FBP (together with its variants) has been commonly used in X-ray CT scanners.

However, if one is mainly interested in the singularities (i.e. jumps or sharp features) of the image, it is possible to use the even simpler backprojection method: just apply the backprojection operator $R^{*}$ to the data $Rf$. Since $R^{*}R$ is an elliptic $\Psi$DO, Theorem 2.4 guarantees that the singularities are recovered:

[TABLE]

Moreover, since $R^{*}R$ is a $\Psi$DO of order $-1$, hence smoothing of order $1$, one expects that $R^{*}Rf$ gives a slightly blurred version of $f$ where the main singularities should still be visible.

3.2. Visible singularities

There are various imaging situations where complete X-ray data (i.e. the function $Rf(s,\omega)$ for all $s$ and $\omega$) is not available. This is the case for limited angle tomography (e.g. in luggage scanners at airports, or dental applications), region of interest tomography, or exterior data tomography. In such cases explicit inversion formulas such as FBP are usually not available, but microlocal analysis (for related normal operators or FIOs) still provides a powerful paradigm for predicting which singularities can be recovered stably from the measurements.

We will try to explain this paradigm a little bit more, starting with an example:

Example 3.5**.**

Let $f$ be the characteristic function of the unit disc $\mathbb{D}$, i.e. $f(x)=1$ if $\lvert x\rvert\leq 1$ and $f(x)=0$ for $\lvert x\rvert>1$. Then $f$ is singular precisely on the unit circle (in normal directions). We have

[TABLE]

Thus $Rf$ is singular precisely at those points $(s,\omega)$ with $\lvert s\rvert=1$, which correspond to those lines that are tangent to the unit circle.

There is a similar relation between the singularities of $f$ and $Rf$ in general, and this is explained by microlocal analysis:

Theorem 3.6**.**

The operator $R$ is an elliptic FIO of order $-1/2$. There is a precise relationship between the singularities of $f$ and singularities of $Rf$.

We will not spell out the precise relationship here, but only give some consequences. It will be useful to think of the Radon transform as defined on the set of (non-oriented) lines in $\mathbb{R}^{2}$. If $\mathcal{A}$ is an open subset of lines in $\mathbb{R}^{2}$, we consider the Radon transform $Rf|_{\mathcal{A}}$ restricted to lines in $\mathcal{A}$. Recovering $f$ (or some properties of $f$) from $Rf|_{\mathcal{A}}$ is a limited data tomography problem. Examples:

• •

If \mathcal{A}=\{\text{lines not meeting \overline{\mathbb{D}}}\}, then $Rf|_{\mathcal{A}}$ is called exterior data.

• •

If $0<a<\pi/2$ and \mathcal{A}=\{\text{lines whose angle with x$-axis is$<a}\}, then $Rf|_{\mathcal{A}}$ is called limited angle data.

It is known that any $f\in C^{\infty}_{c}(\mathbb{R}^{2}\setminus\overline{D})$ is uniquely determined by exterior data (Helgason support theorem), and any $f\in C^{\infty}_{c}(\mathbb{R}^{2})$ is uniquely determined by limited angle data (Fourier slice and Paley-Wiener theorems). However, both inverse problems are very unstable (inversion is not Lipschitz continuous in any Sobolev norms, but one has conditional logarithmic stability).

Definition**.**

A singularity at $(x_{0},\xi_{0})$ is called visible from $\mathcal{A}$ if the line through $x_{0}$ in direction $\xi_{0}^{\perp}$ is in $\mathcal{A}$.

One has the following dichotomy:

• •

If $(x_{0},\xi_{0})$ is visible from $\mathcal{A}$, then from the singularities of $Rf|_{\mathcal{A}}$ one can determine for any $\alpha$ whether or not $(x_{0},\xi_{0})\in WF^{\alpha}(f)$. If $Rf|_{\mathcal{A}}$ uniquely determines $f$, one expects the reconstruction of visible singularities to be stable.

• •

If $(x_{0},\xi_{0})$ is not visible from $\mathcal{A}$, then this singularity is smoothed out in the measurement $Rf|_{\mathcal{A}}$. Even if $Rf|_{\mathcal{A}}$ would determine $f$ uniquely, the inversion is not Lipschitz stable in any Sobolev norms.

4. Gel’fand problem

Seismic imaging gives rise to various inverse problems related to determining interior properties, e.g. oil deposits or deep structure, of the Earth. Often this is done by using acoustic or elastic waves. We will consider the following problem, also known as the inverse boundary spectral problem (see the monograph [KKL01]):

Gel’fand problem: Is it possible to determine the interior structure of Earth by controlling acoustic waves and measuring vibrations at the surface?

In seismic imaging one often tries to recover an unknown sound speed. However, in this presentation we consider the simpler case where the sound speed is constant (equal to one) and one attempts to recover an unknown potential $q\in C^{\infty}_{c}(\Omega)$ at each point $x\in\Omega$, where $\Omega$ is a ball in $\mathbb{R}^{n}$.

Consider the free wave operator

[TABLE]

We assume that the medium is at rest at time $t=0$ and that we take measurements until time $T>0$. If we prescribe the amplitude of the wave to be $f(x,t)$ on $\partial\Omega\times(0,T)$, this leads to a solution $u$ of the wave equation

[TABLE]

Given any $f\in C^{\infty}_{c}(\partial\Omega\times(0,T))$, there is a unique solution $u\in C^{\infty}(\Omega\times(0,T))$ (see [Ev10, Theorem 7 in §7.2.3]). We assume that we can measure the normal derivative $\partial_{\nu}u|_{\partial\Omega\times(0,T)}$, where $\partial_{\nu}u(x,t)=\nabla_{x}u(x,t)\cdot\nu(x)$ and $\nu$ is the outer unit normal to $\partial\Omega$. Doing such measurements for many different functions $f$, the ideal boundary measurements are encoded by the hyperbolic Dirichlet-to-Neumann map (DN map for short)

[TABLE]

The Gel’fand problem for this model amounts to recovering $q(x)$ from the knowledge of the map $\Lambda_{q}$. We will prove the following result due to [RS88].

Theorem 4.1** (Recovering the X-ray transform).**

Let $T>0$ and assume that $q_{1},q_{2}\in C^{\infty}_{c}(\Omega)$. If

[TABLE]

then $q_{1}$ and $q_{2}$ satisfy

[TABLE]

whenever $\gamma$ is a maximal line segment in $\overline{\Omega}$ with length $<T$.

It is natural that the region where one can recover information depends on $T$. By finite propagation speed the map $\Lambda_{q}$ is unaffected if one changes $q$ outside the set 333If $u$ and $\tilde{u}$ solve (4.1) for potentials $q$ and $\tilde{q}$ with the same Dirichlet data $f$, and if $q=\tilde{q}$ in $U:=\{x\in\Omega\,;\,\mathrm{dist}(x,\partial\Omega)<T/2\}$, then $w:=u-\tilde{u}$ solves $(\Box+q)w=F$ where $F:=-(q-\tilde{q})\tilde{u}$ vanishes in $U\times(0,T)$ and in $(\Omega\setminus U)\times(0,T/2)$. Moreover, $w=\partial_{t}w=0$ on $\{t=0\}$ and $w|_{\partial\Omega\times(0,T)}=0$. By finite speed of propagation $\partial_{\nu}w|_{\partial\Omega\times(0,T)}=0$. This proves that $\Lambda_{q}=\Lambda_{\tilde{q}}$.

[TABLE]

For $T$ large enough, one can recover everything:

Corollary 4.2**.**

If $T>\mathrm{diam}(\Omega)$, then $\Lambda_{q_{1}}=\Lambda_{q_{2}}$ implies $q_{1}\equiv q_{2}$.

Proof.

If $T>\mathrm{diam}(\Omega)$, then by Theorem 4.1 one has

[TABLE]

for any maximal line segment $\gamma$ in $\overline{\Omega}$. Thus $q_{1}$ and $q_{2}$ have the same X-ray transform in $\mathbb{R}^{n}$. This transform is injective by Corollary 3.2 when $n=2$. Tiling $\mathbb{R}^{n}$ by two-planes gives injectivity when $n\geq 3$. Thus $q_{1}=q_{2}$. ∎

Theorem 4.1 could be proved based on the following facts, see e.g. [SY18]:

The map $\Lambda_{q}$ is an FIO of order $1$ on $\partial\Omega\times(0,T)$. 2. 2.

The X-ray transform of $q$ can be read off from the symbol of $\Lambda_{q}$ (more precisely, from the principal symbol of $\Lambda_{q}-\Lambda_{0}$).

We will give an elementary proof that is based on testing $\Lambda_{q}$ against highly oscillatory boundary data (compare with (1.8)).

The first step is an integral identity.

Lemma 4.3** (Integral identity).**

Assume that $q_{1},q_{2}\in C^{\infty}_{c}(\Omega)$. For any $f_{1},f_{2}\in C^{\infty}_{c}(\partial\Omega\times(0,T))$, one has

[TABLE]

where $u_{1}$ solves (4.1) with $q=q_{1}$ and $f=f_{1}$, and $u_{2}$ solves an analogous problem with vanishing Cauchy data on $\{t=T\}$:

[TABLE]

Proof.

We first compute the adjoint of the DN map: one has

[TABLE]

where $\Lambda_{q}^{T}g=\partial_{\nu}v|_{\partial\Omega\times(0,T)}$ with $v$ solving $(\Box+q)v=0$ so that $v|_{\partial\Omega\times(0,T)}=g$ and $v=\partial_{t}v=0$ on $\{t=T\}$. To prove this, we let $u$ be the solution of (4.1) and integrate by parts:

[TABLE]

Now, if $u_{1}$ and $u_{2}$ are as stated, the computation above gives

[TABLE]

The result follows by subtracting these two identities. ∎

The second step is to construct special solutions to the wave equation that concentrate near curves $s\mapsto(\gamma(s),s)$ where $\gamma$ is a line segment. These curves are projections to the $(x,t)$ variables of null bicharacteristic curves for $\Box$ (see Example 2.8). Thus the following result is in line with Theorem 2.7 concerning propagation of singularities. The proof is based on a standard geometrical optics / WKB quasimode construction.

Proposition 4.4** (Concentrating solutions).**

Assume that $q\in C^{\infty}_{c}(\Omega)$, and let $\gamma:[\delta,L]\to\overline{\Omega}$ be a maximal line segment in $\overline{\Omega}$ with $0<\delta<L<T$. For any $\lambda\geq 1$ there is a solution $u=u_{\lambda}$ of $(\Box+q)u=0$ in $\Omega\times(0,T)$ with $u=\partial_{t}u=0$ on $\{t=0\}$, such that for any $\psi\in C^{\infty}_{c}(\Omega\times[0,T])$ one has

[TABLE]

Moreover, if $\tilde{q}\in C^{\infty}_{c}(\Omega)$, there is a solution $\tilde{u}=\tilde{u}_{\lambda}$ of $(\Box+\tilde{q})\tilde{u}=0$ in $\Omega\times(0,T)$ with $\tilde{u}=\partial_{t}\tilde{u}=0$ on $\{t=T\}$, such that for any $\psi\in C^{\infty}_{c}(\Omega\times[0,T])$ one has

[TABLE]

At this point it is easy to prove the main result:

Proof of Theorem 4.1.

Using the assumption $\Lambda_{q_{1}}=\Lambda_{q_{2}}$ and Lemma 4.3, we have

[TABLE]

for any solutions $u_{j}$ of $(\Box+q_{j})u_{j}=0$ in $\Omega\times(0,T)$ so that $u_{1}=\partial_{t}u_{1}=0$ on $\{t=0\}$, and $u_{2}=\partial_{t}u_{2}=0$ on $\{t=T\}$.

Let $\gamma:[\delta,L]\to\overline{\Omega}$ be a maximal unit speed line segment in $\overline{\Omega}$ with $L<T$, and let $u_{1}=u_{1,\lambda}$ be the solution constructed in Proposition 4.4 for the potential $q_{1}$ with $u_{1}=\partial_{t}u_{1}=0$ on $\{t=0\}$. Moreover, let $u_{2}=u_{2,\lambda}$ be the solution constructed in the end of Proposition 4.4 for the potential $q_{2}$ with $u_{1}=\partial_{t}u_{1}=0$ on $\{t=T\}$. Taking the limit as $\lambda\to\infty$ in (4.5) and using (4.4) with $\psi(x,t)=(q_{1}-q_{2})(x)$, we obtain that

[TABLE]

Thus the integrals of $q_{1}$ and $q_{2}$ over maximal line segments of length $<T$ in $\overline{\Omega}$ are the same. ∎

Proof of Proposition 4.4.

Let $\gamma:[\delta,L]\to\overline{\Omega}$ be a maximal unit speed line segment in $\overline{\Omega}$ with $L<T$, and let $\eta:\mathbb{R}\to\mathbb{R}^{n}$ be the unit speed line so that $\eta(s)=\gamma(s)$ for $s\in[\delta,L]$. Write $x_{0}:=\eta(0)$ and $\xi_{0}:=\dot{\eta}(0)$, so that $x_{0}\notin\overline{\Omega}$ and $\gamma(s)=x_{0}+(s+\delta)\xi_{0}$. After a translation and rotation, we may assume that $x_{0}=0$ and $\xi_{0}=e_{n}$.

We first construct an approximate solution $v=v_{\lambda}$ for the operator $\Box+q$, having the form

[TABLE]

where $\varphi$ is a real phase function, and $a$ is an amplitude supported near the curve $s\mapsto(\eta(s),s)$. Note that

[TABLE]

Using a similar expression for $\partial_{x_{j}}^{2}$, we compute

[TABLE]

We would like to have $(\Box+q)(e^{i\lambda\varphi}a)=O(\lambda^{-1})$. To this end, we first choose $\varphi$ so that the $\lambda^{2}$ term in (4.6) vanishes. This will be true if $\varphi$ solves the eikonal equation

[TABLE]

There are many possible solutions, but we make the simple choice

[TABLE]

With this choice, (4.6) becomes

[TABLE]

where $L$ is the constant vector field

[TABLE]

It is convenient to consider new coordinates $(x^{\prime},z,w)$ in $\mathbb{R}^{n+1}$, where

[TABLE]

Then $L$ corresponds to $2\partial_{z}$ in the sense that

[TABLE]

where $\breve{F}$ corresponds to $F$ in the new coordinates:

[TABLE]

We next look for the amplitude $a$ in the form

[TABLE]

Inserting this to (4.6) and equating like powers of $\lambda$, we get

[TABLE]

We would like the last expression to be $O(\lambda^{-1})$. This will hold if $a_{0}$ and $a_{-1}$ satisfy the transport equations

[TABLE]

Let $\chi\in C^{\infty}_{c}(\mathbb{R}^{n})$ be supported near [math], and choose

[TABLE]

We will later choose $\chi$ to depend on $\lambda$. Next we choose

[TABLE]

These functions satisfy (4.12), and they vanish unless $w$ is small (i.e. $x_{n}$ is close to $t$). Then (4.9) becomes

[TABLE]

where

[TABLE]

Using the Cauchy-Schwartz inequality, one can check that

[TABLE]

uniformly over $\lambda\geq 1$. This concludes the construction of the approximate solution $v=e^{i\lambda\varphi}a$.

We next find an exact solution $u=u_{\lambda}$ of (4.1) having the form

[TABLE]

where $r$ is a correction term. Note that for $t$ close to [math], $v(\,\cdot\,,t)$ is supported near $x_{0}\notin\overline{\Omega}$ and hence $v=\partial_{t}v=0$ on $\{t=0\}$. Note also that $(\Box+q)v=F_{\lambda}$. Thus $u$ will solve (4.1) for $f=v|_{\partial\Omega\times(0,T)}$ if $r$ solves

[TABLE]

By the wellposedness of this problem [Ev10, Theorem 5 in §7.2.3], there is a unique solution $r$ with

[TABLE]

We now fix the choice of $\chi$ so that (4.3) will hold. Let $\zeta\in C^{\infty}_{c}(\mathbb{R}^{n})$ satisfy $\zeta=1$ near [math] and $\lVert\zeta\rVert_{L^{2}(\mathbb{R}^{n})}=1$, and choose

[TABLE]

where

[TABLE]

With this choice

[TABLE]

It follows that

[TABLE]

Since $u=v+r$, the integral in (4.3) has the form

[TABLE]

Using that $\psi\lvert a_{0}\rvert^{2}$ is compactly supported in $\Omega\times(0,T)$, we have

[TABLE]

by changing variables as in (4.8). Finally, changing $x^{\prime}$ to $\varepsilon x^{\prime}$ and $w$ to $\varepsilon w$ and letting $\lambda\to\infty$ (so $\varepsilon\to 0$) yields

[TABLE]

by the normalization $\lVert\zeta\rVert_{L^{2}(\mathbb{R}^{n})}=1$ and the fact that $\psi\in C^{\infty}_{c}(\Omega\times[0,T])$. This proves (4.3).

It remains to prove (4.4). Since $\eta(T)\notin\overline{\Omega}$, we have $v=\partial_{t}v=0$ on $\{t=T\}$, and we may alternatively arrange that $r$ solves (4.13) with $r=\partial_{t}r=0$ on $\{t=T\}$ instead of $\{t=0\}$. We can do such a construction for the potential $\tilde{q}$ instead of $q$. Since $\varphi$ and $a_{0}$ are independent of $q$, the same argument as above proves (4.4). ∎

5. Calderón problem: boundary determination

Electrical Impedance Tomography (EIT) is an imaging method with potential applications in medical imaging and nondestructive testing. The method is based on the following important inverse problem.

Calderón problem: Is it possible to determine the electrical conductivity of a medium by making voltage and current measurements on its boundary?

The treatment in this section follows [FSU].

Let us begin by recalling the mathematical model of EIT. The purpose is to determine the electrical conductivity $\gamma(x)$ at each point $x\in\Omega$, where $\Omega\subset\mathbb{R}^{n}$ represents the body which is imaged (in practice $n=3$). We assume that $\Omega\subset\mathbb{R}^{n}$ is a bounded open set with $C^{\infty}$ boundary, and that $\gamma\in C^{\infty}(\overline{\Omega})$ is positive.

Under the assumption of no sources or sinks of current in $\Omega$, a voltage potential $f$ at the boundary $\partial\Omega$ induces a voltage potential $u$ in $\Omega$, which solves the Dirichlet problem for the conductivity equation,

[TABLE]

Since $\gamma\in C^{\infty}(\overline{\Omega})$ is positive, the equation is uniformly elliptic, and there is a unique solution $u\in C^{\infty}(\overline{\Omega})$ for any boundary value $f\in C^{\infty}(\partial\Omega)$. One can define the Dirichlet-to-Neumann map (DN map) as

[TABLE]

Here $\nu$ is the outer unit normal to $\partial\Omega$ and $\partial_{\nu}u|_{\partial\Omega}=\nabla u\cdot\nu|_{\partial\Omega}$ is the normal derivative of $u$. Physically, $\Lambda_{\gamma}f$ is the current flowing through the boundary.

The Calderón problem (also called the inverse conductivity problem) is to determine the conductivity function $\gamma$ from the knowledge of the map $\Lambda_{\gamma}$. That is, if the measured current $\Lambda_{\gamma}f$ is known for all boundary voltages $f\in C^{\infty}(\partial\Omega)$, one would like to determine the conductivity $\gamma$.

We will prove the following theorem.

Theorem 5.1** (Boundary determination).**

Let $\gamma_{1},\gamma_{2}\in C^{\infty}(\overline{\Omega})$ be positive. If

[TABLE]

then the Taylor series of $\gamma_{1}$ and $\gamma_{2}$ coincide at any point of $\partial\Omega$.

This result was proved in [KV84], and it in particular implies that any real-analytic conductivity is uniquely determined by the DN map. The argument extends to piecewise real-analytic conductivities. A different proof was given in [SU88], based on two facts:

The DN map $\Lambda_{\gamma}$ is an elliptic $\Psi$DO of order $1$ on $\partial\Omega$. 2. 2.

The Taylor series of $\gamma$ at a boundary point can be read off from the symbol of $\Lambda_{\gamma}$ computed in suitable coordinates. The symbol of $\Lambda_{\gamma}$ can be computed by testing against highly oscillatory boundary data (compare with (1.8)).

Remark 5.2**.**

The above argument is based on studying the singularities of the integral kernel of the DN map, and it only determines the Taylor series of the conductivity at the boundary. The values of the conductivity in the interior are encoded in the $C^{\infty}$ part of the kernel, and different methods (based on complex geometrical optics solutions) are required for interior determination.

Let us start with a simple example:

Example 5.3** (DN map in half space is a $\Psi$DO).**

Let $\Omega=\mathbb{R}^{n}_{+}=\{x_{n}>0\}$, so $\partial\Omega=\mathbb{R}^{n-1}=\{x_{n}=0\}$. We wish to compute the DN map for the Laplace equation (i.e. $\gamma\equiv 1$) in $\Omega$. Consider

[TABLE]

Writing $x=(x^{\prime},x_{n})$ and taking Fourier transforms in $x^{\prime}$ gives

[TABLE]

Solving this ODE for fixed $\xi^{\prime}$ and choosing the solution that decays for $x_{n}>0$ gives

[TABLE]

We may now compute the DN map:

[TABLE]

Thus the DN map on the boundary $\partial\Omega=\mathbb{R}^{n-1}$ is just $\Lambda_{1}=\lvert D_{x^{\prime}}\rvert$ corresponding to the Fourier multiplier $\lvert\xi^{\prime}\rvert$. This shows that at least in this simple case, the DN map is an elliptic $\Psi$DO of order $1$.

We will now prove Theorem 5.1 by an argument that avoids showing that the DN map is a $\Psi$DO, but is rather based on directly testing the DN map against oscillatory boundary data. The first step is a basic integral identity (sometimes called Alessandrini identity) for the DN map.

Lemma 5.4** (Integral identity).**

Let $\gamma_{1},\gamma_{2}\in C^{\infty}(\overline{\Omega})$. If $f_{1},f_{2}\in C^{\infty}(\partial\Omega)$, then

[TABLE]

where $u_{j}\in C^{\infty}(\overline{\Omega})$ solves $\mathrm{div}(\gamma_{j}\nabla u_{j})=0$ in $\Omega$ with $u_{j}|_{\partial\Omega}=f_{j}$.

Proof.

We first observe that the DN map is symmetric: if $\gamma\in C^{\infty}(\overline{\Omega})$ is positive and if $u_{f}$ solves $\nabla\cdot(\gamma\nabla u_{f})=0$ in $\Omega$ with $u_{f}|_{\partial\Omega}=f$, then an integration by parts shows that

[TABLE]

Thus

[TABLE]

The result follows by subtracting the above two identities. ∎

Next we show that if $x_{0}$ is a boundary point, there is an approximate solution of the conductivity equation that concentrates near $x_{0}$, has highly oscillatory boundary data, and decays exponentially in the interior. As a simple example, the solution of

[TABLE]

that decays for $x_{n}>0$ is given by $u=e^{-\lambda x_{n}}e^{i\lambda x^{\prime}\cdot\xi^{\prime}}$, which concentrates near $\{x_{n}=0\}$ and decays exponentially when $x_{n}>0$ if $\lambda$ is large. Roughly, this means that the solution of a Laplace type equation with highly oscillatory boundary data concentrates near the boundary.

Proposition 5.5**.**

(Concentrating approximate solutions) Let $\gamma\in C^{\infty}(\overline{\Omega})$ be positive, let $x_{0}\in\partial\Omega$, let $\xi_{0}$ be a unit tangent vector to $\partial\Omega$ at $x_{0}$, and let $\chi\in C^{\infty}_{c}(\partial\Omega)$ be supported near $x_{0}$. Let also $N\geq 1$. For any $\lambda\geq 1$ there exists $v=v_{\lambda}\in C^{\infty}(\overline{\Omega})$ having the form

[TABLE]

such that

[TABLE]

and as $\lambda\to\infty$

[TABLE]

Moreover, if $\tilde{\gamma}\in C^{\infty}(\overline{\Omega})$ is positive and $\tilde{v}=\tilde{v}_{\lambda}$ is the corresponding approximate solution constructed for $\tilde{\gamma}$, then for any $f\in C(\overline{\Omega})$ and $k\geq 0$ one has

[TABLE]

for some $c_{k}\neq 0$.

We can now give the proof of the boundary determination result.

Proof of Theorem 5.1.

Using the assumption that $\Lambda_{\gamma_{1}}=\Lambda_{\gamma_{2}}$ together with the integral identity in Lemma 5.4, we have that

[TABLE]

whenever $u_{j}$ solves $\mathrm{div}(\gamma_{j}\nabla u_{j})=0$ in $\Omega$.

Let $x_{0}\in\partial\Omega$, let $\xi_{0}$ be a unit tangent vector to $\partial\Omega$ at $x_{0}$, and let $\chi\in C^{\infty}_{c}(\partial\Omega)$ satisfy $\chi=1$ near $x_{0}$. We use Proposition 5.5 to construct functions

[TABLE]

so that

[TABLE]

We obtain exact solutions $u_{j}$ of $\mathrm{div}(\gamma_{j}\nabla u_{j})=0$ by setting

[TABLE]

where the correction terms $r_{j}$ are the unique solutions of

[TABLE]

By standard energy estimates [Ev10, Section 6.2] and by (5.4), the solutions $r_{j}$ satisfy

[TABLE]

We now insert the solutions $u_{j}=v_{j}+r_{j}$ into (5.3). Using (5.5) and (5.4), it follows that

[TABLE]

as $\lambda\to\infty$. Letting $\lambda\to\infty$, the formula (5.2) yields

[TABLE]

In particular, $\gamma_{1}(x_{0})=\gamma_{2}(x_{0})$.

We will prove by induction that

[TABLE]

The case $j=0$ was proved above (here we may vary $x_{0}$ slightly). We make the induction hypothesis that (5.7) holds for $j\leq k-1$. Let $(x^{\prime},x_{n})$ be boundary normal coordinates so that $x_{0}$ corresponds to [math], and $\partial\Omega$ near $x_{0}$ corresponds to $\{x_{n}=0\}$. The induction hypothesis states that

[TABLE]

Considering the Taylor expansion of $(\gamma_{1}-\gamma_{2})(x^{\prime},x_{n})$ with respect to $x_{n}$ gives that

[TABLE]

for some smooth function $f$ with $f(x^{\prime},0)=\frac{\partial_{n}^{k}(\gamma_{1}-\gamma_{2})(x^{\prime},0)}{k!}$. Inserting this formula in (5.6), we obtain that

[TABLE]

Now $x_{n}=\mathrm{dist}(x,\partial\Omega)$ in boundary normal coordinates. Assuming that $N$ was chosen larger than $k$, we may take the limit as $\lambda\to\infty$ and use (5.2) to obtain that

[TABLE]

This shows that $\partial_{n}^{k}(\gamma_{1}-\gamma_{2})(x^{\prime},0)=0$ for $x^{\prime}$ near [math], which concludes the induction. ∎

It remains to prove Proposition 5.5, which constructs approximate solutions (also called quasimodes) concentrating near a boundary point. This is a typical geometrical optics / WKB type construction for quasimodes with complex phase. The proof is elementary, although a bit long. The argument is simplified slightly by using the Borel summation lemma, which is used frequently in microlocal analysis in various different forms.

Lemma 5.6** (Borel summation, [Hö85, Theorem 1.2.6]).**

Let $f_{j}\in C^{\infty}_{c}(\mathbb{R}^{n-1})$ for $j=0,1,2,\ldots$. There exists $f\in C^{\infty}_{c}(\mathbb{R}^{n})$ such that

[TABLE]

Proof of Proposition 5.5.

We will first carry out the proof in the case where $x_{0}=0$ and $\partial\Omega$ is flat near [math], i.e. $\Omega\cap B(0,r)=\{x_{n}>0\}\cap B(0,r)$ for some $r>0$ (the general case will be considered in the end of the proof). We also assume $\xi_{0}=(\xi_{0}^{\prime},0)$ where $\lvert\xi_{0}^{\prime}\rvert=1$.

We look for $v$ in the form

[TABLE]

Write $Pu=D\cdot(\gamma Du)=\gamma D^{2}u+D\gamma\cdot Du$. The principal symbol of $P$ is

[TABLE]

Since $e^{-i\lambda\Phi}D_{j}(e^{i\lambda\Phi}b)=(D_{j}+\lambda\partial_{j}\Phi)b$, we compute

[TABLE]

We want to choose $\Phi$ and $b$ so that $P(e^{i\lambda\Phi}b)=O_{L^{2}(\Omega)}(\lambda^{-N})$. Looking at the $\lambda^{2}$ term in (5.9), we first choose $\Phi$ so that

[TABLE]

We additionally want that $\Phi(x^{\prime},0)=x^{\prime}\cdot\xi_{0}^{\prime}$ and $\partial_{n}\Phi(x^{\prime},0)=i$ (this will imply that $\nabla\Phi(0)=\xi_{0}+ie_{n}$). In fact, using (5.8) we can just choose

[TABLE]

and then $p_{2}(x,\nabla\Phi)=\gamma(\xi_{0}+ie_{n})\cdot(\xi_{0}+ie_{n})\equiv 0$ in $\Omega$.

We next look for $b$ in the form

[TABLE]

Since $p_{2}(x,\nabla\Phi)\equiv 0$, (5.9) implies that

[TABLE]

We will choose the functions $b_{-j}$ so that

[TABLE]

We will additionally arrange that

[TABLE]

and that each $b_{-j}$ is compactly supported so that

[TABLE]

for some fixed $\varepsilon>0$.

To find $b_{0}$, we prescribe $b_{0}(x^{\prime},0),\partial_{n}b_{0}(x^{\prime},0),\partial_{n}^{2}b_{0}(x^{\prime},0)$, $\ldots$ successively and use the Borel summation lemma to construct $b_{0}$ with this Taylor series at $\{x_{n}=0\}$. We first set $b_{0}(x^{\prime},0)=\chi(x^{\prime})$. Writing $\eta:=\nabla\cdot(\gamma\nabla\Phi)$, we observe that

[TABLE]

Thus, in order to have $Lb_{0}|_{x_{n}=0}=0$ we must have

[TABLE]

We prescribe $\partial_{n}b(x^{\prime},0)$ to have the above value (which depends on the already prescribed quantity $b(x^{\prime},0)$). Next we compute

[TABLE]

where $Q$ depends on the already prescribed quantities $b_{0}(x^{\prime},0)$ and $\partial_{n}b_{0}(x^{\prime},0)$. We thus set

[TABLE]

which ensures that $\partial_{n}(Lb_{0})|_{x_{n}=0}=0$. Continuing in this way and using Borel summation, we obtain a function $b_{0}$ so that $Lb_{0}=0$ to infinite order at $\{x_{n}=0\}$. The other equations in (5.16) are solved in a similar way, which gives the required functions $b_{-1},\ldots,b_{-N}$. In the construction, we may arrange so that (5.19) and (5.20) are valid.

If $\Phi$ and $b_{-j}$ are chosen in the above way, then (5.11) implies that

[TABLE]

where each $q_{j}(x)$ vanishes to infinite order at $x_{n}=0$ and is compactly supported in $Q_{\varepsilon}$. Thus, for any $k\geq 0$ there is $C_{k}>0$ so that $\lvert q_{j}\rvert\leq C_{k}x_{n}^{k}$ in $Q_{\varepsilon}$, and consequently

[TABLE]

Since $\mathrm{Im}(\Phi)=x_{n}$ in $Q_{\varepsilon}$ we have

[TABLE]

Choosing $k=N+1$ and computing the integrals over $x_{n}$, we get that

[TABLE]

It is also easy to compute that

[TABLE]

Thus, choosing $a=\lambda^{-1/2}b$, we have proved all the claims except (5.2).

To show (5.2), we observe that

[TABLE]

Using a similar formula for $\tilde{v}=e^{i\lambda\Phi}\tilde{a}$ (where $\Phi$ is independent of the conductivity), we have

[TABLE]

Now $\lvert\nabla\Phi\rvert^{2}=2$ and $a=\lambda^{-1/2}b$ where $\lvert b\rvert\lesssim 1$, and similarly for $\tilde{a}$. Hence

[TABLE]

We can change variables $x_{n}\to x_{n}/\lambda$ and use dominated convergence to take the limit as $\lambda\to\infty$. The limit is

[TABLE]

where $c_{k}=2\int_{0}^{\infty}x_{n}^{k}e^{-2x_{n}}\,dx_{n}\neq 0$.

The proof is complete in the case when $x_{0}=0$ and $\partial\Omega$ is flat near [math]. In the general case, we choose boundary normal coordinates $(x^{\prime},x_{n})$ so that $x_{0}$ corresponds to [math] and $\Omega$ near $x_{0}$ locally corresponds to $\{x_{n}>0\}$. The equation $\nabla\cdot(\gamma\nabla u)=0$ in the new coordinates becomes an equation

[TABLE]

where $A$ is a smooth positive matrix only depending on the geometry of $\Omega$ near $x_{0}$. The construction of $v$ now proceeds in a similary way as above, except that the equation (5.10) for the phase function $\Phi$ can only be solved to infinite order on $\{x_{n}=0\}$ instead of solving it globally in $\Omega$. ∎