1 Introduction

Let (Mg) be a compact Riemannian manifold with piecewise smooth boundary \(\partial M\). The Steklov problem is given by

$$\begin{aligned} {\left\{ \begin{array}{ll} -\Delta _g \phi _{\sigma }=0&{}\hbox {in }M\\ \partial _\nu \phi _{\sigma }=\sigma \phi _{\sigma }&{}\hbox {on }\partial M. \end{array}\right. } \end{aligned}$$
(1.1)

There is a discrete sequence \(0=\sigma _0<\sigma _1\le \sigma _2\le \cdots \) of values of \(\sigma \), with \(\sigma _j\rightarrow \infty \) as \(j\rightarrow \infty \), for which non-trivial solutions satisfying (1.1) exist [9]. These are the Steklov eigenvalues and the corresponding functions \(\phi _{\sigma _j}\) are the Steklov eigenfunctions. This paper studies the asymptotic character of the nodal set of the eigenfunctions of the Steklov eigenvalue problem in the case M equals a bounded open set \(\Omega \in {\mathbb {R}}^2\). In particular the results in this paper show that the nodal set of the eigenfunction \(\phi _{\sigma _j}\) is not dense at scale \(\sigma _j^{-1}\) for some such sets \(\Omega \)—or, more precisely, that there is a dense family \(\mathcal {A}\) of simply-connected two-dimensional domains with analytic boundaries such that, for each \(\Omega \in \mathcal {A}\), the eigenfunction \(\phi _{\sigma _j}\) in the domain \(\Omega \) remains nonzero on a j-dependent ball of j-independent radius. This result addresses a question put forth under “Open Problem 10” in [7].

The behavior of both the Steklov eigenvalues (see e.g. [7, 8, 11]) and eigenfunctions (see e.g. [3, 5, 9, 13,14,15, 18, 20]) have been a topic of recent interest. When M has smooth boundary, the Steklov eigenfunctions \(\phi _{\sigma _j}|_{\partial M}\) behave much like high energy Laplace eigenfunctions with eigenvalue \(\sigma _j^2\). In particular, they oscillate at frequency \(\sigma _j\). References [3, 6, 13, 15, 17,18,19,20] study the nodal sets of \(\phi _{\sigma _j}|_M\), giving both upper and lower bounds on its Hausdorff measure similar to those for Laplace eigenfunctions. In fact, most results regarding Steklov eigenfunctions in the interior of M extract behavior similar to that of high energy Laplace eigenfunctions.

The purpose of this article is to show that, away from the boundary of M, Steklov eigenfunctions behave very differently than high energy Laplace eigenfunctions. Not only do they decay rapidly (see [5, 9]) but, at least for a dense class of analytic domains, they oscillate slowly over certain portions of the domain. Girouard–Polterovich [7, Open Problem 10(i)] raise the question of whether nodal sets of Steklov eigenfunctions are dense at scale \(\sigma _j^{-1}\) in M. One consequence of the results in the present paper is a negative answer to this question. We show that arbitrarily close to any simply-connected domain with analytic boundary \(\Omega _0\subset {\mathbb {R}}^2\), there is a domain \(\Omega _1\) for which the nodal sets are not\(\sigma _j^{-1}\) dense and, indeed, that there is a region within \(\Omega _1\) where the nodal set density does not increase as \(\sigma _j\rightarrow \infty \). Moreover, the Steklov eigenfunctions oscillate no faster than a fixed frequency in this region. These results are summarized in the following theorem.

Theorem 1

Let \(\Omega _0\subset {\mathbb {R}}^2\) be a bounded simply-connected domain with analytic boundary, and let \(k>0\) and \(\varepsilon >0\) be given. Then there exist a set \(\Omega _1\subset {\mathbb {R}}^2\) with analytic boundary given by

$$\begin{aligned} \partial \Omega _1=\{ x+{\nu } g(x)\mid x\in \partial \Omega _0\},\qquad \qquad \Vert g\Vert _{C^k(\partial \Omega _0)}<\varepsilon \end{aligned}$$
(1.2)

(where \({\nu }\) denotes the outward unit normal to \(\partial \Omega _0\) and where g is an analytic function defined on \(\partial \Omega _0\)), a point \(x_0\in \Omega _1\) and numbers \(0<r_1<r_0\), (\(B(x_0, r_0)\subset \Omega _1\)) such that: for each Steklov eigenvalue \(\sigma \) for the domain \(\Omega _1\) there exists a point \(x_\sigma \in B(x_0,r_0)\) such that \(B(x_\sigma , r_1)\subset B(x_0, r_0)\) and each Steklov eigenfunction \(\phi _\sigma \) of eigenvalue \(\sigma \) for the domain \(\Omega _1\) satisfies

$$\begin{aligned} |\phi _{\sigma }|>0\hbox { on }B(x_\sigma , r_1)\subset \Omega _1. \end{aligned}$$

Additionally, “\(\phi _\sigma \) has bounded frequency on \(B(x_0,r_0)\)” (a precise statement follows in Theorem 2).

Theorem 1 is a consequence of the more precise results presented in Theorems 2 and 3 and Corollary 2.2. In particular, these results establish that, for each domain \(\Omega \) in a dense class \({\mathcal {A}}\) of two-dimensional domains, an estimate holds for the truncation error in certain “mapped Fourier expansions” of the eigenfunctions \(\phi _\sigma \) (i.e., Fourier expansions of \(\phi _\sigma \) under a change of variables). This estimate is uniformly valid over a subdomain of \(\Omega \) for all eigenfunctions \(\phi _\sigma \) with \(\sigma \) large enough. To state these results we first introduce certain conventions and notations, and we review known facts and results from complex analysis.

In what follows, and throughout the remainder of this article, \({\mathbb {R}}^2\) is identified with the complex plane \({\mathbb {C}}\), \(\Omega \subset {\mathbb {C}}\) denotes a bounded, simply-connected open set with analytic boundary, and \(D:=\{z\in {\mathbb {C}}\mid |z|<1\}\) denotes the open unit disc in the complex plane. Under these assumptions it follows from the Riemann mapping theorem [2] that there is a smooth map \(f:{\overline{D}}\rightarrow {\mathbb {C}}\) such that \(f|_{D}:D\rightarrow \Omega \) is a biholomorphism and \(|\partial _zf|>0\) on \({\overline{D}}\)—that is to say, \(f|_D:D\rightarrow \Omega \) is a biholomorphic conformal mapping of \(\Omega \) up to and including \(\partial \Omega \). We call such a function f a mapping function for\(\Omega \). Note that, denoting by \(\partial _r\) and \(\partial _\nu \) the radial derivative on the boundary of D and the normal derivative on the boundary of \(\Omega \), respectively, we have \(\partial _r=|\partial _z f|\partial _\nu \) and \(|\partial _zf|>0\). Thus, for \(z\in \partial D\) the function

$$\begin{aligned} u_{\sigma _j}:=\phi _{\sigma _j}\circ f \end{aligned}$$
(1.3)

satisfies,

$$\begin{aligned} \partial _r u_{\sigma _j}(z)=|\partial _zf(z)|\partial _\nu \phi _{\sigma _j}(f(z))=|\partial _zf(z)|\sigma _j \phi _{\sigma _j}(f(z)), \end{aligned}$$

and, hence, the generalized Steklov eigenvalue problem

$$\begin{aligned} {\left\{ \begin{array}{ll} -\Delta u_{\sigma _j}=0&{}\hbox {in }D\\ \partial _r u_{\sigma _j}=\sigma _j |\partial _z f|u_{\sigma _j}&{}\hbox {on }\partial D. \end{array}\right. } \end{aligned}$$
(1.4)

Finally we introduce notation for the relevant Fourier analysis. For \(v\in C({\overline{D}})\) we let

$$\begin{aligned} {\hat{v}}(k)=\frac{1}{2\pi }\int _0^{2\pi } v(\cos \theta ,\sin \theta )e^{-ik\theta }d\theta \end{aligned}$$
(1.5)

denote the “boundary Fourier coefficients”, namely, the Fourier coefficients of the restriction \(v|_{_{\partial D}}\) of v to \(\partial D\). Where notationally useful, we write \({\mathcal {F}}[v] = {\hat{v}}\).

Definition 1.1

We say that the Steklov problem on \(\Omega \) satisfies the tunneling condition if there is \(m_0>0\) and a mapping function f for \(\Omega \), such that for all \(K>0\) there is \(C_0>0\) satisfying for any m

$$\begin{aligned} |{\hat{u}}_{\sigma }(k)|\le C_0^{|k-m|}\left( \sum _{\ell =m-m_0}^{m+m_0}|{\hat{u}}_{\sigma }(\ell )|^2\right) ^{\frac{1}{2}},\qquad |k|\le K\sigma . \end{aligned}$$

Lemma 4.1 shows that any tunneling Steklov problem there exist \(\sigma _0>0\) so that for each \(m\in {\mathbf {Z}}\) there is a constant \(C>0\) such that for \(\sigma >\sigma _0\),

$$\begin{aligned} e^{-C\sigma }\Vert {\hat{u}}_{{{\sigma }}}\Vert _{\ell ^2}\le \left( \sum _{k= m-m_0}^{{{m+}}m_0}|{\hat{u}}_{{{\sigma }}}(k)|^2\right) ^{\frac{1}{2}}. \end{aligned}$$
(1.6)

This estimate and its connections with similar results in quantum mechanics motivate the “tunneling” terminology introduced in Definition 1.1. To explain this, we first review the notion of tunneling from quantum mechanics. In that setting, we consider quantum wave-functions \(\psi _h\) which are solutions to

$$\begin{aligned} (-h^2\Delta +V-E)\psi _h=0,\qquad \Vert \psi _h\Vert _{L^2({\mathbb {R}}^n)}=1 \end{aligned}$$

where \(V\in C^\infty ({\mathbb {R}}^n;{\mathbb {R}})\) satisfies \(V(x)\underset{|x|\rightarrow \infty }{\longrightarrow } \infty \) and \(E\in {\mathbb {R}}\) is the energy of \(\psi _h\). Since \(\psi _h\) is the wave-function of a quantum particle with energy E under the influence of the potential V, the corresponding classical Hamilatonian is \(p=|\xi |^2+V-E\) and the corresponding classical particle lies on the surface \(\{p=0\}\). The classical tunneling result from quantum mechanics then states that, despite the fact that \(\psi _h\) is classically localized to \(\{p=0\}\), for anyU open and bounded, there is \(c>0\) such that

$$\begin{aligned} \Vert \psi _h\Vert _{L^2(U)}\ge ce^{-c/h}. \end{aligned}$$

That is, even for \(U\subset \{V(x)>E\}\), the classically forbidden region, there is some positive (albeit exponentially small in h) probability of finding the particle U. (We refer the reader to e.g. [21, Chapters 4, 7] for a description of the phase space perspective for studying Schrödinger operators and this type of tunneling estimate.)

Returning to the setting of Steklov eigenfunctions, recall that \(u_{{{\sigma }}}\) is an eigenfunction of the Dirichlet to Neumann map which is a pseudodifferential operator on \(\partial \Omega \) with symbol \(|\xi |_g\) where g is the metric on \(\partial \Omega \) [16, Sec. 7.11, Vol 2]. Therefore, the classical problem corresponding to the Steklov problem is the Hamiltonian flow for the Hamiltonian \(|\xi |_g\) on \(T^*\partial \Omega \) at energy \(\sigma ^{-1}|\xi |_g=1\)—which describes the motion of a free particle on \(\partial \Omega \). The allowable energies for this classical problem are given by \(\{\sigma ^{-1}|\xi |_g=1\}\) which, in the Fourier series representation correspond to \(\sigma = |\xi |_g\sim |k|\). Thus, the classically forbidden region is \(\big |\sigma ^{-1}|k|-1\big |>c>0\). Equation (1.6) tells us that, in cases for which the Steklov problem on \(\Omega \) is tunneling, Steklov eigenfunctions carry positive energy even in the classically forbidden region \(\sigma ^{-1}|k|\ll 1\), with an energy value that is no smaller than exponentially decaying in \(\sigma \). (Using the estimates of [5] one can also see that Steklov eigenfunctions carry at most exponentially small energy in the forbidden region.)

Theorem 2

Assume that the Steklov problem on \(\Omega \) is tunneling and let \(\sigma \) denote a Steklov eigenvalue for the set \(\Omega \). Let

$$\begin{aligned} {\tilde{u}}_{\sigma ,\delta }:=u_{\sigma }|_{B(0,\delta )} = \sum _{k=-\infty }^\infty {\hat{u}}_{\sigma }(k) r^{|k|}e^{ik\theta },\qquad {\tilde{u}}_{\sigma ,\delta ,m}:=\sum _{|k|< m} {\hat{u}}_{\sigma }(k) r^{|k|}e^{ik\theta }.\quad \end{aligned}$$
(1.7)

Then, there exist a constant \(c>0\) such that, for each integer \(N>0\), there are constants \(C_N\), \(\sigma _0\), \(\delta _0\), and \(m_0>0\) so that for all \(0<\delta <\delta _0\), \(m>m_0\), and \( \sigma >\sigma _0\) the inequality

$$\begin{aligned} \frac{\Vert {\tilde{u}}_{\sigma ,\delta } -{\tilde{u}}_{\sigma ,\delta ,m}\Vert _{C^N(B(0,\delta ))}}{\Vert {\tilde{u}}_{\sigma ,\delta }\Vert _{L^2(B(0,\delta ))}} \le C_N (\delta ^{m-N-m_0-1}+e^{-c\sigma }) \end{aligned}$$
(1.8)

holds.

Letting \(\{\phi _{\sigma _j}\}_{j=1}^\infty \) denote an orthonormal basis of Steklov eigenfunctions and calling \(u_{\sigma _j}=\phi _{\sigma _j}\circ f\), Theorem 2 shows in particular that

$$\begin{aligned} u_{\sigma _j}=\sum _{|k|<m}{\hat{u}}_{\sigma _j}(k)r^{|k|}e^{ik\sigma }+O\left( (r^{m-m_0-1}+e^{-c\sigma _j})\sqrt{\sum _{|k|<m} |{\hat{u}}_{\sigma _j}(k)|^2 \frac{r^{2k+1}}{2k+1}}\right) .\nonumber \\ \end{aligned}$$
(1.9)

In other words, for r small, \(u_{\sigma _j}\) is well approximated by a function with finitely many Fourier modes. If there is \(c>0\) such that

$$\begin{aligned} |{\hat{u}}_{\sigma _j}(0)|\ge c\sqrt{\sum _{0<|k|<m}|{\hat{u}}_{\sigma _j}(k)|^2}, \end{aligned}$$

then we obtain

$$\begin{aligned} u_{\sigma _j}={\hat{u}}_{\sigma _j}(0) +O((r+e^{-c\sigma _j})|{\hat{u}}_{\sigma _j}(0)|) \end{aligned}$$

and \(u_{\sigma _j}\) is nearly constant on small balls centered around 0. In general, however, finitely many Fourier modes are necessary to capture the lowest-order asymptotics, as indicated in Eq. (1.9).

One of the main components of the proof of Theorem 1, in addition to Theorem 2, is the construction of a large class of domains \(\Omega \) for which the Steklov problem is tunneling. To this end, we introduce some additional definitions. A function \(v\in C(D)\) will be said to be boundary-band-limited provided \({\hat{v}}(k) = 0\) except for a finite number of values of \(k\in {\mathbb {Z}}\). We say that a mapping function f is boundary band limited conformal (BBLC) if \(|\partial _z f|\) is boundary band-limited. If in addition, \(|\partial _z f||_{\partial D}\) is non-constant, we will write that \(\Omega \) is BBLCN. Finally, we say the domain \(\Omega \) is BBLC (BBLCN) if and only if a BBLC (BBLCN) mapping function, \(f:D\rightarrow \Omega \) exists. We now present the main theorem of this paper.

Theorem 3

Assume \(\Omega \) is BBLCN. Then the Steklov problem on \(\Omega \) is tunneling.

Fig. 1
figure 1

Fixed-sign sets for Steklov eigenfunctions over the elliptical domain \(\Omega = x^2+\frac{y^2}{1.01^2}=1\). The yellow and blue regions indicate the subsets over which the eigenfunctions are positive and negative, respectively. The left and right images correspond to the eigenvalues \(\sigma _{20}=9.9502\) and \(\sigma _{30}=14.9253\), respectively. For a circle the nodal lines coincide with a set of j uniformly arranged radial lines from the center to the boundary: they are dense at scale \(\sigma _j^{-1}=j^{-1}\) over the complete domain, including the origin. Under the barely-visible perturbation of the unit disc into the slightly elliptical domain \(\Omega \), regions of asymptotically fixed size on which the eigenfunction does not change sign open-up within \(\Omega \). Indeed, the nodal set corresponding to \(\sigma _{30}\) (right image) shows such an opening, whereas the nodal set corresponding to \(\sigma _{20}\) (left image) does not; cf. also Remark 1.2

Remark 1.2

It is not clear whether the elliptical and kite-shaped domains (Eqs, (6.1) and (6.2)) considered in Figs. 14 and 5 satisfy the BBLCN condition or, more generally, whether they have tunneling Steklov problems (we have not as yet been able to establish that the tunneling condition holds for domains that are not BBLCN). However, domain-opening observations such as those displayed in Fig. 1 and Sect. 6, suggest that these domains may nevertheless be tunneling. This and other domain-opening observations provide support for Conjecture 1.3 below. (Steklov eigenfunctions on a domain which satisfies the BBLCN condition, and, therefore, in view of Theorem 3, is known to be tunneling, are displayed in Fig. 2).

In view of Remark 1.2 we conjecture that every Steklov problem on an analytic domain is tunneling unless the Steklov domain \(\Omega \) is a disc:

Conjecture 1.3

Let \(\Omega \subset {\mathbb {R}}^2\) be a bounded, simply-connected domain with real analytic boundary that is not equal to B(xr) for any \(x\in {\mathbb {R}}^2\), \(r>0\). Then the Steklov problem on \(\Omega \) is tunneling.

1.1 Outline of the Paper

This paper is organized as follows. Section 2 shows that arbitrary analytic, bounded, simply-connected domains can be approximated arbitrarily closely by BBLCN domains (Fig. 3). Then, Sects. 3 and 4 provide proofs for Theorems 3 and 2 , respectively. The numerical methods used in this paper to produce accurate Steklov eigenvalues, eigenfunctions, and associated nodal sets are presented in Sect. 5. Section 6, finally, illustrates the methods with numerical results for elliptical and kite-shaped domains.

Remark 1.4

Throughout this article we abuse notation slightly by allowing C to denote a positive constant that may change from line to line but does not depend on any of the parameters in the problem. In addition \(C_N\) is a positive constant that may change from line to line and depends only on the parameter N.

Fig. 2
figure 2

Steklov eigenfunctions on the domain \(\Omega \) whose mapping function, which is given by Eq. (2.3), maps the center of the disk to the point \(z_0 = (0.8,0)\) (marked by red asterisks in the figures). The corresponding Steklov eigenvalues are given by \(\sigma _{16}=7.9642\) (top left), \(\sigma _{40}=19.8173\) (top right), and \(\sigma _{60}=29.8197\) (bottom left). Note that, according to Corollary 2.2 the set \(\Omega \) is a BBLCN approximation to the disk. As predicted by Theorem 2, oscillations avoid a region around \(z_0\) for high \(\sigma \). The bottom-right image displays a typical eigenfunction on the exact disc. Note the dramatic change that arises in the Steklov eigenfunctions from a barely visible boundary perturbation of the disc

Fig. 3
figure 3

The function \(\lambda \) for an ellipse (left) and a kite-shaped domain (right)

2 Approximation by Tunneling Domains

This section shows that any analytic domain can be approximated arbitrarily closely (in a sense made precise in Corollary 2.2) by a BBLCN domain. To do this, first let \(M\ge 0\), \(\alpha _i\in {\mathbb {C}}\setminus {\overline{D}}\) for \(i=1,\ldots ,N\), and let \(N_i\ge 1\), \(i=1,\ldots M\), and let us seek approximating BBLCN domains whose mappings \(f:{\overline{D}}\rightarrow {\mathbb {C}}\) take the form

$$\begin{aligned} f(z)=\int _0^z p^2(w)dw,\qquad p(z)= \prod _{i=1}^M (z-\alpha _i)^{N_i}. \end{aligned}$$

In words: f is the integral of the square of a polynomial with roots outside \({\overline{D}}\). It follows that

$$\begin{aligned} \partial _z f=\prod _{i=1}^M(z-\alpha _i)^{2N_i},\qquad |\partial _zf |=\prod _{i=1}^M (|z-\alpha _i|^{2})^{N_{{{i}}}}. \end{aligned}$$

In particular,

$$\begin{aligned} |\partial _z f|(e^{i\theta })=\prod _{i=1}^M(1-e^{i\theta }\overline{\alpha _i}-e^{-i\theta }\alpha _i+|\alpha _i|^2)^{N_i} \end{aligned}$$

which manifestly shows that \({|\partial _z f|}\) is boundary-band-limited.

We next show that an arbitrary non-vanishing analytic function on \({\overline{D}}\) can be approximated by the square of a polynomial.

Lemma 2.1

Let \(g:{\overline{D}}\rightarrow {\mathbb {C}}\) smooth with \(g|_{D}\) analytic and \(|g|>0\) on \({\overline{D}}\). Then, for any \(\varepsilon _0>0\) and \(k>0\), there are \(M>0\), \(\alpha _0\), \(\{(\alpha _i,N_i)\}_{i=1}^M\) with \(|\alpha _i|>1\), \(i=1,\ldots ,M\) such that

$$\begin{aligned} \Vert g- \alpha _0\prod _{i=1}^M (z-\alpha _i)^{2N_i}\Vert _{C^k({\overline{D}})}<\varepsilon _0. \end{aligned}$$

Proof

Define \(h:{\overline{D}}\rightarrow {\mathbb {C}}\) by

$$\begin{aligned} h(z)=\int _0^z\frac{g'(w)}{g(w)}dw +\log (g(0)) \end{aligned}$$

Then, since \({{D}}\) is simply-connected and \(|g|>0\) on \({\overline{D}}\), h is analytic in D with smooth extension to \({\overline{D}}\). In addition,

$$\begin{aligned} w(z)=e^{\frac{1}{2}h(z)} \end{aligned}$$

is an analytic function on D such that \(w^2(z)=g(z)\) and w extends smoothly to \({\overline{D}}\). Then, for all \(\varepsilon >0\), there is a polynomial \(p_\varepsilon \) such that

$$\begin{aligned} \Vert w(z)-p_\varepsilon (z)\Vert _{C^{k}({\overline{D}})}<\varepsilon \min (\Vert w(z)\Vert _{C^{k}({\overline{D}})},1) \end{aligned}$$

In particular, since \(|g|>c>0\) on \({\overline{D}}\), for \(0<\varepsilon \) small enough, \(p_\varepsilon \) has no zeros in \({\overline{D}}\). Hence,

$$\begin{aligned} p_\varepsilon =\beta _0\prod _{i=1}^M(z-\beta _i)^{N_i} \end{aligned}$$

for some \(|\beta _0|>0\), \(|\beta _i|>1\), \(i=1,\ldots , M\). Multiplying by \(w+p_\varepsilon \), we have

$$\begin{aligned} \Vert g(z)-p^2_\varepsilon (z)\Vert _{C^k({\overline{D}})}&=\Vert (w-p_\varepsilon )(w+p_\varepsilon )\Vert _{C^k({\overline{D}})}\\&\le C_k\Vert (w-p_\varepsilon )\Vert _{C^{k}({\overline{D}})}\Vert (w+p_\varepsilon )\Vert _{C^{k}({\overline{D}})}\\&\le C_k\varepsilon (2+\varepsilon )\Vert w\Vert _{C^{k}({\overline{D}})} \end{aligned}$$

Choosing \(\varepsilon =\frac{\varepsilon _0}{C_k}\min (\frac{1}{3\Vert w\Vert _{C^{k}({\overline{D}})}},1)\) proves the result with \(\alpha _0=\beta _0^2\) and \(\alpha _i=\beta _i\). \(\square \)

This result can be used to approximate any analytic domain by a BBLCN domain:

Corollary 2.2

For any analytic, bounded, simply-connected domain \(\Omega \), \(k>0\), and \(\varepsilon _0>0\) there is a BBLCN domain \(\Omega _{\varepsilon _0}\) and \(g_{\varepsilon _0}\in C^\infty (\partial \Omega )\) such that with \(\nu \) the outward unit normal to \(\Omega \),

$$\begin{aligned} \partial \Omega _{\varepsilon _0}=\{ x+\nu g_{\varepsilon _0}(x)\mid x\in \partial \Omega \},\qquad \Vert g_{\varepsilon _0}\Vert _{C^k(\partial \Omega )}<\varepsilon _0. \end{aligned}$$
(2.1)

Proof

Since \(\Omega \) is analytic, there is \(f:{\overline{D}}\rightarrow {\mathbb {C}}\) analytic such that \(f|_{D}:D\rightarrow \Omega \) is a biholomorphism and \(|\partial _z f|>0\) on D. Moreover, by [2], \(\partial _zf\) has a smooth extension to \({\overline{D}}\). Then, applying Lemma 2.1 with \(g=\partial _zf(z)\) gives

$$\begin{aligned} p_\varepsilon =\alpha _0\prod _{i=1}^M(z-\alpha _i)^{N_i} \end{aligned}$$

a polynomial with no roots in \({\overline{D}}\) such that

$$\begin{aligned} \Vert \partial _zf(z)-p^2_\varepsilon (z)\Vert _{C^{\max (k,1)}({\overline{D}})}<\varepsilon . \end{aligned}$$

Note also that adjusting p if necessary we may assume that the restriction of \(|p_\varepsilon |\) to \(\partial D\) is not constant. Then, defining

$$\begin{aligned} f_\varepsilon :=\int _0^z p^2_\varepsilon (w)dw+f(0) \end{aligned}$$
(2.2)

we have

$$\begin{aligned} \Vert f_\varepsilon -f\Vert _{C^{\max (k+1,2)}({\overline{D}})}<\varepsilon ,\qquad \partial _zf_\varepsilon =p^2_\varepsilon , \end{aligned}$$

so that \(\big |\partial _z f_\varepsilon \big ||_{_{\partial D}}\) is non-constant and band limited. Moreover, since f is a biholomorphism, for \(\varepsilon >0\) small enough, \(f_\varepsilon \) is also a biholomorphism. We next show that since \(\Vert f_\varepsilon -f\Vert _{C^{\max (k+1,2)}({\overline{D}})}<\varepsilon \), for \(\varepsilon >0\) small enough the curve

$$\begin{aligned} \partial \Omega _\varepsilon =\{f_\varepsilon (z)\mid |z|=1\} \end{aligned}$$

can be expressed in the form (2.1). To do this let

$$\begin{aligned} F(t,\theta ,\omega ,s)=f(e^{i\theta })-tf_\varepsilon (e^{i(\omega +\theta )})-(1-t)f(e^{i(\omega +\theta )})-sf'(e^{i\theta })e^{i\theta } \end{aligned}$$

and note that \(F(1,\theta ,\omega ,s)=0\) if and only if

$$\begin{aligned} f_\varepsilon (e^{i(\omega +\theta )})=f(e^{i\theta })\pm s\nu (\theta ). \end{aligned}$$

Therefore, we aim to find \(s=s(\theta )\) and \(\omega =\omega (\theta )\) such that \(F(1,\theta ,\omega (\theta ),s(\theta ))=0\). Note that

$$\begin{aligned} \begin{aligned} \partial _sF&=-f'(e^{i\theta })e^{i\theta }\\ \partial _{\omega }F&= -ie^{i(\omega +\theta )}\Big (f'(e^{i(\omega +\theta )})+t\big (f_\varepsilon '(e^{i(\omega +\theta )})-f'(e^{i(\omega +\theta )})\big ){{\Big )}} \end{aligned} \end{aligned}$$

In particular,

$$\begin{aligned} \partial _{\omega }F=i\partial _sF+O(\varepsilon )+O(|\omega |). \end{aligned}$$

Writing \({\tilde{F}}=\begin{pmatrix}{{\text {Re}}\,}F\\ {{\text {Im}}\,}F\end{pmatrix}\), we have

$$\begin{aligned} \partial _{s,\omega }{\tilde{F}}=\begin{pmatrix} -{{\text {Re}}\,}(f'(e^{i\theta })e^{i\theta })&{}{{\text {Im}}\,}(f'(e^{i\theta })e^{i\theta })\\ -{{\text {Im}}\,}(f'(e^{i\theta })e^{i\theta })&{}-{{\text {Re}}\,}(f'(e^{i\theta })e^{i\theta })\end{pmatrix}+O(|\omega |+\varepsilon ) \end{aligned}$$

Therefore, there are \(C>0\), \(\delta >0\), \(\varepsilon _0>0\) such that for \(0<\varepsilon <\varepsilon _0\), \(|\omega _0|<\delta \), \(t_0\in (-1,2)\), and \(s_0\in [-1,1]\)

$$\begin{aligned} {{\big \Vert \big (\partial _{s,\omega }{\tilde{F}}\big )^{-1}\big \Vert \le C.}} \end{aligned}$$

In particular, by the implicit function theorem, there is \(c>0\) such that if \(|\omega _0|<\delta \), \(t_0\in (-1,2)\), and \(F(t_0,\theta _0,\omega _0,s_0)=0\), then for \(|t-t_0|<c\), \(|\theta -\theta _0|<c\), \(|\omega -\omega _0|<c\) and \(|s-s_0|<c\), \(\omega =\omega (t, \theta )\) and \(s={s}(t,{{\theta }})\) are the unique solutions of \(F(t,\theta ,\omega ,s)=0\). In particular, since \(F(0,\theta ,0,0)=0\), the solutions \(s=s(t,\theta )\) and \(\omega =\omega (t,\theta )\) can be continued as functions of t as long as \(|\omega (t,\theta )|\) remains small.

We next note that

$$\begin{aligned} \begin{pmatrix} \partial _t \omega \\ \partial _ts\end{pmatrix}=\begin{pmatrix}\partial _\omega {{{\tilde{F}}}}&\partial _s{{{\tilde{F}}}}\end{pmatrix}^{-1}\partial _t{{{\tilde{F}}}}=O(\Vert f_\varepsilon -f\Vert _{L^\infty })=O(\varepsilon ), \end{aligned}$$

and, therefore,

$$\begin{aligned} |\omega (t,\theta )|+|s(t,\theta )|\le \int _0^{t}|\partial _t\omega (r,\theta )|+|\partial _ts(r,\theta )|dr\le C t\varepsilon . \end{aligned}$$

Hence for \(\varepsilon \) small enough the solutions \(\omega (t,\theta )\) and \(s(t,\theta )\) continue to \(t=1\) and satisfy

$$\begin{aligned} |\omega (1,\theta )|+|s(1,\theta )|\le C\varepsilon . \end{aligned}$$

Again, using the implicit function theorem, this implies that \(\omega (\theta ):=\omega (1,\theta )\) and \(s(\theta ):=s(1,\theta )\) are \(2\pi \)-periodic. Differentiating k times now yields

$$\begin{aligned} |\partial _\theta ^ks|\le C_k\varepsilon , \end{aligned}$$

finishing the proof by setting \(g_{\varepsilon _0}=\pm s\) and shrinking \(\varepsilon >0\) as necessary. (Here the ± corresponds to whether \(f(e^{i\theta })\) is positively (−) or negatively (\(+\)) oriented.) \(\square \)

Remark 2.3

Since the map \(f_\varepsilon \) in Eq. (2.2) may send 0 to a point \(z_0\) close to the boundary, it is interesting to see how the Steklov eigenfunctions rearrange their nodal sets in such a way that Theorems 1 and 2 are satisfied on the image of \(f_\varepsilon \). To demonstrate this let \(|a|<1\), consider the biholomorphic function \(f(z):=\frac{z-a}{{\bar{a}}z-1}\), and let \(f_\varepsilon \) denote the approximant of f given by Eq. (2.2) with

$$\begin{aligned} p_\varepsilon ({{w}})=i \sqrt{1- |a|^2}\sum _{j=0}^N({\bar{a}}w)^j\quad \text{ with } N=20 \text{ and } a=0.8. \end{aligned}$$
(2.3)

(This polynomial was obtained as the N-th order Taylor polynomial of \(\sqrt{\partial _z f}\).) In this case, according to Theorems 1 and 2 , the Steklov eigenfunctions should be slowly oscillating in a \(\sigma \) independent neighborhood of \(z_0\). Figure 2 displays corresponding Steklov eignfunction or various orders as well as a typical eigenfunction for the exact disc. Note the dramatic change that arises in the Steklov eigenfunctions from a barely visible boundary perturbation of the disc.

3 BBLCN Domains and Tunneling Steklov Problems

This section presents a proof of Theorem 3. In preparation for that proof, let \(\Omega \subset {\mathbb {C}}\) be a BBLCN domain, and denote by f the corresponding mapping function. Define

$$\begin{aligned} {\mathcal {F}}\left[ \big |\partial _zf\big |\right] (n):=a_n,\quad n\in {\mathbb {Z}} \quad \left( a_0:=\frac{1}{2\pi }\int _0^{2\pi }|\partial _z f(e^{i\theta })|\,d\theta >0\right) . \end{aligned}$$

Since \(\Omega \) is a BBLCN domain, the function \(\big |\partial _z f\big ||_{_{\partial D}}\) is band limited and \(\big |\partial _z f\big ||_{_{\partial D}}\) is not identically constant. It follows that

$$\begin{aligned} m_0:=\sup \{|n|: |a_n|\ne 0\} \end{aligned}$$

satisfies \(1\le m_0<\infty \).

Denoting by \({\hat{u}}(n)\) the boundary Fourier coefficients of an eigenfunction u, the corresponding boundary Fourier coefficients of \(\partial _r u\) are given by \(|n|{\hat{u}}(n)\). Thus, a solution to (1.4) is uniquely determined as an \(\ell ^2\) solution to the equation

$$\begin{aligned} |n|{\hat{u}}(n)=\sigma {\mathcal {F}} \left[ u\big |\partial _z f\big |\right] (n)\quad n\in {\mathbb {Z}}. \end{aligned}$$
(3.1)

In what follows we may, and do, assume that solutions \({\hat{u}}\) have \(\ell ^2\)-norm equal to one.

3.1 Proof of Theorem 3

Since

$$\begin{aligned} {\mathcal {F}}\left[ u\big |\partial _z f\big |\right] {{(n)}}=\sum _{m} a_{m}{\hat{u}}(n-m), \end{aligned}$$

it follows that (3.1) can be re-expressed in the form

$$\begin{aligned} |n|{\hat{u}}(n)=\sum _m \sigma a_m{\hat{u}}(n-m). \end{aligned}$$
(3.2)

From (3.2) we obtain

$$\begin{aligned} a_{-m_0}{\hat{u}}(n+m_0)=\frac{|n|}{\sigma }{\hat{u}}(n)-\sum _{m\ne -m_0}a_m{\hat{u}}(n-m), \end{aligned}$$

and, then, for all \(|n|\le K\sigma \),

$$\begin{aligned} |{\hat{u}}(n+m_0)|&\le |a_{-m_0}|^{-1}\left( \frac{||n|-\sigma a_0|}{\sigma }|{\hat{u}}(n)|+\sum _{m\ne 0,-m_0}|a_m||{\hat{u}}(n-m)|\right) \\&\le |a_{-m_0}|^{-1}\left( \frac{||n|-\sigma a_0|}{\sigma }|{\hat{u}}(n)|+\sum _{\begin{array}{c} m=-m_0+1\\ m\ne 0 \end{array}}^{m_0}|a_m||{\hat{u}}(n-m)|\right) \\&\le |a_{-m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })\sum _{k=n-m_0}^{n+m_0-1}|{\hat{u}}(k)|. \end{aligned}$$

The second inequality follows from the fact that \(a_n\equiv 0\) for \(|n|\ge m_0\), while the third one results from the relation \(a_0>0\) and the positivity, \(\sigma > 0\), of all nontrivial eigenvalues \(\sigma \), which imply that

$$\begin{aligned} ||n|-\sigma a_0|\le \max (|n|,\sigma |a_0|){\le \sigma (\max (K,\Vert a_m\Vert _{\ell ^\infty }))}. \end{aligned}$$

Making an identical argument, but solving for \({\hat{u}}(n-m_0)\), and using that \(|a_{m_0}|=|a_{-m_0}|\ne 0\), we have for all \(|n|\le K\sigma \),

$$\begin{aligned} \begin{aligned} |{\hat{u}}(n+m_0)|&\le |a_{m_0}|^{-1}\max ({{K}},\Vert a_m\Vert _{\ell ^\infty })\sum _{k=n-m_0}^{n+m_0-1}|{\hat{u}}(k)|,\\ |{\hat{u}}(n-m_0)|&\le |a_{m_0}|^{-1}\max ({{K}},\Vert a_m\Vert _{\ell ^\infty })\sum _{k=n-m_0+1}^{n+m_0}\!\!\!\!\!\!|{\hat{u}}(k)|. \end{aligned} \end{aligned}$$
(3.3)

We now use Eq. (3.3) to prove the first half of our tunneling estimate.

Lemma 3.1

Let \(m\in {\mathbb {Z}}\), \(K>0\), and

$$\begin{aligned} A_m:=\left( \sum _{k=m-m_0}^{m+m_0}|{\hat{u}}(k)|^2\right) ^{\frac{1}{2}}. \end{aligned}$$

Then, there exists \(C_0>0\) so that for all \(\sigma >0\) and for \(-K\sigma \le n+m\le K\sigma \) we have

$$\begin{aligned} |{\hat{u}}(n+m)|\le C_0^{|n|}A_m. \end{aligned}$$
(3.4)

Proof

We will assume \(m\ge 0\) since the other case follows similarly. The cases of \(n=-m_0,\ldots ,m_0\) are clear if we take \(C_0\ge 1\). Suppose (3.4) holds for \(-m_0\le n\le \ell \) with \(m_0\le \ell \). Then, by (3.3),

$$\begin{aligned} |{\hat{u}}(m+\ell +1)|&\le |a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })\sum _{k=\ell -2m_0+1}^{\ell }|{\hat{u}}(k+m)|\\&\le |a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })\sum _{k=\ell -2m_0+1}^{\ell }C_0^{|k|}A \end{aligned}$$

Now, if \(m_0\le \ell <2m_0\), then

$$\begin{aligned} |{\hat{u}}(m+\ell +1)|&\le |a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })\left( \sum _{k=0}^{\ell }C_0^{k}+\sum _{k=1}^{2m_0-\ell -1}C_0^k\right) A\\&\le |a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })\Big (\frac{C_0^{\ell +1}-1+C_0^{2m_0-\ell +1}-C_0}{C_0-1}\Big )A \end{aligned}$$

In particular, taking

$$\begin{aligned} C_0\ge 2|a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })+1 \end{aligned}$$

we have

$$\begin{aligned} |{\hat{u}}(m+\ell +1)|\le C_0^{\ell +1}A. \end{aligned}$$

Next, if \(2m_0\le \ell \), then

$$\begin{aligned} |{\hat{u}}(\ell +m+1)|&\le |a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })A\frac{ C_0^{\ell +1}-C_0^{\ell -2m_0+1}}{C_0-1} \end{aligned}$$

Taking \({C_0\ge 2 |a_{m_0}|^{-1}\max (K,\Vert a_m\Vert _{\ell ^\infty })+1}\) completes the proof for \(-m_0\le n\le K\sigma -m\).

An almost identical argument gives the \(-K\sigma -m \le n\le 0\) case. \(\square \)

4 Analysis of Tunneling Steklov Problems

The proof of Theorem 2 now follows in two steps. First, we show that, for eigenfunctions of any tunneling Steklov problem, the boundary Fourier coefficients of low frequency contain a mass no smaller than exponential in \(\sigma \). To finish the proof, we use the fact that the harmonic extension of \(e^{in\theta }\) decays exactly as \(r^{|n|}\). Examining the solution on the ball of radius \(\delta >0\) for some \(\delta \) small enough, it will be shown that the low frequencies dominate the behavior of u.

Lemma 4.1

Suppose that \(\Omega \) has tunneling Steklov problem. Then there exist \(\sigma _0>0\) so that for all \(m>0\) there is \(C>0\) such that for \(\sigma >\sigma _0\),

$$\begin{aligned} e^{-C\sigma }\Vert {\hat{u}}\Vert _{\ell ^2}\le \left( \sum _{k= m-m_0}^{{{m+}}m_0}|{\hat{u}}(k)|^2\right) ^{\frac{1}{2}}=:A_m. \end{aligned}$$

Proof

First, note that by e.g. [5, Corollary 1.3], for \(\sigma >3m\) there is \(C>0\) so that

$$\begin{aligned} \sum _{|k-m|\le 2\sigma }|{\hat{u}}(k)|^2\ge \Vert {\hat{u}}\Vert _{\ell ^2}^2(1-Ce^{-\sigma /C}). \end{aligned}$$

By Lemma 3.1

$$\begin{aligned} \sum _{|k-m|\le 2\sigma }|{\hat{u}}(k)|^2\le \sum _{0\le k-m\le 2\sigma } C_0^{2k} A_m^2+\sum _{-2\sigma \le k-m<0}C_0^{2|k|}A_m^2\le 2 \frac{2C_0^{4\sigma +2}-1}{C^2_0-1}A_m^2 \end{aligned}$$

In particular,

$$\begin{aligned} \frac{C^2_0-1}{2(2C_0^{4\sigma +2}-1)}\Vert {\hat{u}}\Vert _{\ell ^2}^2(1-Ce^{-C\sigma })\le A_m^2=\sum _{k=-m_0}^{m_0}|{\hat{u}}(k)|^2. \end{aligned}$$

Taking \(\sigma _0\) large enough so that \(Ce^{-C\sigma }\le \frac{1}{2}\), finishes the proof. \(\square \)

Proof of Theorem 2

In what follows we utilize the definitions (1.7) for a given eigenvalue \(\sigma _j = \sigma \), and, for that eigenvalue we denote \({\hat{u}}(k) = {\hat{u}}_{\sigma _j}(k) = {\hat{u}}_{\sigma }(k)\). Then, applying the relation

$$\begin{aligned} \int _{B(0,\delta )} \left| \sum _k b_k r^{|k|}e^{ik\theta }\right| ^2= \sum _k |b_k|^2\frac{2\pi \delta ^{2|k|+{2}}}{2|k|+{2}}, \end{aligned}$$
(4.1)

which is valid for all sequences \(\{b_k\}_{k\in {\mathbb {Z}}}\subset {\mathbb {C}}\), to the right-hand equation in (1.7), for \(m\ge m_0\) we obtain

$$\begin{aligned} \Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert ^2_{L^2}= & {} \sum _{|k|\le m}\frac{2\pi \delta ^{2|k|+{2}}}{2{{|}}k{{|}}+{2}}|{\hat{u}}(k)|^2\nonumber \\\ge & {} 2\pi \frac{\delta ^{2m_0+{2}}}{2m_0+{2}}\sum _{|k|\le m_0}|{\hat{u}}(k)|^2=2\pi \frac{\delta ^{2m_0+{2}}}{2m_0+{2}}A^2. \end{aligned}$$
(4.2)

To estimate the error in approximating \(u_{\sigma ,\delta }\) by \({\tilde{u}}_{\sigma ,\delta ,m}\), first note that

$$\begin{aligned} \Big \Vert \sum _{|k|\ge 2\sigma } {\hat{u}}(k)r^{|k|}e^{ik\theta }\Big \Vert _{C^N(B(0,\delta ))}&\le \sum _{|k|\ge 2\sigma } |{\hat{u}}(k)|\cdot \Vert r^{|k|}e^{ik\theta }\Vert _{C^N(B(0,\delta ))}\\&\le \Big (\sum _{|k|\ge 2\sigma } |{\hat{u}}(k)|^2\Big )^{\frac{1}{2}}\Big (\sum _{|k|\ge 2\sigma } \Vert r^{|k|}e^{ik\theta }\Vert ^2_{C^N(B(0,\delta ))}\Big )^{\frac{1}{2}}\\&\le \Big (\sum _{|k|\ge 2\sigma } |{\hat{u}}(k)|^2\Big )^{\frac{1}{2}}\Big (\sum _{|k|\ge 2\sigma } k^{2N}\delta ^{2k-2N}\Big )^{\frac{1}{2}}\\&\le C_N\Vert {\hat{u}}\Vert _{\ell ^2}\delta ^{-N}{\sigma ^N}\delta ^{2\sigma }. \end{aligned}$$

Applying Lemma 4.1 with \(m=0\), and absorbing the \(\sigma ^N\) into the exponential factor we then obtain

$$\begin{aligned} \Big \Vert \sum _{|k|\ge 2\sigma } {\hat{u}}(k)r^{|k|}e^{ik\theta }\Big \Vert _{C^N(B(0,\delta ))}\le C_N\delta ^{-N}\delta ^{2\sigma }e^{C\sigma }A \end{aligned}$$

where \(A=A_0\) is given by

$$\begin{aligned} A:=\Big (\sum _{k=-m_0}^{m_0}|{\hat{u}}(k)|^2\Big )^{\frac{1}{2}}, \end{aligned}$$

We can now estimate

$$\begin{aligned}&\Big \Vert \sum _{|k|\ge m} {\hat{u}}(k) r^{|k|}e^{ik\theta }\Big \Vert _{C^N(B(0,\delta ))}\\&\quad \le \sum _{m\le |k|< 2\sigma } \Vert {\hat{u}}(k) r^{|k|}e^{ik\theta }\Vert _{C^N(B(0,\delta ))} +\Big \Vert \sum _{|k|\ge 2\sigma } {\hat{u}}(k) r^{|k|}e^{ik\theta }\Vert _{C^N(B(0,\delta ))}\\&\quad \le \sum _{m\le |k|< 2\sigma } |{\hat{u}}(k)|\cdot \Big \Vert r^{|k|}e^{ik\theta }\Big \Vert _{C^N(B(0,\delta ))} +C_N\delta ^{-N}\delta ^{2\sigma }e^{C\sigma }A\\ \end{aligned}$$

Thus, using the definition of tunneling (Definition 1.1), we obtain

$$\begin{aligned}&\Big \Vert \sum _{|k|\ge m} {\hat{u}}(k) r^{|k|}e^{ik\theta }\Big \Vert _{C^N(B(0,\delta ))}\\&\quad \le C_N\delta ^{m-N}A\sum _{m\le |k|< 2\sigma } C_0^{|k|}|k|^N\delta ^{|k|-m}+C_N\delta ^{-N}\delta ^{2\sigma }e^{C\sigma }A\\&\quad \le C_N\delta ^{m-N} A+C_N\delta ^{-N}\delta ^{2\sigma }e^{C\sigma }A \end{aligned}$$

provided that \(\delta <\frac{1}{2}C_0^{-1}\). Therefore, using (4.2),

$$\begin{aligned} \frac{\Vert {\tilde{u}}_{\sigma ,\delta } -{\tilde{u}}_{\sigma ,\delta ,m}\Vert _{C^N(B(0,\delta ))}}{\Vert {\tilde{u}}_{\sigma ,\delta }\Vert _{L^2(B(0,\delta ))}} \le C_N\delta ^{m-N-m_0-1}+C_N\delta ^{2\sigma -N-m_0-1}e^{C\sigma }. \end{aligned}$$

Thus, choosing \(\delta >0\) such that \(\delta < e^{-2C}\) and taking \(\sigma _0>N{+m_0+1}\) the claim follows. \(\square \)

We can now present a proof of Theorem 1.

Proof of Theorem 1

From Corollary 2.2 we know that there exists a tunneling domain \(\Omega _1\subset {\mathbb {C}}\) satisfying (1.2) for the given value \(\varepsilon >0\). Let \(\sigma _0\) be as in Theorem 2. Clearly, it suffices to prove the statement of the theorem for \(\sigma > \sigma _0\), since for \(\sigma \le \sigma _0\) the statement follows from the fact that there are finitely many Steklov eigenvalues below \(\sigma _0\) and that \({{\phi }}_{\sigma }\) cannot vanish in any open set. Therefore, we may and do assume \(\sigma > \sigma _0\) along with the other assumptions in Theorem 2, so that, in particular, inequality (1.8) holds. In what follows we write

$$\begin{aligned} L^2(B(0,\delta )) = L^2_\delta \quad \text{ and }\quad L^\infty (B(0,\delta )) = L^\infty _\delta \end{aligned}$$
(4.3)

Fixing \(m\ge m_0+2\), and letting \({\tilde{u}}_{\sigma ,\delta }\) and \({\tilde{u}}_{\sigma ,\delta ,m}\) be given by (1.7) (with \(u_{\sigma }\) related to \(\phi _\sigma \) via  (1.3)) we note that

$$\begin{aligned} \Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^\infty _\delta }\ge \frac{\Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta }}{\sqrt{\pi } \delta }. \end{aligned}$$

It follows that there exists \(x_0\in B(0,\delta )\) such that

$$\begin{aligned} |{\tilde{u}}_{\sigma ,\delta ,m}(x_0)|\ge \frac{\Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta }}{\sqrt{\pi } \delta }. \end{aligned}$$
(4.4)

Now, since \(\Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{C^1}\le C_{m,\delta }\Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2}\), it follows from (4.4) that there is \(r_{m,\delta }\in {\mathbb {R}}\), \(0< r_{m,\delta }<\delta \) (in particular, independent of \(\sigma \)) such that

$$\begin{aligned} |{\tilde{u}}_{\sigma ,\delta ,m}(x)|>\frac{\Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta }}{2\sqrt{\pi } \delta },\qquad x\in B(x_0,r_{m,\delta }). \end{aligned}$$
(4.5)

But, since \(m\ge m_0+2\), the estimate (1.8) with \(N=0\) yields

$$\begin{aligned} |{\tilde{u}}_{\sigma ,\delta ,m}(x)|\le & {} |{\tilde{u}}_{\sigma ,\delta }(x)| + |{\tilde{u}}_{\sigma ,\delta }(x)-{\tilde{u}}_{\sigma ,\delta ,m}(x)|\nonumber \\\le & {} C_0(\delta +e^{-c\sigma }) \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta } + |{\tilde{u}}_{\sigma ,\delta }(x)| \end{aligned}$$
(4.6)

and

$$\begin{aligned} \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta }\le & {} \Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta } + \Vert {\tilde{u}}_{\sigma ,\delta } - {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta } \nonumber \\\le & {} \Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta } + \sqrt{\pi }\delta C_0(\delta +e^{-c\sigma }) \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta }. \end{aligned}$$
(4.7)

(To establish the rightmost inequality in (4.7) the relation \(\Vert {\tilde{u}}_{\sigma ,\delta } - {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta }\le \sqrt{\pi }\delta \Vert {\tilde{u}}_{\sigma ,\delta } - {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^\infty _\delta }\) was used before the inequality (1.8) was applied.) From (4.7) we obtain

$$\begin{aligned} \Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta }\ge \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta } - \sqrt{\pi }\delta C_0(\delta +e^{-c\sigma }) \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta }. \end{aligned}$$
(4.8)

It follows from (4.5), (4.6) and (4.8) that

$$\begin{aligned}&C_0(\delta +e^{-c\sigma }) \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta } + |{\tilde{u}}_{\sigma ,\delta }(x)|\nonumber \\&\quad >\frac{\Vert {\tilde{u}}_{\sigma ,\delta ,m}\Vert _{L^2_\delta }}{2\sqrt{\pi } \delta } \nonumber \\&\quad \ge \frac{\Vert {\tilde{u}}_{\sigma ,\delta }\Vert _{L^2_\delta }}{2\sqrt{\pi } \delta }- \frac{C_0}{2}(\delta +e^{-c\sigma }) \Vert {\tilde{u}}_{\sigma ,\delta } \Vert _{L^2_\delta },\quad x\in B(x_0,r_{m,\delta }), \end{aligned}$$
(4.9)

and, therefore

$$\begin{aligned} |{\tilde{u}}_{\sigma ,\delta }(x)|> \Vert {\tilde{u}}_{\sigma ,\delta }\Vert _{L^2_\delta }\left( \frac{1}{2\sqrt{\pi } \delta }- \frac{3C_0}{2}(\delta +e^{-c\sigma }) \right) , \qquad x\in B(x_0,r_{m,\delta }).\qquad \end{aligned}$$
(4.10)

Taking \(\delta _1\) sufficiently small and \(\delta \le \delta _1\) the inequality

$$\begin{aligned} \frac{3C_0}{2}(\delta +e^{-c\sigma _0})<\frac{1}{2\sqrt{\pi }\delta } \end{aligned}$$

holds, and it therefore follows that for a certain constant \(D>0\) we have

$$\begin{aligned} |{{{\tilde{u}}}}_{\sigma ,\delta }(x)|>\frac{D\Vert {\tilde{u}}_{\sigma ,\delta }\Vert _{L^2_\delta }}{\delta }\quad \text{ for }\quad x\in B(x_0,r_{m,\delta }) \end{aligned}$$

provided \(\delta <\delta _1\). In particular,

$$\begin{aligned} |\phi _{\sigma }(x)|>0,\qquad x\in f(B(x_0,r_{m,\delta })). \end{aligned}$$

Since the derivative of f never vanishes, for \(\delta <\delta _1\) and for a certain \(E>0\) there is a ball \(\mathcal {B}\) of radius \(Er_{m,\delta }\) such that \(\phi _{\sigma }\) does not vanish on \(\mathcal {B}\). The proof is now complete.

\(\square \)

5 Numerical Formulation

5.1 Integral Representation

Let \(\Omega \subset {\mathbb {R}}^2\) denote a domain with, say, a \(C^2\) boundary, and let

$$\begin{aligned} S[\phi ](x) := \int \limits _{\partial \Omega } G(x,y) \phi (y) d{{s}}(y),\quad x\in {\mathbb {R}}^2 ,\quad G(x,y) = -\frac{1}{2\pi }\log |x-y|, \end{aligned}$$

denote the Single Layer Potential (SLP) for a given density \(\phi :\partial \Omega \rightarrow {\mathbb {R}}\) in a certain Banach space H of functions. Both Sobolev and continuous spaces H of functions lead to well developed Fredholm theories in this context [10, 12]. It is useful to recall that, as shown e.g. in the aforementioned references, the limiting values of the potential S and its normal derivative on \(\partial \Omega \) can be expressed in terms of well known “jump conditions” that involve the single and double layer boundary integral operators

$$\begin{aligned} {\mathcal {S}}[\phi ](x):= & {} \int _{\partial \Omega } G(x,y)\phi (y) ds(y) \quad \quad \text{ and }\quad {\mathcal {T}}[\phi ](x)\\:= & {} \int _{\partial \Omega } \frac{\partial G(x,y)}{\partial {\nu }(x)}\phi (y)ds(y),\quad x\in \partial \Omega , \end{aligned}$$

respectively.

In view of the jump conditions for the SLP [10], use of the representation

$$\begin{aligned} u(x) = S[\phi ](x),\quad x\in \Omega , \end{aligned}$$
(5.1)

for the eigenfunction u, the Steklov boundary condition in Eq. (1.1) gives rise to the generalized eigenvalue problem

$$\begin{aligned} \left( \frac{1}{2}{\mathcal {I}}+{\mathcal {T}}\right) \left[ \phi \right] = {\sigma }\, {\mathcal {S}}[\phi ] \quad \text{ for }\quad x \in \partial \Omega . \end{aligned}$$
(5.2)

Unfortunately, however, the single layer operator \({\mathcal {S}}\) on the right side of this equation is not always invertible. In order to avoid singular right-hand sides and the associated potential sensitivity to round-off errors, in what follows we utilize the Kress potential

$$\begin{aligned} u(x) = S_0[\phi ](x) = \int \limits _{\partial \Omega } G(x,y) \left( \phi (y)-{\overline{\phi }}\right) d{{s}}(y) + {\overline{\phi }},\quad x \in \Omega \end{aligned}$$
(5.3)

(where \({\overline{\phi }}\) denotes the average of \(\phi \) over \(\partial \Omega \)), which leads to the modified eigenvalue equation [1]

$$\begin{aligned} \left( \frac{1}{2}{\mathcal {I}}+{\mathcal {T}}\right) \left[ \phi - {\overline{\phi }}\right] = {\sigma }\left( {\mathcal {S}}[\phi -{\overline{\phi }}]+ {\overline{\phi }} \right) \quad \text{ for }\quad x \in \partial \Omega . \end{aligned}$$
(5.4)

The right-hand operator in this equation is invertible [10, Theorem 7.41], as desired. For either formulation, the evaluation of a given eigenfunction u requires evaluation of the SLP, in accordance with either (5.1) or (5.3), for the solution \(\phi \) of the corresponding generalized eigenvalue problem (5.2) or (5.4), respectively, at all required points \(x\in \Omega \).

Remark 5.1

Note that for a given harmonic function u in \(\Omega \), \(\phi \) in (5.2) and that in (5.4) are not the same.

5.2 Fourier Expansion and Exponential Decay

In terms of a given \(2\pi \)-periodic parametrization C(t) of \(\partial \Omega \), the Steklov eigenfunction u corresponding to a given solution \((\phi ,{\sigma })\) of the regularized eigenvalue problem (5.4), which is given by the single layer expression (5.3), can be expressed, for a given point \(x = (x_1,x_2)\in \Omega \),

$$\begin{aligned} u(x_1,x_2)= & {} {\overline{\phi }}{+}\frac{1}{{4}\pi } \int \limits _0^{2\pi } \log \left[ \left( x_1 - C_1(t) \right) ^2 + \left( x_2 -C_2(t) \right) ^2 \right] \left[ \phi \left( C(t)\right) - {\overline{\phi }} \right] \nonumber \\&\times \left| {\dot{C}}(t) \right| dt, \end{aligned}$$
(5.5)

where \(C(t)=(C_1(t),C_2(t))\) and where \({\overline{\phi }}\) denotes the average of \(\phi \) over the curve \(\partial \Omega \). Unfortunately, a direct use of this expression does not capture important elements in the eigenfunction within \(\Omega \), such as the nodal sets, since, for analytic domains, the eigenfunctions decay exponentially fast within \(\Omega \) as the frequency increases [5, 13]. In regions where the actual values of the eigenfunction may be significantly below machine precision the expression (5.5) must be inaccurate: this expression can only achieve the exponentially small values via the cancellations that occur as the solution \(\phi \) becomes more and more oscillatory. But such cancellations cannot take place numerically below the level of machine precision. In order to capture the decay explicitly within the numerical algorithm we proceed in a manner related to the construction used in [13].

To accurately obtain the exponentially decaying values of the Steklov eigenfunction we proceed as follows. We first consider the Fourier expansion

$$\begin{aligned} \left[ \phi \left( C(t)\right) - {\overline{\phi }}\, \right] \left| {\dot{C}}(t) \right| = \sum _{\begin{array}{c} n\in {\mathbb {Z}}\\ n\ne 0 \end{array}} A_n e^{int}. \end{aligned}$$
(5.6)

of the product \(\left[ \phi (C(t))-{\overline{\phi }}\,\right] \left| {\dot{C}}(t) \right| \); note that, as is easily checked, the \(n=0\) term in the Fourier expansion (5.6) is indeed equal to zero. Inserting this expansion in (5.5) we obtain

$$\begin{aligned} u(x_1,x_2) = {\overline{\phi }}+\sum _{\begin{array}{c} n\in {\mathbb {Z}}\\ n\ne 0 \end{array}} A_n B_n^0(x_1,x_2),\quad \text{ where } \\ B_n^0(x_1,x_2) = -\frac{1}{{4}\pi } \int \limits _0^{2\pi } \log \left[ \left( x_1 - C_1(t) \right) ^2 + \left( x_2 -C_2(t) \right) ^2 \right] e^{int}dt. \end{aligned}$$

Then, assuming an analytic boundary, as is relevant in the context of this paper, and further assuming, for simplicity, that C(t) is in fact an entire function of t (as are, for example, all parametrizations C(t) given by vector Fourier series containing finitely many terms), we introduce, for \(x=(x_1,x_2)\in \Omega \), the quantities

$$\begin{aligned} {\lambda }(x) = \sup \left\{ s\ge 0 : x\ne C(t+ ir) \text{ for } \text{ all } r\hbox { with } |r|\le s \hbox { and for all } t \in [0,2\pi ]\right\} \end{aligned}$$

and

$$\begin{aligned}&B_n(x_1,x_2,s)\nonumber \\&\quad = -\frac{1}{{4}\pi } \int \limits _0^{2\pi } \log \nonumber \\&\qquad \times \left[ \left( x_1 - C_1(t+is{\text {sgn}}(ns) ) \right) ^2 + \left( x_2 -C_2(t+is{\text {sgn}}(ns) ) \right) ^2 \right] e^{int}dt.\qquad \end{aligned}$$
(5.7)

Using Cauchy’s Theorem for \(x=(x_1,x_2)\in \Omega \) and any \(s\in {\mathbb {R}}\) satisfying \(|s| \le {\lambda }(x)\), we obtain

$$\begin{aligned} B_n^0(x_1,x_2)=e^{-|ns|} B_n(x_1,x_2,s), \end{aligned}$$
(5.8)

and, thus, letting \(s = \alpha {\lambda }(x)\) for any \(\alpha \in {\mathbb {R}}\) satisfying \(|\alpha |\le 1\), the eigenfunction u is given by

$$\begin{aligned} u(x_1,x_2) = {\overline{\phi }}+\sum _{\begin{array}{c} n\in {\mathbb {Z}}\\ n\ne 0 \end{array}} A_n e^{-|n\alpha | {\lambda }(x_1,x_2)} B_n(x_1,x_2,\alpha {\lambda }(x_1,x_2)) \end{aligned}$$
(5.9)

Lemma 5.2

There is \(C>0\) such that for all \(n>0\),

$$\begin{aligned} |B_n(x_1,x_2,\lambda (x_1,x_2))|\le \frac{C}{1+|n|}. \end{aligned}$$

Moreover, there is \(c>0\) and a sequence \(\{n_k\}_{k=1}^\infty \) with \(|n_k|\rightarrow \infty \) such that

$$\begin{aligned} |B_{n_k} ( x_1,x_2,{\lambda }(x_1,x_2) ) |\ge \frac{c}{n_k}. \end{aligned}$$
(5.10)

A proof of Lemma 5.2 is given in Appendix B. It follows from Lemma 5.2 that Eq. (5.8) optimally captures the exponential decay of the \(B_n\) terms as \(\sigma \rightarrow \infty \). Note that this setup does not capture the exponential decay of the coefficients \(A_n\) below machine precision away from \(|n|\sim \sigma \), and, therefore, the accuracy of the resulting interior eigenfunction reconstructions does not exceed that accuracy level. But the function \(\lambda (x_1,x_2)\) does capture the exponential decay and the geometrical character of the eigenfunction as long as the (spatially constant) coefficients \(A_n\) for low n remain above machine precision.

For general curves C(t) no closed form expressions exist for the function \({\lambda }(x)\), and a numerical algorithm must be used for the evaluation of this quantity, as part of a numerical implementation of the eigenfunction expression (5.9). In our implementation the function \({\lambda }\) was evaluated via an application of Newton’s method to the nonlinear equation

$$\begin{aligned} h(z)=(x_1-C_1(z))^2+(x_2-C_2(z))^2=0. \end{aligned}$$

Explicit expressions can be obtained for circles and ellipses, however:

  1. (1)

    For a circle of radius 1:

    $$\begin{aligned} {\lambda }(x_1,x_2) = -\log \left( \sqrt{x_1^2 + x_2^2}\right) . \end{aligned}$$
  2. (2)

    For an ellipse of semiaxes \(a>b\):

    $$\begin{aligned} {\lambda }(x_1,x_2) = \mathrm {arcosh}\left( \frac{a}{\sqrt{a^2-b^2} }\right) - \mathrm {Re}\left\{ \mathrm {arcosh} \left( \frac{x_1+ix_2}{\sqrt{a^2-b^2}}\right) \right\} . \end{aligned}$$
    (5.11)

The derivation of the expression (5.11) is outlined in Appendix A.

5.3 Exponential Decay and Verification of Cauchy’s Theorem

Table 1 Verification of the Cauchy-theorem-based identity (5.8) for the domain \(\Omega \) bounded by the elliptical curve (6.1) with \(a=2\) and \(b=1\)
Table 2 Same as Fig. (1) but for the kite-shaped domain \(\Omega \) bounded by the curve (6.2)

Tables 1 and 2 demonstrate the validity of Eq. (5.8) (since in both cases the results in the second and third columns closely agree with each other for \(n\le 50\)), as well as the exponential decay of the exact coefficients \(B_n^0\)—as born by the results in the third column of these tables. The disagreement observed for \(n>50\) is caused by the lack of precision of the results in the second column beyond machine accuracy, a problem that is eliminated in the third column via an application of the relation (5.8).

6 Numerical Results

Figures 4 and 5 present density plots and fixed-sign sets for Steklov eigenfunctions over domains bounded by the elliptical and kite-shaped curves parametrized by the vector functions

$$\begin{aligned} C(t) = ((a\cos (t),b\sin (t))\quad (0\le t <2\pi ) \end{aligned}$$
(6.1)

with \(a=2\) and \(b=1\), and

$$\begin{aligned} C(t) = (\cos (t) + 0.65 \cos (2t) - 0.65,1.5 \sin (t))\quad (0\le t <2\pi ), \end{aligned}$$
(6.2)

respectively. These figures demonstrate, in particular, domain-opening and non-density of nodal sets as discussed in Remark 1.2.

Fig. 4
figure 4

Density-plots (first and third rows) and fixed-sign sets (second and forth rows) for Steklov eigenfunctions over the elliptical domain (6.1). The eigenfunctions of orders 57 and 81 demonstrate the onset of the asymptotic character. In particular, regions of asymptotically fixed size open up. (As indicated in Sect. 6, the Steklov problem for the ellipse does not reduce to separated variables. Assuming Conjecture 1.3 is valid, further, the Steklov problem only reduces to separation of variables for a bounded analytic and simply-connected domain \(\Omega \) in case \(\Omega \) equals a disc)

Fig. 5
figure 5

Density-plots (first row) and fixed-sign sets (second row) for Steklov eigenfunctions over the kite-shaped domain (6.2)

As suggested in the caption to Fig. 4, the question might be considered as to whether the Steklov problem for the ellipse can be reduced to solution of decoupled separated-variables one-dimensional problems in elliptical coordinates. As mentioned in that caption, however, a simple calculation shows that the Steklov problem does not separate in those variables. Note also that the observed domain-opening and non-density of nodal sets are inconsistent with separation of variables—since, if variables separate, a simple eigenfunction can be written as \(u_\sigma =G(\mu )T(\tau )\), \(\mu \in [0,\mu _0]\)\(\tau \in [0,2\pi )\) where \(\mu \) and \(\tau \) are as in Appendix A. But, since \(T(\tau )\) must have \(1/\sigma \) dense zeros set and \(|G(\mu _0)|>0\), we obtain that the nodal set is \(1/\sigma \) dense throughout the ellipse which is indeed inconsistent with the fixed sign domains presented in the figure. An analogous argument can be used for high energy eigenfunctions, and for an arbitrary domain for which reducibility to separation of variables is assumed, even in cases involving multiple eigenvalues, once the work of Shamma [14] describing the asymptotic structure of Steklov eigenfunctions is incorporated. Thus, any bounded analytic and simply-connected domain \(\Omega \) for which the Steklov eigenproblem can be reduced to problems of separated variables cannot exhibit asymptotic domain-openings. It follows that any such domain \(\Omega \) cannot be a tunneling domain and must thus equal a disc—at least if Conjecture 1.3 holds.