*Article* **The Dirichlet-to-Neumann Map in a Disk with a One-Step Radial Potential: An Analytical and Numerical Study**

**Sagrario Lantarón and Susana Merchán \***

Departamento de Matemática e Informática Aplicadas a las Ingenierías Civil y Naval, Escuela de Caminos, Canales y Puertos, Universidad Politécnica de Madrid, Calle del Profesor Aranguren, 3, 28040 Madrid, Spain; sagrario.lantaron@upm.es

**\*** Correspondence: susana.merchan@upm.es

**Abstract:** Herein, we considered the Schrödinger operator with a potential *q* on a disk and the map that associates to *q* the corresponding Dirichlet-to-Neumann (DtN) map. We provide some numerical and analytical results on the range of this map and its stability for the particular class of one-step radial potentials.

**Keywords:** Dirichlet-to-Neumann map; Schrödinger operator; stability

#### **1. Introduction**

Let <sup>Ω</sup> <sup>⊂</sup> <sup>R</sup><sup>2</sup> be a bounded domain with smooth boundary *<sup>∂</sup>*Ω. For each *<sup>q</sup>* <sup>∈</sup> *<sup>L</sup>*∞(Ω), consider the so called Dirichlet-to-Neumann map (DtN) given by:

$$\begin{array}{rcl} \Lambda\_q: H^{1/2}(\partial \Omega) & \to & H^{-1/2}(\partial \Omega) \\ & f & \to & \frac{\partial \mu}{\partial n}|\_{\partial \Omega} \end{array} \tag{1}$$

where *u* is the solution of the following problem:

$$\begin{cases} \Delta u + q(\mathbf{x})u = 0, & \mathbf{x} \in \Omega, \\ u = f\_{\prime} & \partial \Omega\_{\prime} \end{cases} \tag{2}$$

and *<sup>∂</sup><sup>u</sup> <sup>∂</sup><sup>n</sup>* |*∂*<sup>Ω</sup> denotes the normal derivative of *u* on the boundary *∂*Ω.

Note that the uniqueness of *u* as solution of (2) requires that 0 is not a Dirichlet eigenvalue of <sup>Δ</sup> <sup>+</sup> *<sup>q</sup>*. A sufficient condition to guarantee that <sup>Λ</sup>*<sup>q</sup>* is well defined is to assume *<sup>q</sup>*(*x*) <sup>&</sup>lt; *<sup>λ</sup>*1, the first Dirichlet eigenvalue of the Laplace operator in <sup>Ω</sup>, since, in this case, the solution in (2) is unique. We assume that this condition holds and lets us define the space

$$L\_{<\lambda\_1}^{\infty}(\Omega) = \{ q \in L^{\infty}(\Omega) \text{ . s.t. } q(\mathbf{x}) < \lambda\_1 \text{ a.e. } \}.$$

In this work, we were interested in the following map:

$$\begin{array}{rcl}\Lambda : L^{\infty}\_{<\lambda\_1}(\Omega) & \to & \mathcal{L}(H^{1/2}(\partial\Omega); H^{-1/2}(\partial\Omega)) \\ \neq & \Lambda\_q. \end{array} \tag{3}$$

This has an important role in inverse problems, where the aim is to recover the potential *q* from boundary measurements. In practice, these boundary measurements correspond to the associated DtN map, and therefore, the mathematical statement of the classical inverse problem consists of the inversion of Λ.

It is known that <sup>Λ</sup> is one-to-one as long as *<sup>q</sup>* <sup>∈</sup> *<sup>L</sup><sup>p</sup>* with *<sup>p</sup>* <sup>&</sup>gt; 2 (see [1]). Therefore, the inverse map Λ−<sup>1</sup> can be defined in the range of Λ. There are, however, two related important and difficult questions that are not well understood: A characterization of the

**Citation:** Lantarón, S.; Merchán, S. The Dirichlet-to-Neumann Map in a Disk with a One-Step Radial Potential: An Analytical and Numerical Study. *Mathematics* **2021**, *9*, 794. https:// doi.org/10.3390/math9080794

Academic Editor: Lucas Jódar

Received: 2 March 2021 Accepted: 25 March 2021 Published: 7 April 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

range of Λ and its stability, i.e., a quantification of the difference of two potentials, in the *L*<sup>∞</sup> topology in terms of the distance of their associated DtN maps. Obviously, this stability will affect the efficiency of any inversion or reconstruction algorithm to recover the potential from the DtN map (see [2] and [3].

The first question, i.e., the characterization of the range of Λ is widely open. Tothe best of our knowledge, the further result is due to [4], where a characterization is obtained for the adherence, with respect to the weak topology in <sup>2</sup> <sup>−</sup>1, of the sequence of eigenvalues associated with the orthogonal basis of eigenvectors {*eink*}*k*∈Z. Here, <sup>2</sup> *<sup>α</sup>* is the space of sequences {*cn*}*n*∈Z, such that ∑*n*∈<sup>Z</sup> |*n*| <sup>2</sup>*α*|*cn*<sup>|</sup> <sup>2</sup> < ∞. This topology is not the usual one in <sup>L</sup>(*H*1/2(*∂*Ω); *<sup>H</sup>*−1/2(*∂*Ω)) and it is not easy to interpret how the adherence enlarges this set. Furthermore, the characterization does not give much practical information on the range, as, for instance, the convexity or accurate bounds on the eigenvalues. In fact, characterizing such properties is one of the main motivations of this work, since we could establish easily a priori if a desired linear map in <sup>L</sup>(*H*1/2(*∂*Ω); *<sup>H</sup>*−1/2(*∂*Ω)) can be associated with a DtN map. On the contrary, we have to take into account that in practice, the DtN map is estimated from physical measurements, which are subject to errors and may provide only partial information. A precise knowledge of the range of Λ is useful to find the best DtN map that fits the measurements and to design an inversion algorithm in such situations.

Concerning the stability, it is well known that the problem is ill posed and that the most we can expect is logarithmic stability in general (see [5]). There are more explicit results when we assume that the potential *<sup>q</sup>* has some smoothness. In particular, if *<sup>q</sup>* <sup>∈</sup> *<sup>H</sup>s*(Ω) with *s* > 0, the following log −stability condition is known (see [1]):

$$\|\|q\_1 - q\_2\|\|\_{L^\infty} \le V(\|\Lambda\_{\mathbb{q}\_1} - \Lambda\_{\mathbb{q}\_2}\|\_{\mathcal{L}(H^{1/2}; H^{-1/2})}),\tag{4}$$

where *<sup>V</sup>*(*t*) = *<sup>C</sup>* log(1/*t*)−*<sup>α</sup>* for some constants *<sup>C</sup>*, *<sup>α</sup>* <sup>&</sup>gt; 0. Stronger stability conditions are known in some particular cases. For example, in [6], it was shown that when *q* is piecewise constant and the components where it takes a constant value are fixed and satisfy some technical conditions, the stability is Lipschitz, i.e., there exists a constant *C* > 0, such that:

$$\|\|q\_1 - q\_2\|\|\_{L^{\mathfrak{su}}} \le \mathbb{C} \|\Lambda\_{\mathfrak{q}1} - \Lambda\_{\mathfrak{q}2}\|\_{\mathcal{L}(H^{1/2}; H^{-1/2})}.\tag{5}$$

In this work, we tried to understand better the situation by considering the simplest case of a disk with one-step radial potentials *q*. More precisely, we provide some results on the range of <sup>Λ</sup> and its stability when we restrict to the particular case <sup>Ω</sup> = *<sup>B</sup>*(0, 1) = *<sup>x</sup>* <sup>∈</sup> <sup>R</sup><sup>2</sup> : *<sup>r</sup>* <sup>=</sup> <sup>|</sup>*x*<sup>|</sup> <sup>&</sup>lt; <sup>1</sup> and *<sup>q</sup>* <sup>∈</sup> *<sup>F</sup>* <sup>⊂</sup> *<sup>L</sup>*∞(Ω) given by:

$$F = \{ q \in L^{\infty}(\Omega) \; : \; q(r) = \gamma \chi\_{(0,b)}(r), \; r = |\mathbf{x}|, \; b \in (0,1), \; \gamma \in [0,1] \},\tag{6}$$

where *<sup>χ</sup>*(0,*b*)(*r*) is the characteristic function of the interval (0, *<sup>b</sup>*). Note that *<sup>F</sup>* is a twoparametric family depending on *γ* and *b*.

It is worth mentioning that, as we show below, the solution of (2) is unique for all *<sup>b</sup>* <sup>∈</sup> (0, 1) and *<sup>γ</sup>* <sup>≥</sup> 0, and therefore, the DtN map is well defined for all of these one-step potentials. In other words, 0 is not an eigenvalue of the operator <sup>Δ</sup> + *<sup>q</sup>* and, in particular, we do not need to restrict ourselves to the constraint *<sup>q</sup>*(*x*) <sup>&</sup>lt; *<sup>λ</sup>*1. However, we still restrict ourselves to the bounded set *F* to simplify.

Even in this simple case, a complete analytic answer to the previous questions (range of the DtN map and sharp stability conditions) is unknown. Therefore, we also considered a numerical approach based on a discrete sampling of the set *F*. Given an integer *N* > 0, we define *<sup>h</sup>* = 1/*<sup>N</sup>* and:

$$F\_h \quad = \{ q \in L^\infty(\Omega) \, : \, q(r) = \gamma \chi\_{(0,b)}(r), \, b = h i, \, \gamma = h j, \, \} \tag{7}$$

$$i = 1, \dots, N - 1, \, j = 0, \dots, N \}. \tag{7}$$

Note that *Fh* has *<sup>N</sup>*(*<sup>N</sup>* <sup>−</sup> <sup>1</sup>) + 1 functions from *<sup>F</sup>*. As *<sup>h</sup>* <sup>→</sup> 0, we can obtain a better description of *F* and, in particular, we should recover the properties for *q* ∈ *F*.

The main contributions of this paper are given below:


$$G\_{\gamma} = \{ q \in L^{\infty}(\Omega) \; : \; q(r) = \gamma \chi\_{(0,b)}(r), \; b \in (0,1) \}. \tag{8}$$

We prove that if *b* ≥ *b*<sup>0</sup> > 0, there is stability of the DtN map with respect to the position of the discontinuity *b* for potentials in *Gγ*. More precisely, we obtain a stability constant depending on *γ*−1*b*−3, which is uniformly bounded for *b* > *b*<sup>0</sup> and fixed *γ* (see Theorem 3 below). Note, however, that this constant blows up as *γ* → 0. This stability result does not give information about the stability with respect to the *L*<sup>∞</sup> norm of the potentials, but it provides stability with respect to the *L*<sup>1</sup> norm, which is sensitive to the position of the discontinuities, when *b* > *b*<sup>0</sup> > 0 and *γ* > *γ*<sup>0</sup> > 0. In fact, we show numerical evidence of such stability when considering potentials in *F*.

4. For the range of Λ, we give a characterization in terms of the first two eigenvalues of the DtN map. We also analyze the region where the stability constant is larger, and, therefore, the potentials for which any recovering algorithm for *q* from the DtN map will have more difficulties.

We mention that a similar analysis can be conducted for the closely related—and more classical—conductivity problem, where (2) is replaced by:

$$\begin{cases} -\operatorname{div} a(\mathbf{x}) \nabla v = 0, & \mathbf{x} \in \Omega, \\ v = f, & \partial \Omega, \end{cases} \tag{9}$$

and the Dirichlet-to-Neumann map, or voltage-to-current map, is given by:

$$\begin{array}{rcl} \Lambda\_{\mathfrak{a}} : H^{1/2}(\partial \Omega) & \to & H^{-1/2}(\partial \Omega) \\ & f & \to & a \frac{\partial v}{\partial n}|\_{\partial \Omega} \end{array} \tag{10}$$

In this case, the relationship between piecewise constant radial conductivities and the eigenvalues of the DtN map is known [8] through a suitable recurrence formula. However, there is not a direct transformation that relates both problems, and the analysis must be done specifically for this case.

We refer to the review paper [9] and the references therein for theoretical results on the DtN map in this case.

The rest of this paper is divided as follows: In Section 2 below, we characterize the DtN map in terms of its eigenvalues using polar coordinates. In Sections 3 and 4, we analyze the stability and range results, respectively. In Section 5, we briefly describe the main conclusions, and finally, Section 5 contains the proofs of the theorems stated in the previous sections.

#### **2. The Dirichlet-to-Neumann Map**

In this section, we characterize the Dirichlet-to-Neumann map in the case of a disk. System (2) in polar coordinates reads:

$$\begin{cases} \begin{aligned} r^2 \frac{\partial^2 v}{\partial r^2} + r \frac{\partial v}{\partial r} + \frac{\partial^2 v}{\partial \theta^2} + r^2 q(r) v = 0, & (r, \theta) \in (0, 1) \times [0, 2\pi), \\ \lim\_{r \to 0, r > 0} v(r, \theta) < \infty, & \end{aligned} \end{cases} \tag{11}$$

where *<sup>v</sup>*(*r*, *<sup>θ</sup>*) = *<sup>u</sup>*(*<sup>r</sup>* cos *<sup>θ</sup>*,*<sup>r</sup>* sin *<sup>θ</sup>*) and *<sup>g</sup>*(*θ*) = *<sup>f</sup>*(cos *<sup>θ</sup>*, sin *<sup>θ</sup>*) is a periodic function.

An orthonormal basis in *<sup>L</sup>*2(0, 2*π*) is given by {*einθ*}*n*∈Z. Here, we use this complex basis to simplify the notation, but in the analysis below, we only consider the subspace of real valued functions. Therefore, any function *<sup>g</sup>* <sup>∈</sup> *<sup>L</sup>*2(0, 2*π*) can be written as:

$$\log(\theta) = \sum\_{n \in \mathbb{Z}} g\_n e^{in\theta}, \quad g\_n = \frac{1}{2\pi} \int\_0^{2\pi} \mathfrak{g}(t) e^{-int} dt,\tag{12}$$

and *<sup>g</sup>* <sup>2</sup> *<sup>L</sup>*2(0,2*π*) <sup>=</sup> <sup>∑</sup>*n*∈<sup>Z</sup> <sup>|</sup>*gn*<sup>|</sup> 2. Associated with this basis, we define the usual Hilbert spaces: *H<sup>α</sup>* # , for *α* > 0, as

$$H^{\mathfrak{a}}\_{\mathfrak{b}} = \{ \mathfrak{g} \; : \; ||\mathfrak{g}||\_{\mathfrak{a}}^2 = \sum\_{n \in \mathbb{Z}} (1 + n^2)^{\mathfrak{a}} |\mathfrak{g}\_n|^2 < \infty \}.$$

The Dirichlet-to-Neumann map in this case is defined as:

$$\begin{array}{rcl} \Lambda\_{\emptyset}: H\_{\emptyset}^{1/2}(0, 2\pi) & \to & H\_{\emptyset}^{-1/2}(0, 2\pi) \\ & \text{g} & \to & \frac{\partial v}{\partial r}(1, \cdot), \end{array} \tag{13}$$

where *v* is the unique solution of (11).

In the above basis, the Dirichlet-to-Neumann map turns out to be diagonal. In fact, we have the following result:

**Theorem 1.** *Let* Ω *be the unit disk and q* ∈ *F. Then:*

$$
\Lambda\_{\emptyset} \left( \epsilon^{in\theta} \right)\_{\prime} = \quad c\_n e^{in\theta}, \quad n \in \mathbb{Z}, \tag{14}
$$

*where:*

$$\begin{array}{rcl} c\_0 &=& \frac{-b\sqrt{\gamma}l\_1(\sqrt{\gamma}b)}{b\log b\sqrt{\gamma}l\_1(\sqrt{\gamma}b) + l\_0(\sqrt{\gamma}b)} \\ &=& \frac{-b\sqrt{\gamma}l\_1(\sqrt{\gamma}b) + l\_0(\sqrt{\gamma}b)}{b\log b} \end{array} \tag{15}$$

$$\mathcal{L}\_{n} = \quad \mathcal{L}\_{-n} = n \frac{J\_{n-1}(\sqrt{\gamma}b) - b^{2n} J\_{n+1}(\sqrt{\gamma}b)}{J\_{n-1}(\sqrt{\gamma}b) + b^{2n} J\_{n+1}(\sqrt{\gamma}b)}, \quad n \in \mathbb{N},\tag{16}$$

*and Jn*(*r*) *are the Bessel functions of the first kind.*

Note that the range of Λ, when restricted to *F*, is characterized by the set of sequences {*cn*}*n*≥<sup>0</sup> of the form (15) and (16) for all possible *<sup>b</sup>*, *<sup>γ</sup>*. In particular, when *<sup>q</sup>* <sup>=</sup> 0, we have:

$$\mathfrak{c}\_{\mathfrak{n}} = \mathfrak{n}, \quad \mathfrak{n} = 0, 1, 2, \dots, \tag{17}$$

and this sequence must be in the range of Λ.

The norm of Λ*q*, when restricted to *F*, is given by:

$$\|\Lambda\_{\mathfrak{q}}\|\_{\mathcal{L}(H^{1/2}\_{\mathfrak{q}};H^{-1/2}\_{\mathfrak{q}})} = \sup\_{n\geq 0} \frac{|c\_n|}{1+n^2}.\tag{18}$$

**Proof of Theorem 1.** We first compute *<sup>c</sup>*<sup>0</sup> in (14). As the boundary data at *<sup>r</sup>* <sup>=</sup> 1 in (11) is the constant *<sup>g</sup>*(*θ*) = 1, we assume that *<sup>v</sup>*(*r*, *<sup>θ</sup>*) is radial, i.e., *<sup>v</sup>*(*r*, *<sup>θ</sup>*) = *<sup>a</sup>*0(*r*). Then, *<sup>a</sup>*<sup>0</sup> should satisfy:

$$\begin{cases} \, \, ^2a\_0^{\prime\prime} + ra\_0^{\prime} + r^2 q(r) a\_0 = 0, \quad 0 < r < 1, \\\, \, a\_0(1) = 1, \quad \lim\_{r \to \pm 0} a\_0(r) < \infty. \end{cases} \tag{19}$$

For *<sup>r</sup>* <sup>∈</sup> (0, *<sup>b</sup>*), we solve the ODE with the boundary data at *<sup>r</sup>* = 0, while for *<sup>r</sup>* <sup>∈</sup> (*b*, 1), we use the boundary data at *<sup>r</sup>* = 1. In the first case, the ODE is the Bessel ODE of order 0, and therefore, we have:

$$a\_0(r) = \begin{cases} A\_0 l\_0(\sqrt{\gamma}r), & r \in (0, b), \\ 1 + \mathcal{C}\_0 \log r, & r \in (b, 1), \end{cases}$$

where *J*<sup>0</sup> is the Bessel function of the first kind and *A*<sup>0</sup> and *C*<sup>0</sup> are constants. These are computed by imposing continuity of *a*<sup>0</sup> and *a* <sup>0</sup> at *<sup>r</sup>* <sup>=</sup> *<sup>b</sup>*. In this way, we obtain:

$$\begin{cases} A\_0 l\_0(\sqrt{\gamma}b) = 1 + \mathcal{C}\_0 \log|b| \\ A\_0 \sqrt{\gamma} l\_0'(\sqrt{\gamma}b) = \mathcal{C}\_0 \frac{1}{b} . \end{cases}$$

Solving the system for *<sup>A</sup>*<sup>0</sup> and *<sup>C</sup>*<sup>0</sup> and taking into account that <sup>Λ</sup>*q*(1) = *<sup>∂</sup><sup>v</sup> <sup>∂</sup><sup>r</sup>* (1, *<sup>θ</sup>*) = *a* <sup>0</sup>(1) = *<sup>C</sup>*0, we easily obtain (14).

Similarly, to compute *cn* in (14), we have to consider *<sup>g</sup>*(*θ*) = *<sup>e</sup>in<sup>θ</sup>* in (11), and therefore, we assume that the solution *<sup>v</sup>*(*r*, *<sup>θ</sup>*) can be written in separate variables, i.e., *<sup>v</sup>*(*r*, *<sup>θ</sup>*) = *an*(*r*)*einθ*. Then, *an* must satisfy:

$$\begin{cases} r^2 a\_n^{\prime\prime} + r a\_n^{\prime} + (r^2 q(r) - n^2) a\_n = 0, & 0 < r < 1, \\\ a\_n(1) = 1, & \lim\_{r \to 0, r > 0} a\_n(r) < \infty, \quad n \ge 1. \end{cases} \tag{20}$$

As in the case of *<sup>c</sup>*0, for *<sup>r</sup>* <sup>∈</sup> (0, *<sup>b</sup>*), we solve the ODE with the boundary data at *<sup>r</sup>* = 0, while for *<sup>r</sup>* <sup>∈</sup> (*b*, 1), we use the boundary data at *<sup>r</sup>* = 1. We have:

$$a\_n(r) = \begin{cases} A\_n \mathcal{J}\_n(\sqrt{\gamma}r), & r \in (0, b), \\ C\_n(r^n - r^{-n}) + r^n, & r \in (b, 1), \end{cases}$$

where *An* and *Cn* are constants. These are computed by imposing continuity of *an* and *a n* at *<sup>r</sup>* = *<sup>b</sup>*. In this way, we obtain:

$$\begin{cases} A\_n l\_n(\sqrt{\gamma}b) = \mathbb{C}\_n(b^n - b^{-n}) + b^n\\ A\_n \sqrt{\gamma} l\_n'(\sqrt{\gamma}b) = n \mathbb{C}\_n(b^{n-1} + b^{-n-1}) + nb^{n-1} \dots \end{cases}$$

Solving the system for *An* and *Cn*, we obtain, in particular:

$$\mathbb{C}\_{n} = \frac{-b^{n}l\_{n}'(\sqrt{\gamma}b) + n\frac{b^{n-1}}{\sqrt{\gamma}}l\_{n}(\sqrt{\gamma}b)}{-(b^{-n-1} + b^{n-1})\frac{n}{\sqrt{\gamma}}l\_{n}(\sqrt{\gamma}b) - (b^{-n} - b^{n})l\_{n}'(\sqrt{\gamma}b)}.$$

We simplify this expression using the well-known identity:

$$2f\_n'(r) = f\_{n-1}(r) - f\_{n+1}(r),$$

and we obtain:

$$\mathbb{C}\_{\mathfrak{n}} = \frac{-J\_{\mathfrak{n}+1}(\sqrt{\gamma}b)}{b^{-2\mathfrak{n}}f\_{\mathfrak{n}-1}(\sqrt{\gamma}b) + f\_{\mathfrak{n}+1}(\sqrt{\gamma}b)}.$$

Now, taking into account that <sup>Λ</sup>*q*(*ein<sup>θ</sup>* ) = *<sup>∂</sup><sup>v</sup> <sup>∂</sup><sup>r</sup>* (1, *<sup>θ</sup>*) = *<sup>a</sup> <sup>n</sup>*(1)*ein<sup>θ</sup>* = (2*nCn* <sup>+</sup> *<sup>n</sup>*)*einθ*, we easily obtain (14).

**Remark 1.** *In this proof of Theorem 1, we do not use the restriction γ* ≤ 1 *that satisfies the potentials in F. In fact, the statement of the theorem still holds for any step potential, as in F, but with any arbitrary large γ* ≥ 0*.*

#### **3. Stability**

In this section, we focus on the stability results for the map Λ. Some results are analytical and they are stated as theorems. The proofs are given in Appendix A. We divided this section in three subsections, where we consider the negative stability result for *q* ∈ *F* norm, and some partial results when we consider the subsets *Fb* defined by:

$$F\_b = \left\{ q \in L^{\infty}(\Omega) \; : \; q(r) = \mathfrak{z}\_{(0,b)}(r), \quad \gamma \in [0,1] \right\},$$

and *Gγ* defined in (8).
