Subwavelength guided modes for acoustic waves in bubbly crystals with a line defect

TL;DR

This paper demonstrates that a line defect in a bubbly crystal acts as a waveguide for specific acoustic frequencies, supported by a new formula for defect mode frequencies and numerical validation.

Contribution

It introduces an original formula for defect mode frequencies in bubbly crystal waveguides and numerically verifies wave localization along the line defect.

Findings

Line defect creates localized waveguiding for certain frequencies.

New analytical formula for defect mode frequencies.

Numerical simulations confirm waveguiding behavior.

Abstract

The recent development of subwavelength photonic and phononic crystals shows the possibility of controlling wave propagation at deep subwavelength scales. Subwavelength bandgap phononic crystals are typically created using a periodic arrangement of subwavelength resonators, in our case small gas bubbles in a liquid. In this work, a waveguide is created by modifying the sizes of the bubbles along a line in a dilute two-dimensional bubbly crystal, thereby creating a line defect. Our aim is to prove that the line defect indeed acts as a waveguide; waves of certain frequencies will be localized to, and guided along, the line defect. The key result is an original formula for the frequencies of the defect modes. Moreover, these frequencies are numerically computed using the multipole method, which numerically illustrates our main results.

Equations

$$\{\begin{array}{ll}{\Gamma }^k(x,y)=-\frac{i}{4}{H}_0^{(1)}(k|x-y|),& k>0,\\ \text{vskip}6.0ptplus2.0ptminus2.0pt& \\ {\Gamma }^0(x,y)=\frac{1}{2\pi }\ln |x-y|,& k=0,\end{array}$$

$${\mathcal{S}}_D^k[\varphi ](x)={\int }_{\partial D}{\Gamma }^k(x,y)\varphi (y)\mathrm{d}\sigma (y),x\in {\mathbb{R}}^2.$$

$${\mathcal{K}}_D^{k,\ast }[\varphi ](x)={\int }_{\partial D}\frac{\partial }{\partial {\nu }_x}{\Gamma }^k(x,y)\varphi (y)\mathrm{d}\sigma (y),x\in \partial D.$$

$$\mathcal{S}_{D}^{k}[\phi]\big{|}_{+}=\mathcal{S}_{D}^{k}[\phi]\big{|}_{-},$$

$$\frac{\partial}{\partial\nu}\mathcal{S}_{D}^{k}[\phi]\bigg{|}_{\pm}=\left(\pm\frac{1}{2}I+\mathcal{K}_{D}^{k,*}\right)[\phi].$$

$$-\frac{i}{4}{H}_0(k|x-y|)=\frac{1}{2\pi }\ln |x-y|+{\eta }_k+\sum _{j=1}^{\infty }\left(b_j\ln (k|x-y|)+c_j\right)(k|x-y|{)}^{2j},$$

$${\eta }_k=\frac{1}{2\pi }(\ln k+\gamma -\ln 2)-\frac{i}{4},b_j=\frac{(-1{)}^j}{2\pi }\frac{1}{2^{2j}(j!{)}^2},c_j=b_j\left(\gamma -\ln 2-\frac{i\pi }{2}-\sum _{n=1}^j\frac{1}{n}\right),$$

$${\hat{S}}_D^k[\varphi ](x)={\mathcal{S}}_D[\varphi ](x)+{\eta }_k{\int }_{\partial D}\varphi \mathrm{d}\sigma .$$

$${\mathcal{S}}_D^k={\hat{\mathcal{S}}}_D^k+O(k^2\ln k).$$

$$({\Delta }_x+k^2){\Gamma }^{\alpha ,k}(x,y)=\sum _{n\in {\mathbb{R}}^2}\delta (x-y-n)e^{in\cdot \alpha },x,y\in Y,$$

$${\mathcal{S}}_D^{\alpha ,k}[\varphi ](x)={\int }_{\partial D}{\Gamma }^{\alpha ,k}(x,y)\varphi (y)\mathrm{d}\sigma (y),x\in {\mathbb{R}}^2.$$

$$\mathcal{S}_{D}^{\alpha,k}[\phi]\big{|}_{+}=\mathcal{S}_{D}^{\alpha,k}[\phi]\big{|}_{-},$$

$$\frac{\partial}{\partial\nu}\Big{|}_{\pm}\mathcal{S}_{D}^{\alpha,k}[\phi]=\left(\pm\frac{1}{2}I+(\mathcal{K}_{D}^{-\alpha,k})^{*}\right)[\phi]\quad\mbox{on}\leavevmode\nobreak\ \partial D,$$

$$({\mathcal{K}}_D^{-\alpha ,k}{)}^{\ast }[\varphi ](x)={\int }_{\partial D}\frac{\partial }{\partial {\nu }_x}{\Gamma }^{\alpha ,k}(y,y)\varphi (y)\mathrm{d}\sigma (y).$$

$$\mathcal{F}[f](x_1,\alpha )=\sum _{m\in \mathbb{Z}}f(x_1-m)e^{i\alpha m},$$

$${\mathcal{F}}^{-1}[g](x_1)=\frac{1}{2\pi }{\int }_{Y^{\ast }}g(x_1,\alpha )\mathrm{d}\alpha .$$

$$v_w=\sqrt{\frac{{\kappa }_w}{{\rho }_w}},v_b=\sqrt{\frac{{\kappa }_b}{{\rho }_b}},k_w=\frac{\omega }{v_w}\text{and}{k}_b=\frac{\omega }{v_b}$$

$$\rho (x)={\rho }_b{\chi }_{\mathcal{C}}(x)+{\rho }_w(1-{\chi }_{\mathcal{C}}(x)),\kappa (x)={\kappa }_b{\chi }_{\mathcal{C}}(x)+{\kappa }_w(1-{\chi }_{\mathcal{C}}(x)),$$

$$\delta =\frac{{\rho }_b}{{\rho }_w}\ll 1.$$

$$v_w=v_b=1.$$

$$\kappa (x)\nabla \cdot \left(\frac{1}{\rho (x)}\nabla v(x)\right)+{\omega }^{2}v(x)=0,x\in {\mathbb{R}}^2.$$

$$-\kappa \nabla \cdot \frac{1}{\rho }\nabla .$$

$$\{\begin{array}{l}\kappa (x)\nabla \cdot \left(\frac{1}{\rho (x)}\nabla v(x)\right)+{\omega }^{2}v(x)=0,x\in Y_s,\\ e^{-i{\alpha }_1{x}_1}u\text{mbox}isperiodicin{x}_1,\end{array}$$

$$\{\begin{array}{l}\kappa (x)\nabla \cdot \left(\frac{1}{\rho (x)}\nabla v(x)\right)+{\omega }^{2}v(x)=0,x\in Y^2,\\ e^{-i\alpha \cdot x}u\text{mbox}isperiodicinx,\end{array}$$

$$0\le {\omega }_1^{\alpha }\le {\omega }_2^{\alpha }\le \cdots ,$$

$${\Lambda }_{0,{\alpha }_1}$$

$${\Lambda }_0$$

$$v=\{\begin{array}{ll}{\mathcal{S}}_D^{\alpha ,k_w}[{\phi }^{\alpha }]& \text{in}\text{leavevmode}{Y}^2\setminus \overline{D},\\ {\mathcal{S}}_D^{k_b}[{\psi }^{\alpha }]& \text{in}\text{leavevmode}D,\end{array}$$

$${\mathcal{A}}^{\alpha }(\omega ,\delta )[{\Phi }^{\alpha }]=0,$$

$${\mathcal{A}}^{\alpha }(\omega ,\delta )=\left(\begin{array}{cc}{\mathcal{S}}_D^{k_b}& -{\mathcal{S}}_D^{\alpha ,k}\\ -\frac{1}{2}+{\mathcal{K}}_D^{k_b,\ast }& -\delta \left(\frac{1}{2}+{\left({\mathcal{K}}_D^{-\alpha ,k}\right)}^{\ast }\right)\end{array}\right),{\Phi }^{\alpha }=\left(\begin{array}{c}{\phi }^{\alpha }\\ {\psi }^{\alpha }\end{array}\right).$$

$${\omega }_1^{\alpha }=\sqrt{\frac{\delta {\mathrm{Cap}}_{D,\alpha }}{|D|}}+O({\delta }^{3/2}),$$

Full text

Mathematics Subject Classification (MSC2000).  35R30, 35C20.

Keywords.  bubble, subwavelength resonance, subwavelength phononic crystal, subwavelength waveguide, line defect, weak localization.

1 Introduction

Line defects in bandgap photonic or phononic bandgap crystals are of interest due to their possible applications in low-loss waveguides. The main mathematical problem of interest is to show that the spectrum of the defect operator has a non-zero overlap with the original bandgap. Moreover, it is also of interest to understand the nature and location of the defect spectrum. For previous works regarding line defects in bandgap crystals we refer to [11, 23, 24, 13, 14, 15, 16, 18].

In this work, we consider a line defect in a phononic bandgap crystal comprised of gas bubbles in a liquid. The gas bubbles are known to resonate at a low frequency, called the Minnaert frequency. The corresponding wavelength is larger than the bubble by several orders of magnitude [1, 30]. Based on this, it is possible to create subwavelength bandgap crystals, which operate at wavelengths much larger than the unit cell size of the microstructured material. One of the main motivations for studying subwavelength bandgap materials is to manipulate wave propagation at subwavelength scales. A second motivation is for their use in devices where conventional bandgap materials, based on Bragg scattering, would create infeasibly large devices [28, 33]. Mathematical properties of bubbly phononic bandgap materials have been studied in, for example, [1, 3, 6, 7, 8, 12], and subwavelength phononic bandgap materials have been experimentally realised in [25, 26, 27].

Wave localization due to a point defect in a bubbly bandgap material was first proven in [3]. In [2], where some additions and minor corrections to [3] were made, it is shown that the mechanism for creating localized modes using small perturbations is quite different depending on the volume fraction of the bubbles. In order to create localized modes in the dilute regime, the defect should be smaller than the surrounding bubbles, while in the non-dilute regime, the defect has to be larger. Based on this, in the case of a line defect, it is natural to expect different behaviour in these two different regimes. This suggests that different methods of analysis are needed in the two regimes. In this paper, we will mainly focus on the dilute regime, taking the radius of the bubbles sufficiently small.

If the defect size is small, i.e. if the size of the perturbed bubble is close to its original size, then the band structure of the defect problem will be a small perturbation of the band structure of the original problem [4, 9]. This way, it is possible to shift the defect band upwards, and a part of the defect band will fall into the subwavelength bandgap. However, because of the curvature of the original band, it is impossible to create a defect band entirely inside the bandgap with this approach.

In order to create defect bands which are entirely located inside the subwavelength bandgap, we have to consider slightly larger perturbations. In this paper, we will show that for arbitrarily small defects, a part of the defect band will lie inside the bandgap. Moreover, we will show that for suitably large perturbation sizes, the entire defect band will fall into the bandgap, and we will explicitly quantify the size of the perturbation needed in order to achieve this. Because of this, our results are more general than previous weak localization results since we explicitly show how the defect band depends on the perturbation size.

In order to have guided waves along the line defect, the defect mode must not only be localized to the line, but also propagating along the line. In other words, we must exclude the case of standing waves in the line defect, i.e. modes which are localized in the direction of the line. As discussed in [22, 23], such modes are associated with the point spectrum of the perturbed operator which appears as a flat band in the dispersion relation. Proving the absence of bound modes in phononic or photonic waveguides is a challenging problem; for example in [32] this was proven by imposing “hard-wall” Dirichlet or Neumann boundary conditions along the waveguide, while in [20] the absence of bound modes was proven in the case of a simpler Helmholtz-type operator. In this paper, we use the explicit formula for the defect band to show that it is nowhere flat, and hence does not correspond to bound modes in the direction of the line.

The paper is structured as follows. In Section 2 we discuss preliminary results on layer potentials, and outline the main results from [7]. In Section 3 we restrict to circular domains and follow the approach of [3, 2] to model the line defect using the fictitious source superposition method, originally introduced in [34]. In Section 4 we prove the existence of a defect resonance frequency, and derive an asymptotic formula in terms of the density contrast in the dilute regime. Using this formula, we show that the defect modes are localized to, and guided along, the line defect. In Section 5 we compute the defect band numerically, in order to verify the formula and also illustrate the behaviour in the non-dilute regime. The paper ends with some concluding remarks in Section 6. In Appendix A, we restrict ourselves to small perturbations to derive an asymptotic formula valid in the non-dilute regime. In Appendix B we outline the fictitious source superposition method in the case of non-circular domains.

2 Preliminaries

2.1 Layer potentials

Let $Y^{2}=[-1/2,1/2)^{2}\subset\mathbb{R}^{2}$ be the unit cell and assume that the bubble occupies a bounded and simply connected domain $D\in Y^{2}$ with $\partial D\in C^{1,s}$ for some $0<s<1$. Let $\Gamma^{0}$ and $\Gamma^{k},k>0$ be the Green’s functions of the Laplace and Helmholtz equations in dimension two, respectively, i.e.,

[TABLE]

where $H_{0}^{(1)}$ is the Hankel function of the first kind and order zero. Here, the outgoing Sommerfeld radiation condition is used for selecting the physical Helmholtz Green’s function [4].

Let $\mathcal{S}_{D}^{k}:L^{2}(\partial D)\rightarrow H_{\textrm{loc}}^{1}(\mathbb{R}^{2})$ be the single layer potential defined by

[TABLE]

Here, $H_{\textrm{loc}}^{1}(\mathbb{R}^{2})$ denotes the space of functions that, on every compact subset of $\mathbb{R}^{2}$, are square integrable and have a weak first derivative that is also square integrable.

We also define the Neumann-Poincaré operator $\mathcal{K}_{D}^{k,*}:L^{2}(\partial D)\rightarrow L^{2}(\partial D)$ by

[TABLE]

The following so-called jump relations of $\mathcal{S}_{D}^{k}$ on the boundary $\partial D$ are well-known (see, for example, [4]):

[TABLE]

and

[TABLE]

Here, $\partial/\partial\nu$ denotes the outward normal derivative, and $|_{\pm}$ denote the limits from outside and inside $D$. In two dimensions, we have the following expansion of the Green’s function for the Helmholtz equation [4]

[TABLE]

where $\ln$ is the principal branch of the logarithm and

[TABLE]

with $\gamma$ being the Euler constant. Define, for $\phi\in L^{2}(\partial D)$,

[TABLE]

Then the following expansion holds:

[TABLE]

We also introduce a quasi-periodic version of the layer potentials. For $\alpha\in[0,2\pi)^{2}$, the quasi-periodic Green’s function $\Gamma^{\alpha,k}$ is defined to satisfy

[TABLE]

where $\delta$ is the Dirac delta function. The function $\Gamma^{\alpha,k}$ is $\alpha$-quasi-periodic in $x$, i.e., $e^{-i\alpha\cdot x}\Gamma^{\alpha,k}(x,y)$ is periodic in $x$ with respect to $Y$.

We define the quasi-periodic single layer potential $\mathcal{S}_{D}^{\alpha,k}$ by

[TABLE]

It satisfies the following jump formulas:

[TABLE]

and

[TABLE]

where $(\mathcal{K}_{D}^{-\alpha,k})^{*}$ is the operator given by

[TABLE]

We recall that $\mathcal{S}_{D}^{\alpha,0}:L^{2}(\partial D)\rightarrow H^{1}(\partial D)$ is invertible for $\alpha\neq 0$ [4].

2.2 Floquet transform

A function $f(x_{1})$ is said to be $\alpha$-quasi-periodic in the variable $x_{1}\in\mathbb{R}$ if $e^{-i\alpha x_{1}}f(x_{1})$ is periodic. Given a function $f\in L^{2}(\mathbb{R})$, the Floquet transform in one dimension is defined as

[TABLE]

which is $\alpha$-quasi-periodic in $x_{1}$ and periodic in $\alpha$. Let $Y=[-1/2,1/2)$ be the unit cell and $Y^{*}:=\mathbb{R}/2\pi\mathbb{Z}\simeq[0,2\pi)$ be the Brillouin zone. The Floquet transform is an invertible map $\mathcal{F}:L^{2}(\mathbb{R})\rightarrow L^{2}(Y\times Y^{*})$, with inverse (see, for instance, [21, 4])

[TABLE]

2.3 Bubbly crystals and subwavelength bandgaps

Here we briefly review the subwavelength bandgap opening of a bubbly crystal from [7].

Assume that a single bubble occupies the region $D$ specified in Section 2.1. We denote by $\rho_{b}$ and $\kappa_{b}$ the density and the bulk modulus inside the bubble, respectively. We let $\rho_{w}$ and $\kappa_{w}$ be the corresponding parameters outside the bubble. We introduce

[TABLE]

as the speed of sound outside and inside the bubbles, and the wavenumber outside and inside the bubbles, respectively. Here, $\omega$ corresponds to the operating frequency of the acoustic waves. Let $\mathcal{C}=\cup_{n\in\mathbb{Z}^{2}}(D+n)$ be the periodic bubbly crystal. Define, for $x\in\mathbb{R}^{2}$,

[TABLE]

where $\chi_{\mathcal{C}}$ is the characteristic function of $\mathcal{C}$.

We assume that there is a large contrast in the density, that is, the density contrast $\delta$ satisfies

[TABLE]

Recall that under (2.3), there exists a subwavelength resonance of the bubble in free space [1].

In the following, we shall also make the assumption stated below.

Assumption 2.1**.**

Without loss of generality, we assume that

[TABLE]

In this case we have $k_{b}=k_{w}=\omega$. Assumption 2.1 only serves to simplify the expressions. The methods presented in this paper indeed apply as long as the wave speeds outside and inside the bubbles are comparable to each other.

The wave propagation problem inside the periodic crystal can be modelled as

[TABLE]

We denote by $\Lambda_{0}$ the set of propagating frequencies, i.e., the set of $\omega$ such that $\omega^{2}$ is in the spectrum of the operator

[TABLE]

Denote by $Y_{s}=Y\times\mathbb{R}$ the unit strip and recall that $Y^{2}=[-1/2,1/2)^{2}$ is the unit cell of the crystal. Applying the Floquet transformation, first in $x_{1}$-direction and then in $x_{2}$-direction, equation (2.4) can be decomposed first as

[TABLE]

where $\alpha_{1}\in Y^{*}$, and then as

[TABLE]

where $\alpha=(\alpha_{1},\alpha_{2})\in Y^{*}\times Y^{*}$. We denote by $\Lambda_{0,\alpha_{1}}$ the set of $\omega$ such that $\omega^{2}$ is in the spectrum of the operator implied by (2.5) and by $\Lambda^{ess}_{0,\alpha_{1}}$ the essential part of this spectrum. It is known that (2.6) has non-trivial solutions for discrete values of $\omega$:

[TABLE]

and we have the following band structure of propagating frequencies for the periodic bubbly crystal $\mathcal{C}$:

[TABLE]

In [7], it is proved that there exists a subwavelength spectral gap opening in the band structure. Let us briefly review this result. We look for a solution $v$ of (2.6) which has the following form:

[TABLE]

for some densities $\varphi^{\alpha},\psi^{\alpha}\in L^{2}(\partial D)$. Using the jump relations for the single layer potentials, one can show that (2.6) is equivalent to the boundary integral equation

[TABLE]

where

[TABLE]

Since it can be shown that $\omega=0$ is a characteristic value for the operator-valued analytic function $\mathcal{A}(\omega,0)$, we can conclude the following result by the Gohberg-Sigal theory [4, 19].

Lemma 2.1**.**

For any $\delta$ sufficiently small, there exists a characteristic value $\omega_{1}^{\alpha}=\omega_{1}^{\alpha}(\delta)$ to the operator-valued analytic function $\mathcal{A}^{\alpha}(\omega,\delta)$ such that $\omega_{1}^{\alpha}(0)=0$ and $\omega_{1}^{\alpha}$ depends on $\delta$ continuously.

The next theorem gives the asymptotic expansion of $\omega_{1}^{\alpha}$ as $\delta\rightarrow 0$.

Theorem 2.1**.**

[7]* For $\alpha\neq 0$ and sufficiently small $\delta$, we have*

[TABLE]

where the constant ${\mathrm{Cap}}_{D,\alpha}$ is given by

[TABLE]

Here, $\langle\,\cdot\,,\,\cdot\,\rangle$ stands for the standard inner product of $L^{2}(\partial D)$ and $\chi_{\partial D}$ denotes the characteristic function of $\partial D$.

Let $\omega_{1}^{*}=\max_{\alpha}\omega_{1}^{\alpha}$. The following theorem expresses the fact that a subwavelength bandgap opens in the band structure of the bubbly crystal.

Theorem 2.2**.**

[7]* For every $\varepsilon>0$, there exists $\delta_{0}>0$ and $\tilde{\omega}>\omega_{1}^{*}$ such that*

[TABLE]

for $\delta<\delta_{0}$.

3 Integral representation for bubbly crystals with a defect

3.1 Formulation of the line defect problem

In the following, we will consider the case when all the bubbles are circular disks. This gives a convenient presentation, and makes the problem similar to the point defect problem studied in [3, 2]. In Appendix B, we will outline the analysis in the case of non-circular bubbles.

Consider a perturbed crystal, where all the disks along the $x_{1}$-axis are replaced by defect disks of radius $R_{d}$ with $0<R_{d}<R$. Denote the centre defect disk by $D_{d}$ and let

[TABLE]

be the perturbed crystal, depicted in Figure 1. Moreover, let $\varepsilon=R_{d}-R<0,\varepsilon\in(-R,0)$ be the perturbation of the radius. Define

[TABLE]

The wave propagation problem inside the periodic crystal can be modelled as

[TABLE]

We denote by $\Lambda_{d}$ the set of propagating frequencies in the line defect crystal, i.e. the set of $\omega$ such that $\omega^{2}$ is in the spectrum of the operator

[TABLE]

Since the defect crystal is periodic in the $x_{1}$-direction, we can use the Floquet transformation to decompose (3.1) as

[TABLE]

where $\alpha_{1}\in Y^{*}$ and $Y_{s}$ again denotes the strip $Y_{s}=[-1/2,1/2)\times\mathbb{R}$. We will denote by $\Lambda_{d,\alpha_{1}}$ the set of $\omega$ such that $\omega^{2}$ is in the spectrum of the operator implied by (3.2) and by $\Lambda^{ess}_{d,\alpha_{1}}$ the corresponding essential part of the spectrum.

In the strip $Y_{s}$, the perturbations $\rho_{d}-\rho$ and $\kappa_{d}-\kappa$ have compact support. Since the essential spectrum is stable under compact perturbations [17, 31], it can be shown that the essential spectra $\Lambda^{ess}_{0,\alpha_{1}}$ and $\Lambda^{ess}_{d,\alpha_{1}}$ coincide.

In this paper, we want to show that introducing the line defect creates a defect band $\omega^{\varepsilon}(\alpha_{1})\notin\Lambda_{0,\alpha_{1}}$. Moreover, we want to show that $\varepsilon$ can be chosen such that $\omega^{\varepsilon}(\alpha_{1})\notin\Lambda_{0}$ for all $\alpha_{1}\in Y^{*}$, which means that any Bloch mode is localized to the line defect. We also want to show that $\omega^{\varepsilon}(\alpha_{1})$ is not contained in the pure point part of $\Lambda_{d,\alpha_{1}}$, which means that there are no bound modes in the defect direction.

3.2 Effective sources for the defect

Here we describe an effective sources approach to the solution of (3.2) in the strip. The idea is to model the defect bubble $D_{d}$ as an unperturbed bubble $D$ with additional fictitious monopole and dipole sources $f$ and $g$. This method was originally introduced in [34] and then it was applied in [3, 2] for a point defect in a bubbly crystal.

Let us consider the following problem:

[TABLE]

where $f$ and $g$ are the source terms and $\delta_{m,n}$ is the Kronecker delta function. Note that the sources are present only on the boundary of the central bubble $D$.

We denote the solution to the original problem (3.2) by $u$ and the effective source solution (3.3) by $\widetilde{u}$. We want to find appropriate conditions on $f$ and $g$ in order to achieve

[TABLE]

Then $u$ can be recovered by extending $\widetilde{u}$ to the whole region including $D\setminus D_{d}$ with boundary conditions on $\partial D$ and $\partial D_{d}$. The conditions for the effective sources $f$ and $g$, which are necessary in order to correctly model the defect, will be characterized in the next subsection.

3.3 Characterization of the effective sources

Here we clarify the relation between the effective source pair $(f,g)$ and the layer density pair $(\varphi,\psi)$ defined in equation (3.5) below.

First, we observe that away from the central unit cell $Y^{2}$, the equations (3.2) and (3.3) satisfy the same geometric and quasi-periodic conditions. Thus, in order for (3.4) to hold, it is sufficient for $u$ and $\widetilde{u}$ to coincide inside the central unit cell $Y^{2}$.

Inside $Y^{2}$, the solution $\widetilde{u}$ can be represented as

[TABLE]

for some pair $(\varphi,\psi)\in L^{2}(\partial D)^{2}$, with $H$ satisfying the homogeneous equation $(\Delta+k_{w}^{2})H=0$ in $Y^{2}$. In (3.5), the local properties of $\widetilde{u}$ around $\partial D$ are given by the single-layer potentials, while $H$ can be chosen to make $\widetilde{u}$ satisfy the quasi-periodic condition. From the jump conditions given in Section 2.1, the pair $(\varphi,\psi)$ satisfies

[TABLE]

Similarly, inside $Y^{2}$, the solution $u$ can be represented as

[TABLE]

where

[TABLE]

Now, having the two solutions coincide inside $(Y^{2}\setminus D)\cup D_{d}$ is equivalent to the conditions

[TABLE]

and

[TABLE]

Assuming $D$ is a disk, the above equations were solved in [3, 2], and we state the results in Proposition 3.1 below. First, we introduce some notation. Since $D$ and $D_{d}$ are circular disks, we can use a Fourier basis for functions in $L^{2}(\partial D)$ or $L^{2}(\partial D_{d})$. For $n\in\mathbb{Z}$, define the subspace $V_{n}$ of $L^{2}(\partial D)$ as $V_{n}:=\mbox{span}\{e^{im\theta}\}$. Then define the subspace $V_{mn}$ of $L^{2}(\partial D)^{2}$ as

[TABLE]

Similarly, let ${V}_{mn}^{d}$ be the subspace of $L^{2}(\partial D_{d})^{2}$ with the same Fourier basis. Then it can be shown that the operator $\mathcal{A}_{D}$ in (3.6) has the following matrix representation as an operator from $V_{mn}$ to $V_{m^{\prime}n^{\prime}}$:

[TABLE]

Similarly, the operator $\mathcal{A}_{D_{d}}$ in (3.7) is represented as follows:

[TABLE]

In [3, 2], the following proposition was shown.

Proposition 3.1**.**

The density pair $(\varphi,\psi)$ and the effective sources $(f,g)$ satisfy the following relation

[TABLE]

where the operators $\mathcal{P}_{1}:L^{2}(\partial D)^{2}\rightarrow L^{2}(\partial D_{d})^{2}$ and $\mathcal{P}_{2}:L^{2}(\partial D)^{2}\rightarrow L^{2}(\partial D_{d})^{2}$ are defined by

[TABLE]

and $\mathcal{A}_{D}^{\varepsilon}$ is defined as

[TABLE]

3.4 Floquet transform of the solution

In view of Proposition 3.1, we can identify the solutions $u$ and $\widetilde{u}$. In this section, we derive an integral equation for the effective source problem (3.3). This problem is already quasi-periodically reduced in the $x_{1}$-direction, with quasi-periodicity $\alpha_{1}$. For some quasi-periodicity $\alpha_{2}\in Y^{*}$, we set $\alpha=(\alpha_{1},\alpha_{2})$ and apply the Floquet transform to the solution $u$ in the $x_{2}$-direction as follows:

[TABLE]

The transformed solution $u^{\alpha}$ satisfies

[TABLE]

The solution $u^{\alpha}$ is $\alpha$-quasi-periodic in the two-dimensional cell $Y^{2}$, and can be represented using quasi-periodic layer potentials as

[TABLE]

where, similarly as in equation (2.7), the pair $(\varphi^{\alpha},\psi^{\alpha})\in L^{2}(\partial D)^{2}$ is the solution to

[TABLE]

Since the operator $\mathcal{A}^{\alpha}$ is invertible for small enough $\delta$ and for $\omega$ inside the bandgap [7], we have

[TABLE]

The solution $u$ to problem (3.3) can be recovered by the inversion formula as

[TABLE]

Now, by the same arguments as those in [3, 2], we obtain the following proposition.

Proposition 3.2**.**

The density pair $(\varphi,\psi)$ and the effective source pair $(f,g)$ satisfy

[TABLE]

for small enough $\delta$ and for $\omega\notin\Lambda_{0,\alpha_{1}}$ inside the bandgap.

3.5 The integral equation for the layer densities

Here we state the integral equation for the layer density pair $(\phi,\psi)$. The following result is an immediate consequence of Propositions 3.1 and 3.2.

Proposition 3.3**.**

The layer density pair $(\phi,\psi)\in L^{2}(\partial D)^{2}$ satisfies the following integral equation:

[TABLE]

for small enough $\delta$ and for $\omega\notin\Lambda_{0,\alpha_{1}}$ inside the bandgap.

The expression of this integral equation resembles the one for a point defect found in [3, 2]. However, this similarity is not obvious, and can be seen as a consequence of the cancellation of $H$ in Proposition 3.1.

The significance of Proposition 3.3 is as follows. If we can show that there is a characteristic value $\omega=\omega^{\varepsilon}$ of $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ inside the bandgap, i.e. if there is a non-trivial pair $(\phi,\psi)$ such that $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega^{\varepsilon})\left(\begin{smallmatrix}\phi\\ \psi\end{smallmatrix}\right)=0$, then $\omega^{\varepsilon}$ is a resonance frequency for the defect mode.

4 Subwavelength guided modes in the defect

Here, we will prove the existence of a resonance frequency $\omega=\omega^{\varepsilon}(\alpha_{1})$ inside the bandgap of the unperturbed crystal at $\alpha_{1}$. We will give an asymptotic formula for $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ in terms of $\delta$ in the dilute regime. Moreover, we will show that the defect band is not contained in the pure point spectrum of the defect operator, and for perturbation sizes $\varepsilon$ larger than some critical $\varepsilon_{0}$, the entire defect band is located in the bandgap region of the original operator.

4.1 Asymptotic expansions for small $\delta$

In this section, we will asymptotically expand $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ in the limit as $\delta\rightarrow 0$ and with $\omega$ in the subwavelength regime, i.e., $\omega=O(\sqrt{\delta})$. Throughout this section we assume that $\alpha\neq(0,0)$. We begin by studying the operator $(\mathcal{A}^{\alpha}(\omega,\delta))^{-1}$.

Define $\psi_{\alpha}$ as

[TABLE]

Since we know that \big{(}\frac{1}{2}I+(\mathcal{K}_{D}^{-\alpha,0})^{*}\big{)}[\psi_{\alpha}]=\psi_{\alpha}, we can decompose this operator as

[TABLE]

where

[TABLE]

is a projection on $\psi_{\alpha}$. Then it can be shown that $Q_{\alpha}[\psi_{\alpha}]=0$ and $Q_{\alpha}^{*}[\chi_{\partial D}]=0$, where $Q_{\alpha}^{*}$ is the adjoint of $Q_{\alpha}$.

For small $\delta$ and for $\omega=O(\sqrt{\delta})$ inside the corresponding bandgap, the operator $\mathcal{A}^{\alpha}(\omega,\delta)$ can be decomposed as

[TABLE]

Define the operators

[TABLE]

and

[TABLE]

The motivation for defining these operators is given in Lemmas 4.1 and 4.2. Introducing these operators enables the explicit computation of $(A^{\alpha})^{-1}$. We will compute the asymptotic expansion of these operators for small $\omega$ and $\delta$.

Lemma 4.1**.**

The following results hold for $A_{0}$ and $A_{1}$:

• (i)

For $\omega\neq 0$, $A_{0}:L^{2}(\partial D)\rightarrow L^{2}(\partial D)$ is invertible, and as $\omega\rightarrow 0$ and $\delta\rightarrow 0$,

[TABLE]

• (ii)

For $\omega\neq\omega^{\alpha}$, $A_{1}:L^{2}(\partial D)\rightarrow L^{2}(\partial D)$ is invertible, and as $\omega\rightarrow 0$ and $\delta\rightarrow 0$,

[TABLE]

where $P_{\alpha}^{\perp}=I-\frac{(\omega^{\alpha})^{2}}{\omega^{2}}P_{\alpha}$.

Proof of (i). We easily find that

[TABLE]

which is well-defined since $-\frac{1}{2}I+\mathcal{K}_{D}^{\omega,*}:L^{2}(\partial D)\rightarrow L^{2}(\partial D)$ is invertible for $\omega\neq 0$ [4]. From the low-frequency expansion of $\mathcal{S}_{D}^{\alpha,\omega}$ [4], and using the Neumann series, we have

[TABLE]

Using the Fourier basis, the operator $-\frac{1}{2}I+\mathcal{K}_{D}^{\omega,*}$ can be represented as [7]

[TABLE]

Using standard asymptotics we can compute

[TABLE]

Hence the operator $\left(-\frac{1}{2}I+\mathcal{K}_{D}^{\omega,*}\right)^{-1}$ can be written as

[TABLE]

Moreover, we have from (2.1) that $\mathcal{S}_{D}^{\omega}[\chi_{\partial D}]=(2\pi R\eta_{\omega}+R\ln R)\chi_{\partial D}+O(\omega^{2}\ln\omega)$, and so

[TABLE]

Combining equations (4.1), (4.1), (4.3) and (4.4) proves (i). ∎

Proof of (ii). Using the definition of $A_{1}$, and the expression for $A_{0}$, we can compute

[TABLE]

Recall the asymptotic expression of $\omega^{\alpha}$ given in Theorem 2.1:

[TABLE]

We then find that

[TABLE]

Define $P_{\alpha}^{\perp}=I-\frac{(\omega^{\alpha})^{2}}{\omega^{2}}P_{\alpha}$. For $\omega$ small enough, $A_{1}$ is invertible precisely when $P_{\alpha}^{\perp}$ is invertible, i.e. when $\omega\neq\omega^{\alpha}$. Moreover, we have

[TABLE]

This proves (ii). ∎

Lemma 4.2**.**

For $\omega\neq\omega^{\alpha}$, and as $\omega\rightarrow 0$ and $\delta\rightarrow 0$, we have

[TABLE]

Proof.

We have already established the invertibility of $A_{0}$ and $A_{1}$. Using this fact, we have

[TABLE]

Because $Q_{\alpha}^{*}\chi_{\partial D}=0$, we have that

[TABLE]

We then have that

[TABLE]

where the last step follows using the Neumann series. ∎

Next, we compute the operator $(\mathcal{A}_{D}^{\varepsilon}-\mathcal{A}_{D})$. Using Proposition 3.1 and equation (3.10), we have

[TABLE]

Consequently, the operator $(\mathcal{A}_{D}^{\varepsilon}-\mathcal{A}_{D})$ is given by

[TABLE]

Introduce the notation

[TABLE]

Using asymptotic expansions of the Bessel function $J_{n}(z)$ and the Hankel function $H_{n}^{(1)}(z)$, for small $z$, straightforward computations show that

[TABLE]

Moreover, we have

[TABLE]

We are now ready to compute the full operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$.

Proposition 4.1**.**

The operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)$ has the form

[TABLE]

where the operators $M_{0}(\omega),M_{1}(\omega):L^{2}(\partial D)\rightarrow L^{2}(\partial D)$ depend on $\varepsilon,\delta,\alpha_{1}$. Moreover, as $\omega\rightarrow 0$, $\delta\rightarrow 0$ and $\omega\notin\Lambda_{0,\alpha_{1}}$ we have

[TABLE]

Proof.

The expression of $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ given in (4.8) follows from equations (3.12) and (4.6). Combining Lemmas 4.1 and 4.2, we find that

[TABLE]

Combining this with equations (3.12), (4.6) and (4.7) yields the desired expression for $M_{0}(\omega)$. ∎

Remark 4.1**.**

It is clear that $\omega=\omega^{\varepsilon}$ is a characteristic value of $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ if and only if $\omega^{\varepsilon}$ is a characteristic value for $I+M_{0}$. We have thus reduced the characteristic value problem for the two-dimensional matrix operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ to the scalar operator $I+M_{0}$.

4.2 Defect resonance frequency in the dilute regime

The following theorem is the main result of this paper. Again, we say a frequency $\omega$ is in the subwavelength regime if $\omega=O(\sqrt{\delta})$.

Theorem 4.1**.**

For $\delta$ and $R$ small enough, there is a unique characteristic value $\omega^{\varepsilon}(\alpha_{1})$ of $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)$ such that $\omega^{\varepsilon}(\alpha_{1})\neq\Lambda_{0,\alpha_{1}}$ and $\omega^{\varepsilon}(\alpha_{1})$ is in the subwavelength regime. Moreover,

[TABLE]

where $\hat{\omega}$ is the root of the following equation:

[TABLE]

Proof.

We seek the characteristic values of the operator $I+M_{0}$. We consider the dilute regime, i.e. where $R$ is small. As shown in [7], in this case we have

[TABLE]

where $\mathcal{R}_{\alpha}(x)=\Gamma^{\alpha,0}(x)-\Gamma^{0}(x)$. In particular,

[TABLE]

We will compute $M_{0}$ in the Fourier basis. It is known that [4]

[TABLE]

which gives

[TABLE]

Moreover,

[TABLE]

In total, we have on the subspace $V_{0}$,

[TABLE]

Moreover, if $n\neq 0$, then

[TABLE]

In summary, the operator $I+M_{0}$ can be written as

[TABLE]

where the limiting operator $\hat{M}(\omega)$ is a diagonal operator in the Fourier basis, with non-zero diagonal entries for $n\neq 0$. We conclude that $\omega=\hat{\omega}$ is a characteristic value for $\hat{M}(\omega)$ if and only if one of the diagonal entries vanishes at $\omega=\hat{\omega}$, i.e. if

[TABLE]

Next, we show that equation (4.11) has a zero $\hat{\omega}\notin\Lambda_{0,\alpha_{1}}$ satisfying $\hat{\omega}=O(\sqrt{\delta})$. Introduce the notation

[TABLE]

then equation (4.11) implies

[TABLE]

For a fixed $\alpha_{1}\in Y^{*}$, define $\omega^{*}=\omega^{(\alpha_{1},\pi)},$ which is the edge of the first band in $\Lambda_{0,\alpha_{1}}$. Observe that ${1}/{I(\omega,\alpha_{1})}$ is monotonically increasing in $\omega$, and

[TABLE]

where $\omega_{0}^{2}$ is the average

[TABLE]

In the dilute regime, we can compute

[TABLE]

so as $\omega\rightarrow\omega^{*}$, the right-hand side of equation (4.12) tends to

[TABLE]

Since $R_{d}<R$, the leading-order term is negative. On the other hand, as $\omega\rightarrow\infty$, the right-hand side of (4.12) tends to $\infty$. Since the right-hand side of equation (4.12) is monotonically increasing, this equation has a unique zero $\omega=\hat{\omega}$. It can be verified that this zero has multiplicity one. Moreover, $\hat{\omega}$ satisfies $\hat{\omega}=O(\sqrt{\delta})$.

Now, we turn to the full operator $I+M_{0}(\omega)$. Since $I+M_{0}(\omega)=\hat{M}(\omega)+O\left(R^{2}+\omega\right),$ by the Gohberg-Sigal theory [4, 9, 19], close to $\hat{\omega}$ there is a unique characteristic value $\omega^{\varepsilon}$ of the operator $I+M_{0}(\omega)$, satisfying

[TABLE]

This concludes the proof. ∎

Remark 4.2**.**

In the case of $R_{d}>R$, i.e. larger defect bubbles, similar arguments show that any subwavelength frequency $\omega^{\varepsilon}(\alpha_{1})\in\Lambda_{d,\alpha_{1}}\setminus\Lambda_{0,\alpha_{1}}$ satisfies equation (4.9) in the dilute regime. However, it is easily verified that this equation has no solutions $\hat{\omega}>\omega^{*}$ in the case $R_{d}>R$. The conclusion is that we must reduce the size of the defect bubbles in order to create subwavelength guided modes in the dilute regime.

4.2.1 Absence of bound modes in the line defect direction

In this section, we will show that the defect band is not contained in the pure point spectrum of the defect crystal.

Lemma 4.3**.**

For $(\alpha_{1},\alpha_{2})\in Y^{*}\times Y^{*}$, $\alpha_{2}\neq 0$, the partial derivative of the quasi-periodic Green’s function

[TABLE]

is zero precisely when $\alpha_{1}=0$ or $\alpha_{1}=\pi$.

Proof.

From the spectral form of the Green’s function [4]:

[TABLE]

it can be easily shown that

[TABLE]

By symmetry of the summation, we find that

[TABLE]

if and only if $\alpha_{1}=0$ or $\alpha_{1}=\pi$. ∎

Proposition 4.2**.**

For $\delta$ and $R$ small enough, and for $\alpha_{1}\neq 0,\pi$, the characteristic value $\omega^{\varepsilon}=\omega^{\varepsilon}(\alpha_{1})$ satisfies

[TABLE]

Proof.

To simplify the computations, we introduce the following notation:

[TABLE]

Then equation (4.9) reads

[TABLE]

Denote by $x^{\prime}=\frac{\partial x}{\partial\alpha_{1}}$ and $y^{\prime}=\frac{\partial y}{\partial\alpha_{1}}$, then we have

[TABLE]

or equivalently,

[TABLE]

where

[TABLE]

and

[TABLE]

First, we show that $A\neq 0$ which implies that the zeros of $x^{\prime}$ coincides with the zeros of $B$. We have

[TABLE]

since $y(ay+b)>0$ for all $(\alpha_{1},\alpha_{2})\in Y^{*}\times Y^{*}$.

Next, we show that the leading order of $B$ vanishes exactly at the points $\alpha_{1}=0$ and $\alpha=\pi$. Using equations (4.5) and (4.10), we have

[TABLE]

Since $\mathcal{R}_{\alpha}=\Gamma^{\alpha,0}-\Gamma^{0}$, using Lemma 4.3 we conclude that for $\delta$ and $R$ small enough and $\alpha_{1}\neq 0,\pi$, $y^{\prime}$ is non-zero for any $\alpha_{2}$. Hence $B$ is non-zero, which concludes the proof. ∎

Proposition 4.2 shows that the defect dispersion relation is not flat, except for local extrema at $\alpha_{1}=0$ and $\alpha_{1}=\pi$. Thus, the defect band is not in the pure point spectrum of the defect operator, and corresponding Bloch modes are not bounded in the line defect direction.

4.2.2 Bandgap located defect bands

In this section, we will demonstrate that it is possible to position the entire defect band in the bandgap region with a suitable choice of $\varepsilon$. Recall that $\varepsilon=R_{d}-R$. As before, let

[TABLE]

Lemma 4.4**.**

For a fixed $\omega\notin\Lambda_{0}$, the minimum

[TABLE]

is attained at $\alpha_{1}=\alpha_{0}$ with $\alpha_{0}\rightarrow 0$ as $\delta\rightarrow 0$.

Proof.

We begin by observing that the minima of $I(\omega,\alpha_{1})$ and $\omega^{(\alpha_{1},\alpha_{2})}$ are attained at the same point $\alpha_{1}=\alpha_{0}\in Y^{*}$. Using Lemma 4.2, for every fixed $\alpha_{2}\neq 0$ the minimum of ${\mathrm{Cap}}_{D,\alpha}$ is attained at $\alpha_{1}=0$, so by Theorem 2.1 the minimum of $\omega^{(\alpha_{1},\alpha_{2})}$ is attained at $\alpha_{1}=\alpha_{0}$ with $\alpha_{0}\rightarrow 0$ as $\delta\rightarrow 0$. Since $\omega^{(0,0)}=0$ (see [7]) this is true for all $\alpha_{2}\in Y^{*}$. ∎

Proposition 4.3**.**

For $\delta$ small enough, there exists an $\varepsilon$ such that for all $\alpha_{1}\in Y^{*}$ we have

[TABLE]

Proof.

We want to show that

[TABLE]

Using Lemma 4.4, it is easy to see that

[TABLE]

is attained at $\alpha_{1}=\alpha_{0}$. Moreover, from [10] we know that

[TABLE]

is attained at $\alpha=\alpha^{*}=(\pi,\pi)$. Using Theorem 4.1, we conclude that the lower edge of the defect band coincides with the upper edge of the unperturbed band if

[TABLE]

For small enough $R$ and $\delta$, the right-hand side is positive, while the left-hand side ranges from [math] to $+\infty$ for $\varepsilon\in(-R,0)$. Hence we can find a solution $\varepsilon_{0}$ to this equation, and the statement holds for $\varepsilon>\varepsilon_{0}$. ∎

Remark 4.3**.**

In practice, for $\delta$ small enough, we can approximate $\varepsilon_{0}$ as the root to the equation

[TABLE]

5 Numerical illustrations

5.1 Implementation

5.1.1 Discretization of the operator

The operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)$ was approximated as a matrix $M(\omega)$ using the truncated Fourier basis $e^{-iN\theta}$, $e^{-i(N-1)\theta},\ldots,e^{iN\theta}$. We refer to [3, 2] for the details of the discretization. The integral over $Y^{*}$ in (3.12) was approximated using the trapezoidal rule with $100$ discretization points. The characteristic value problem for $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)$ was formulated as the root-finding problem $\det M(\omega)=0$ and solved using Muller’s method [4].

5.1.2 Evaluation of the asymptotic formula

The integral over $Y^{*}$ in equation (4.9) was approximated using the trapezoidal rule with $100$ discretization points. Again, the equation was numerically solved using Muller’s method.

5.2 Dilute regime

Figure 2 shows the unperturbed band structure and the defect band for $\alpha_{1}$ over the Brillouin zone $[0,2\pi]$. The material parameters were chosen as $\kappa_{b}=\rho_{b}=1$, $\kappa_{w}=\rho_{w}=5000$, $R=0.05$ and $\varepsilon=-0.2R$. It can be seen that the entire defect band is located inside the deep subwavelength regime of the bandgap. Moreover, the defect frequencies computed using the asymptotic formula agree well with the values computed by discretizing the operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$. Also, we see that the defect band is not flat. In summary, these results show that the defect crystal supports guided modes in the subwavelength regime, localized to the line defect.

5.2.1 Computation of $\varepsilon_{0}$

In this section, we numerically compute the critical perturbation size, where the entire defect band is located in the bandgap. The critical perturbation size was computed in two ways: by solving equation (4.13) for the leading order term, and by solving the root-finding problem $\omega^{\varepsilon_{0}}(0)=\omega^{*}$ where $\omega^{\varepsilon}$ was computed by discretizing the operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$.

Figure 3 shows $\varepsilon_{0}$ for different $R$ in the dilute regime. The material parameters were chosen as $\kappa_{b}=\rho_{b}=1$ and $\kappa_{w}=\rho_{w}=10000$. The values obtained from the asymptotic formula and by discretizing the operator agree, with a smaller radius $R$ giving a smaller error. Quantitatively, for $R$ in this regime, we require a decrease of the bubble size by around $14\%$ to $26\%$ in order that the defect band be located inside the bandgap.

5.3 Non-dilute regime

Here we compute the defect band in the non-dilute regime, in both cases $\varepsilon<0$ and $\varepsilon>0$, corresponding to smaller and larger defect bubbles, respectively. Theorem A.1 in Appendix A shows that there is a defect frequency $\omega^{\varepsilon}$ in the bandgap for small $\varepsilon>0$ but not for small $\varepsilon<0$.

5.3.1 Larger defect bubbles

Figure 4 shows the band structure in the non-dilute case with $\varepsilon>0$. The material parameters were chosen as $\kappa_{b}=\rho_{b}=1$, $\kappa_{w}=\rho_{w}=5000$, $R=0.4$ and $\varepsilon=0.45R$. As expected from Theorem A.1, there is a defect band above the first band of the unperturbed crystal. Moreover, it is possible to position the entire band inside the bandgap.

5.3.2 Smaller defect bubbles

Figure 5 shows the band structure in the non-dilute case with $\varepsilon=-0.6R$. The material parameters were chosen as $\kappa_{b}=\rho_{b}=1$, $\kappa_{w}=\rho_{w}=5000$, $R=0.4$ and $\varepsilon=0.45R$. In this case a defect band is present inside the bandgap. Here $\varepsilon$ is quite large, in contrast to Theorem A.1 which is only valid for small $\varepsilon$.

6 Concluding remarks

In this paper, we have for the first time proved the possibility of creating subwavelength guided waves localized to a line defect in a bubbly phononic crystal. We have shown that introducing a defect line, by shrinking the bubbles along the line, creates a defect frequency band inside the bandgap of the original crystal. An arbitrarily small perturbation will create a non-zero overlap between the defect band and the bandgap, and we have explicitly quantified the required defect size in order to position the entire defect band inside the bandgap. Moreover, we have shown for the first time that the defect band is not contained in the pure point spectrum of the perturbed operator. This shows that we can create truly guided modes, which are not localized in the direction of the defect. In the future, we plan to study more sophisticated waveguides, with bends and junctions. Moreover, we also plan to study waveguides in phononic subwavelength bandgap crystals with non-trivial topology, rigorously proving the existence of topologically protected subwavelength states in bubbly crystals.

Appendix A The resonance frequency of the defect mode for small perturbations

Here we derive a formula for the resonance frequency of the defect mode in the case of small $\varepsilon$, following the approach of [3, 2]. The strength of this approach is that it is valid in both the dilute and non-dilute regimes. We begin by reformulating the integral equation (3.12) in terms of the effective sources $(f,g)$ instead of the layer densities $(\phi,\psi)$. The following proposition is a restatement of Proposition 3.3.

Proposition A.1**.**

The effective source pair $(f,g)\in L^{2}(\partial D)^{2}$ satisfies the following integral equation:

[TABLE]

for small enough $\delta$ and for $\omega\notin\Lambda_{0,\alpha_{1}}$ inside the bandgap.

In this section, we derive an expression for the characteristic value $\omega^{\varepsilon}$ of $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)$ located slightly above $\omega^{\alpha}$ for both the dilute and non-dilute regimes.

Let us first analyse the operator $\int_{Y^{*}}(\mathcal{A}^{\alpha})^{-1}\>\mathrm{d}\alpha$. Since $\omega^{\alpha}$ is a simple pole of the mapping $\omega\mapsto\mathcal{A}^{\alpha}(\omega,)^{-1}$ in a neighbourhood of $\omega^{\alpha}$, according to [4], we can write

[TABLE]

where the operator-valued function $\mathcal{R}^{\alpha}(\omega)$ is holomorphic in a neighbourhood of $\omega^{\alpha}$, and the operator $\mathcal{L}^{\alpha}$ maps $L^{2}(\partial D)^{2}$ onto $\ker\mathcal{A}^{\alpha}(\omega^{\alpha},\delta)$. Let us write

[TABLE]

where ∗ denotes the adjoint operator. Then, as in [4, 5], it can be shown that

[TABLE]

where again $\langle\,\cdot\,,\,\cdot\,\rangle$ stands for the standard inner product of $L^{2}(\partial D)^{2}$.

Hence the operator $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}$ can be decomposed as

[TABLE]

Note that the third term in the right-hand side is holomorphic with respect to $\omega$.

Denote by $\alpha^{*}=(\alpha_{1},\pi)$ and $\omega^{*}(\alpha_{1})=\omega^{(\alpha_{1},\pi)}$. Using similar arguments as in [10] and the fact that each bubble is a circular disk, we can prove the following result on the shape of the dispersion relation close to $\alpha^{*}$.

Lemma A.1**.**

For a fixed $\alpha_{1}\in Y^{*}$, the characteristic value $\omega^{\alpha}$ attains its maximum over $\alpha_{2}$ at $\alpha_{2}=\pi$, i.e. at $\alpha=\alpha^{*}$. Moreover, for $\alpha_{2}$ near $\pi$, we have

[TABLE]

Here, $c_{\delta}(\alpha_{1})$ is a positive function of $\alpha_{1}$ and $\delta$.

The operator $\int_{Y^{*}}\frac{\mathcal{L}^{\alpha}}{\omega-\omega^{\alpha}}\>\mathrm{d}\alpha_{2}$ becomes singular when $\omega\rightarrow\omega^{\alpha}$. Moreover, since we want to compute the defect band inside the bandgap of the periodic problem at $\alpha_{1}$, we can assume $\omega$ is inside this bandgap. Consequently, the singularity occurs as $\omega\rightarrow\omega^{*}$. Let us extract its singular part explicitly. Denote by $\mathcal{A}^{*}=\mathcal{A}^{(\alpha_{1},\pi)},\Phi^{*}=\Phi^{(\alpha_{1},\pi)},$ and $\mathcal{L}^{*}=\mathcal{L}^{(\alpha_{1},\pi)}$. Moreover, denote by $B_{j}$ a bounded function with respect to $\omega$ in $V$. Then, by Lemma A.1, we have

[TABLE]

We therefore get

[TABLE]

for some $\mathcal{R}^{\varepsilon}(\omega)=O(\varepsilon)$ which is analytic and bounded for $\omega$ close to $\omega^{*}$. We look for characteristic values $\omega=\omega^{\varepsilon}$ of $\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)$, i.e. values such that there exists some $\Psi^{\varepsilon}\neq 0$ with ${\mathcal{M}^{\varepsilon,\delta,\alpha_{1}}(\omega)}\Psi^{\varepsilon}=0$. Expanding this equation, we have

[TABLE]

Then, multiplying by $\Phi^{*}$, we obtain

[TABLE]

Since $\mathcal{R}^{\varepsilon}(\omega)=O(\varepsilon)$, it follows from the above equation that $\langle\Phi^{*},\Psi^{\varepsilon}\rangle\neq 0$. Therefore, choose $\Psi^{\varepsilon}$ such that $\langle\Phi^{*},\Psi^{\varepsilon}\rangle=1$. This gives

[TABLE]

In order to derive a more explicit expression, we will consider the asymptotic limit of $\delta\rightarrow 0$. As in [3, 2], we have the following lemma.

Lemma A.2**.**

The following results hold:

• (i)

When $\delta\rightarrow 0$, we have

[TABLE]

which is positive for $\delta$ small enough.

• (ii)

For a fixed $\varepsilon$, when $\delta\rightarrow 0$ we have

[TABLE]

where $\psi_{\alpha^{*}}=(\mathcal{S}_{D}^{\alpha^{*},0})^{-1}[\chi_{\partial D}]$ and ${\mathrm{Cap}}_{D,\alpha^{*}}=-\langle\psi_{\alpha^{*}},\chi_{\partial D}\rangle$. For small $\varepsilon$ and $\delta$, this expression is positive for $\varepsilon<0$ if $R$ is small enough, or $\varepsilon>0$ for $R$ close enough to $1/2$.

Combining equation (A.3) and Lemma A.2, we obtain the following result

Theorem A.1**.**

Assume that $\delta$ is small enough and the pair $(R,\varepsilon)$ satisfies one of the two assumptions

• (i)

$R$* small enough and $\varepsilon<0$ small enough in magnitude (Dilute regime),*

• (ii)

$R$* close enough to $1/2$ and $\varepsilon>0$ small enough (Non-dilute regime).*

Then there exists one frequency value $\omega^{\varepsilon}(\alpha_{1})$ such that the problem (3.2) has a non-trivial solution and $\omega^{\varepsilon}(\alpha_{1})$ is slightly above $\omega^{*}(\alpha_{1})$. Moreover, as $\delta\rightarrow 0$ we have

[TABLE]

where $\alpha^{*}=(\alpha_{1},\pi)$.

Remark A.1**.**

It is easily verified that equation (4.9) evaluated for small $\varepsilon$ coincides with equation (A.4) evaluated for small $R$.

Appendix B Characterization of the effective sources for non-circular bubbles

Let now $D$ be a general simply connected domain with $\partial D\in C^{1}$. In this section, we will restrict to the case of small size perturbations $\varepsilon<0$. Define the defect bubble $D_{d}\in D$ as the domain with boundary

[TABLE]

where $\nu_{x}$ is the outward unit normal of $\partial D$ at $x\in\partial D$. We will need some results given in [5]. First, we introduce some notation. Define the mapping $p:\partial D\rightarrow\partial D_{d},\ p(x)=x+\varepsilon\nu_{x}$. Let $x,y\in\partial D$ and let $\widetilde{x}=p(x)\in\partial D_{d}$ and $\widetilde{y}=p(y)\in\partial D_{d}$. Define $q:L^{2}(\partial D)\rightarrow L^{2}(\partial D_{d}),\ q(\phi)(\widetilde{x})=\phi(p^{-1}(\widetilde{x}))$, and for a surface density $\phi$ on $\partial D$, define $\widetilde{\phi}=q(\phi)$ on $\partial D_{d}$.

We also define the signed curvature $\tau=\tau(x),x\in\partial D$ in the following way. Let $x=x(t)$ be a parametrization of $\partial D$ by arc length. Then define $\tau$ by

[TABLE]

Observe that $\tau$ is independent of the orientation of $\partial D$. The following results are given in [5], but adjusted to the case where $\varepsilon<0$.

Proposition B.1**.**

Let $k>0$. Let $\phi\in L^{2}(\partial D)$ and let $x,y,\widetilde{x},\widetilde{y},\widetilde{\phi}$ be as above. Then

[TABLE]

[TABLE]

[TABLE]

Proposition B.2**.**

Let $\phi\in L^{2}(\partial D)$ and let $x,y,\widetilde{x},\widetilde{y},\widetilde{\phi}$ be as above. Then

[TABLE]

where $\mathcal{K}_{1}^{k}$ is given by

[TABLE]

Here $\frac{\partial^{2}}{\partial T^{2}}$ denotes the second tangential derivative, which is independent of the orientation of $\partial D$.

We also state the following result which is given, for example, in [29].

Proposition B.3**.**

For $x\in\partial D$ and $k\geq 0$ we have

[TABLE]

As in Section 4, we consider the defect problem (3.1), modelled by the fictitious sources as in equation (3.3). Observe that Proposition 3.2 is valid even for the case of non-circular bubbles. To derive the analogue of Proposition 3.1, we again study equations (3.8) and (3.9), i.e.,

[TABLE]

Since $\omega$ is in the subwavelength regime, $k_{b}$ is not a Dirichlet eigenvalue. Together with the uniqueness of the exterior Dirichlet problem, we conclude that it is sufficient to consider these equations on the boundaries. Using the notation from above, this means

[TABLE]

Using the expansions (B.1), (B.2) and (B.3), we find

[TABLE]

with $q$ defined as above. From this we find that

[TABLE]

where $Q$ is the bijection $Q:\left(L^{2}(\partial D)\right)^{2}\rightarrow\left(L^{2}(\partial D_{d})\right)^{2},Q=(q,q)$ and $\mathcal{P}_{1}:\left(L^{2}(\partial D)\right)^{2}\rightarrow\left(L^{2}(\partial D_{d})\right)^{2}$.

Using the asymptotic expansions (B.2) and (B.4), we can expand the operator $\mathcal{A}_{D_{d}}$ as

[TABLE]

where

[TABLE]

Using Taylor expansion, we have that

[TABLE]

We use the Laplacian in the curvilinear coordinates defined by $T_{\widetilde{x}},\nu_{\widetilde{x}}$ for $\widetilde{x}\in\partial D_{d}$,

[TABLE]

It is easily verified that the curvatures on the two boundaries satisfy

[TABLE]

Hence we obtain

[TABLE]

In total, we have

[TABLE]

where the operator $\mathcal{P}_{2}^{-1}:L^{2}(\partial D_{d})^{2}\rightarrow L^{2}(\partial D)^{2}$ is given by

[TABLE]

Combining equations (3.6), (3.7), (B) and (B.7) we arrive at

[TABLE]

As before, we define $\mathcal{A}_{D}^{\varepsilon}=\mathcal{P}_{2}^{-1}\mathcal{A}_{D_{d}}\mathcal{P}_{1}$. Finally, we can compute this explicitly using equations (B), (B.6) and (B.7) and Proposition B.3 to obtain the following proposition, which is the analogue of Proposition 3.1 in the case of non-circular bubbles.

Proposition B.4**.**

The density pair $(\varphi,\psi)$ and the effective sources $(f,g)$ satisfy the following relation

[TABLE]

where the operator $\mathcal{A}_{D}^{\varepsilon}$ satisfies

[TABLE]